Fault Current Constraint Transmission Expansion Planning Based on the Inverse Matrix Modification Lemma and a Valid Inequality

In the transmission expansion planning (TEP) problem, it is challenging to consider a fault current level constraint due to the time-consuming update process of the bus impedance matrix, which is required to calculate the fault currents during the search for the optimal solution. In the existing studies, either a nonlinear update equation or its linearized version is used to calculate the updated bus impedance matrix. In the former case, there is a problem in that the mathematical formulation is derived in the form of mixed-integer nonlinear programming. In the latter case, there is a problem in that an error due to the linearization may exist and the change of fault currents in other buses that are not connected to the new transmission lines cannot be detected. In this paper, we use a method to obtain the exact updated bus impedance matrix directly from the inversion of the bus admittance matrix. We propose a novel method based on the inverse matrix modification lemma (IMML) and a valid inequality is proposed to find a better solution to the TEP problem with fault current level constraint. The proposed method is applied to the IEEE two-area reliability test system with 96 buses to verify the performance and effectiveness of the proposed method and we compare the results with the existing methods. Simulation results show that the existing TEP method based on impedance matrix modification method violates the fault current level constraint in some buses, while the proposed method satisfies the constraint in all buses in a reasonable computation time.


Introduction
Because electricity is difficult to store economically, and it can only be delivered through transmission lines, transmission networks must be expanded in proportion to the increase in electricity demand. Even though the expansion of a transmission network provides positive effects, such as reduction of congestion, an increase in system reliability, and a reduction of transmission loss, a side effect of fault current increase (especially near large metropolitan areas) occurs due to the decrease in the impedance of the overall system [1].
Most of the widely used methods solve this problem by splitting buses, installing additional equipment such as a superconducting current limiter, or upgrading the circuit breaker with a higher breaking capacity. However, one of the more fundamental solutions to solve this problem is to consider the fault current level during transmission expansion planning (TEP), which is the main topic of this paper.
In general, the capacity of a circuit breaker is chosen to be greater than the three-phase line-to-ground fault current, which is calculated as the pre-fault voltage divided by the Thévenin equivalent impedance. In power system optimization problems with short time horizons, such as the optimal power flow or security constraint unit commitment, it can be assumed that the network configuration is fixed, and the fault current can be considered as linear constraint or upper bound on the bus voltage magnitude [2][3][4]. However, in a problem with a relatively long time horizon, it becomes mathematically challenging to consider the fault current in the optimization problem, because the bus impedance matrix is continuously changing.
Research on the TEP problem that considers fault current has been published only recently, and there are two reasons why the fault current has not been studied in the TEP problem in the past. The first is that the commercially available circuit breaker has a sufficient breaking capability margin, compared to the fault current possible in the system, and the second is that the computational burden has been too large to account for the fault current in the TEP problem. There is limited research on the fault current in the TEP problem [5][6][7][8]. The most important feature of these studies is summarizing how to reduce the calculation burden due to updating the system matrix (the bus impedance matrix), which is essential for calculation of the fault current.
For example, in order to overcome this problem, previous research used a metaheuristic algorithm such as the genetic algorithm (GA), removed infeasible cases by applying the Benders cut, or simplified the problem by linear approximation. However, the GA algorithm has a problem in that it cannot reproduce the same results due to the inherent randomness of the algorithm, not to mention its slow convergence [5]. Additionally, methods based on Benders decomposition must solve the problem repeatedly until infeasible cases do not occur, and hence, the number of iterations required to find the optimal solution increases exponentially as the number of candidates for the new transmission line increases [6]. In the linear approximation method, even though the method can find the solution in a relatively short period of time, there is a disadvantage in that the calculation result of the fault current may differ significantly from the actual value due to the model error brought about by the linearization process [7,8]. As a result, most of the existing methods have a problem in that the calculation time is too slow or the accuracy is too limited to be applied to an actual large-scale power system. It should be noted that our study is confined to the balanced fault such as three-phase line-to-ground fault. For power systems where the unbalanced faults are dominant, such as distribution systems, special methods are employed to calculate power flow and fault current, as described in [9,10].
Therefore, the main purpose of this paper is to propose a more exact and efficient method to solve the TEP problem with a fault current level constraint based on the inverse matrix modification lemma (IMML) and a valid inequality. The proposed method consists of the following two steps.
The first is to efficiently calculate the diagonal terms of the bus impedance matrix, which is the Thévenin equivalent impedance of the system [11], and this bus impedance matrix should be recalculated repeatedly during the optimization process. In the existing research [7], a linearized and simplified version of the nonlinear impedance update equation is used to formulate the resulting optimization problem into mixed integer linear programming (MILP) form. The authors of [7] concluded that the normalized root mean square (RMS) error of about 2-3% due to linearization is a minor problem from a practical standpoint. However, in some power systems, such as those in South Korea, where the margins between the fault current and the circuit breaking capability are very narrow, more exact calculation of the fault current is required. Another problem with the method used in [7] is that it cannot guarantee the accuracy of the fault current calculation results at other buses to which the new transmission lines are not connected directly.
In this paper, in order to avoid these problems, the updated bus impedance matrix is calculated from the inversion of the bus admittance matrix. As will be explained in the following sections, this method does not require the nonlinear update equation for calculating the bus impedance matrix, and therefore, the overall optimization problem can be formulated as a MILP optimization problem without linearization. Additionally, the fault currents at the other buses to which the added transmission line is not directly connected can also be calculated without error.
The disadvantage of this method is obviously that the inverse of the bus admittance matrix must be calculated repeatedly in the search for the optimal solution to the TEP problem. In order to mitigate the computational burden of computing the inverse of the bus admittance matrix at every step to find the optimal solution to the TEP problem, we apply the inverse matrix modification lemma (IMML) and the valid inequality.
Due to the nature of the TEP problem, only a part of the transmission network is slightly changed during the search for the optimal solution. Therefore, it has a structure that is well suited to applying the well-known IMML, which can help to calculate the inverse of the bus admittance matrix efficiently [12]. A linear update equation for the diagonal terms of the bus impedance matrix will be explained in detail in the later sections.
Even though the time to calculate the updated impedance can be reduced greatly by using the IMML, the overall computation time to find the optimal solution to the TEP problem still becomes intractable as the number of candidates for new transmission lines increase. Thus, in this paper, we have devised a new method based on a valid inequality to reduce the feasible area of the problem. The valid inequalities are generated every search step during the optimization and added as additional inequality constraints to branch-and-cut (B&C) algorithm to reduce the feasible search region. As shown in the simulation section, the combination of the IMML and the proposed valid inequality can find the optimal solution to the TEP problem with consideration of the fault current very accurately in a reasonable calculation time.
The rest of this paper is organized as follows: Section 2 describes the meanings of all the parameters and variables used in the paper. Section 3 explains the MILP formulation of the TEP problem with the fault current level constraint. The method used for the solution in this study, based on the B&C algorithm, and the suggested valid inequality method are explained in Section 4. The proposed method is applied to IEEE Reliability Test System with 96 buses (RTS-96), which is frequently used in studies for the TEP problem. The simulation results and the concluding remarks are summarized in sections 5 and 6, respectively.

Objective Function
The objective function used in this paper has a form similar to that used in the general TEP problem, which is composed of generator operation cost, transmission line construction cost, and salvage value of the transmission lines, as follows: where (⋅) is the linearized generation cost function given as follows: which can be further simplified into the following equation:

Conventional Method for the Bus Impedance Update
Generally, the three-phase line-to-ground fault current can be calculated by the following equation: where and ℎ are the pre-fault voltage and Thévenin equivalent impedance at bus , respectively. Since Thévenin equivalent impedance is the same as the corresponding diagonal element of the bus impedance matrix if the pre-failure voltage is assumed to be 1.0 (pu), the three-phase line-to-ground fault current constraint in the TEP problem can be expressed as follows [13]: where , , is the (i,j)th element of bus impedance matrix in year t. In the existing study, the following update equation is used to calculate each element of the bus impedance matrix after the transmission lines are added [7,8]: , ≥ 1 where and are the coefficients for linearization, which can be calculated using the following equations: In order to improve the computational speed in the process of applying the above equations, the existing studies [7,8] assume that the off-diagonal elements of the bus impedance matrix are negligibly smaller than the diagonal elements (� , � ≪ � , �), and the diagonal elements of the bus impedance matrix before the update are significantly smaller than the impedance of the new transmission lines that will be connected to the buses (� , � ≪ � , �). Hence, the values of and in Equation (8) are all assumed to be 1. This method has an advantage in that the computational burden for calculating the fault currents can be largely reduced, however, there is a problem in that it is impossible to consider the changes in the fault currents occurring on buses to which the new transmission lines are directly connected. Additionally, as the number of new transmission lines increases, the calculation error is cumulative.

Derivation of the Bus Impedance Update Equation Based on the IMML
In this paper, we adopt a method to directly calculate the bus impedance matrix by inversion of the bus admittance matrix of the new system with the added transmission line, as follows: It is obvious that the calculation time will significantly increase due to the curse of dimensionality if the system becomes larger. Therefore, it is essential to calculate the inverse matrix more effectively to incorporate the fault current constraint (or Equation (5)). In order to find the optimal solution to the TEP problem, a new transmission network configuration for the next search step is constructed by adding one new transmission line to the network configured in the previous search step. The bus admittance matrix of the next search step will have almost the same value, and only small number of elements in the admittance matrix is changed from the bus admittance matrix of the previous step. Therefore, the following inverse matrix modification lemma can be applied in a very straightforward manner to acquire the updated bus impedance matrix.
The so-called Woodbury equation of the IMML can be described as the following generic equation [12]: If , the variation of the admittance matrix , can be decomposed into the product of three matrices, or ⋅ ⋅ , the new admittance matrix, , can be expressed as follows: where is the variation of matrix Y (m×m), is the connection matrix (n×m), is the number of total buses and is the number of modified buses. Here, time index is omitted for simplicity. By applying the Woodbury equation of the IMML to Equation (12), a new impedance matrix, , can be obtained as follows: where is the impedance matrix of the previous optimization step. In order to eliminate the unnecessary inverse operation, Equation (14) is replaced with the following equation: Here, we apply an extra assumption that only one new transmission line can be added to the existing bus in order to derive a more compact formulation for Equation (16). In case that more than one transmission line should be added to a bus, the above assumption is still valid since we can insert a dummy bus between the new transmission line and the existing bus with small reactance between the buses. Therefore, most of the elements in matrix are not changed except for the corresponding elements of the pair of the buses to which the new transmission line candidate is connected. Thus, can be expressed as the following equation: where , and , are [(m − 2) × (m − 2)] and [2 × 2] matrices, respectively. Furthermore, the following equations can be obtained by substituting Equation (18) into Equation (16) with appropriate partitioning: , � + , , � + , , , = 0 , ( + , , ) + , , , = 0 (22) Without loss of generality, by appropriate indexing and simple matrix manipulation, we can define subscript A as the set of buses where the transmission line is constructed, and subscript B as the set of buses where the transmission line is not constructed, which means that all elements of , in Equation (19) will be zero. Then the following equations can be easily derived from Equation (20) through Equation (23): , ( + , It should be noted that both Equation (24) and Equation (25) should be regarded as nonlinear, because they contain terms in the production of two variables: and . In this paper, we apply the following algorithm, which has the same effect as the above equations, so that the resulting optimization problem becomes a MILP formulation. If = 1, then: And if = 0, then: In the above equations, and denote the real and imaginary parts of the product of matrix and matrix , respectively. All logic expressions in the above are implemented using the indicator constraint option of the CPLEX solver incorporated in the General Algebraic Modeling System (GAMS) [14]. Once all the elements of are calculated using Equations (27)−(29), they can be entered into Equation (12), and the equation to calculate the diagonal element of the new bus impedance matrix (or the Thévenin equivalent impedance) is derived as follows: Generally, the reactance is more dominant than the resistance of the transmission line, thus, it can be assumed that { , } ≫ { , } without significant loss of generality. Hence, the constraint given by Equation (5) can be rewritten as follows:

Other Constraints
The following constraints are incorporated into the model to account for the technical limitations of the power system and its components.

Power Flow Equation
The following typical DC power flow equations are included in the model [13]: Equation (32) is the DC power flow equation for the existing transmission lines, and Equation (33) is for the candidates to be newly constructed. For the latter, in the search process, constructing ( , = 1) and not constructing ( , = 0) are expressed separately.

Node Balance Equation
The demand and supply condition at each bus is represented as the following equality constraint:

Limit on Generator Output
All the generators should satisfy the following inequality constraints for the limits on the active power output: In general, in the TEP problem, the on/off status of generator is not considered due to the computational burden. In this paper, for slight improvement of the accuracy of the generation cost, the on/off status of generator is included in the model. However, if the overall computation time is crucial during the simulation, Equation (35) can be simplified to 0 ≤ , , ≤ without significant problem.

Transmission Limit
The following inequality constraints should be satisfied to limit the apparent power flow through the transmission lines:

Bus Voltage Angel Limit
The following inequality constraint is included to limit the bus voltage angles:

New Transmission Line Constraint
In this study, we assume that both existing and newly constructed transmission lines will be maintained until the end of the planning period once they are constructed. The following constraints are included in the model for this purpose:

Solution Methods
The TEP problem with the fault current constraint formulated in the previous section is a typical MILP optimization problem. One of the most commonly used methods to find the optimal solution of the MILP problem is the branch-and-cut algorithm. In general, the B&C algorithm is known to be able to obtain the optimal solution for a relatively small-scale MILP problem, but as explained in the simulation section, the curse of dimensionality still occurs in the TEP problem when considering the proposed fault current constraint. As a result, when the scale of the system increases, the optimal solution cannot be searched for within a reasonable period of time.
In this section, we propose a valid inequality method that can help efficiently find the optimal solution to a given TEP problem. The valid inequality method is one of the ways to improve the search speed by removing the search space from the feasible region where there is no possibility of the solution in a given optimization problem. For example, valid inequalities are effectively used to improve the performance of the B&C algorithm to solve MILP problems such as the traveling salesman problem [15], unit commitment [16], and the oil pipeline design problem [17]. Since most of the valid inequalities are very problem-specific in nature, the valid inequalities used in the previous studies cannot be applied to the TEP problem. Hence, this study proposes a new valid inequality that is well suited to the proposed TEP problem, and it shows through simulation results that the proposed valid inequality works efficiently. Next, we briefly describe the B&C algorithm and then explain the valid inequality suggested in this study.

Branch and Cut(B&C) Algorithm
The B&C algorithm is a combination of the branch-and-bound (B&B) algorithm and the cutting plane method. Figure 1 shows a flow diagram of the B&C algorithm used in this study [18].

Valid Inequalities
As mentioned above, it is assumed in this paper that the existing transmission lines are not decommissioned during the planning period. We believe that this assumption is reasonable because the existing transmission lines are rarely decommissioned even in an actual situation, and in most cases, the aged lines are maintained through reinforcement works.
The proposed model includes two binary variables, one related to transmission line construction ( ), and the other to generator operation status ( ). If a combinatorial set of does not satisfy the fault current constraint, a feasible solution does not exist, whatever values are assigned to . Additionally, in general, as more transmission lines are built, the fault current at each bus tends to increase.
Having said that, we can assert that at least one of the new transmission lines included in the expansion plan should be removed from the plan if a given combination of does not satisfy the fault current constraint. Therefore, the following valid inequality constraint can be added to the mathematical model if the fault current constraint or Equation (5) is not satisfied: where round() is the round-off operator. The valid inequality is added as a new constraint to the cutting plane. In this study, the above valid inequality is implemented using the branch-and-cut-and-heuristic (BCH) of GAMS [14].

Simulation Results
In order to verify the performance of the proposed method, several simulations were performed on the IEEE RTS-96 system, which is frequently used for the TEP problem [19]. The studied power system is depicted in Figure 2. Additionally, Figure 3 shows the load duration curve used in the study, where the durations of load blocks are 1752 hours for a 100% load, 5256 hours at 80%, and 1752 hours at 40%. The growth rate of annual electricity demand and the discount rate are arbitrarily set at 10% and 7%, respectively. The base power and the upper limit of the fault current in each bus are set to 100 MVA and 10 kA, respectively. Other technical parameters used in the simulations are listed separately in the Appendix A.  The proposed model was implemented using the CPLEX solver, which is one of the most frequently used MILP solvers available with the General Algebraic Modeling System (GAMS) [20]. In short, the mathematical formulation of the proposed method is modeled using GAMS to find the optimal solution to the TEP problem and the result from the GAMS model is entered into the PSS/e ASCC Psoftware package [21] to calculate the exact fault currents. For comparison, the existing method based on the linearized impedance update equation [7] is also simulated. The simulation study consists of the following three cases: Case 1: TEP without fault current constraint Case 2: TEP with both the fault current constraint and the valid inequality (proposed method) Case 3: TEP with the fault current constraint based on the linearized impedance update equation (existing method) [7] Tables 1-3 show the optimization results for the three cases. Figure 4 shows the fault current at each bus calculated by the PSS/e ASCC software package [21] using the transmission expansion plans acquired from the results of the three cases. In the figure, the Base case is the fault current before the transmission expansion planning starts and is included in the figure to show that the fault current increases monotonically when transmission lines are added to the system. Table 4 shows the calculation results for the fault current for each year on the bus with a fault current exceeding the upper limit (10 kA).   In Figure 4, the result of Case 1 yields fault currents that exceed 10 kA in the buses of 113, 203, 209, 215, 216, and 218 (blue and red circles). In view of the fault current, the biggest problem of Case 1 is that some of the new transmission lines are constructed on the bus 215, 216, and 218, which have the largest increase in fault currents. It is an obvious result because the fault current level constraint is not considered in the optimization model in Case 1.
On the other hand, any of the fault currents of Case 2 does not exceed the limit of the failure level and the new transmission lines are constructed on the buses 112, 211, and 224 with relatively low increase in fault current.
In Case 3, it can be observed that the fault current at bus 209 exceeds the upper limit as shown in Figure 4 (red circle) and Table 4, and the cause of this problem is due to the linearization of the bus impedance update equation explained in Section 2.2.1. This result can be seen to be due to a rather extreme setting, but we believe that there is always the possibility of this problem occurring in a real system considering the recent circumstances of the power system where the margin between the fault current level and the circuit breakers' breaking capacity decreases gradually. One thing that should be noted in the result of Case 3 is that the fault current violation happened at bus 209, which

Fault Current (kA)
is not connected to any of the new transmission lines as shown in Table 3. The reason for this problem is that when the impedance matrix is updated in Case 3 using Equations (6)- (9), the change of the impedance values between other buses that are not connected to the new transmission lines are ignored in order to reduce the calculation time. This problem does not occur in the proposed method (Case 2) because fault current is calculated using exact impedance matrices that reflect the change in the impedance of all buses.
The objective value and the calculation time for each case are shown in Table 5. As shown in the table, Case 2 (the proposed method) takes longer to find the optimal solution than other cases, and the objective function (construction, operation costs, and salvage values) shows that the proposed method is slightly improved compared to Case 3, but the differences are not significant considering the circumstances of the TEP problem. Therefore, the main advantage of the proposed method is that the TEP problem considering the fault current constraint can be solved more precisely by formulating it in MILP form without linearization.  One last analysis of the proposed method is the impact of the valid inequality on the time to convergence at the optimal solution. Figure 5 shows the results of comparing convergence times with and without using the valid inequality. As shown in the figure, the duality gap converges at 0 in about 267 seconds, but it takes longer than 1000 minutes when the valid inequality is not used. The results of the transmission expansion plan in both cases are exactly the same, and therefore, the effect of the valid inequality is clear.

Concluding Remarks
In this paper, we propose a new method to solve the TEP problem with more exact consideration of fault current level constraint. In the existing TEP studies, approximation errors occur due to the linearization of the impedance update equation, which is formulated as a nonlinear equation to calculate the fault current of a system. In this paper, we adopt a method of eliminating the approximation error by obtaining the updated impedance matrix directly from the inversion of the admittance matrix using the IMML.
This study also proposes a novel valid inequality method to reduce the computation time to find the optimal solution to the TEP problem with fault current level constraint. The concept of valid inequality is a way of reducing the overall computation time by adding an inequality constraint at every computation step that removes the infeasible region of the optimization problem from the search space. In this paper, the valid inequality is constructed using the characteristic that the fault current in the system generally increases due to the increase of the parallel circuit when new transmission lines are added. As shown in the simulation results, this valid inequality proves to be capable of reducing the overall computation time significantly.
The contribution of this study can be summarized as follows. First, the TEP problem considering the fault current constraint is formulated as a MILP optimization problem without linearization. Second, the convergence of the MILP optimization problem is greatly improved by reducing the feasible area through the valid inequality.
In addition to the TEP problem considering the fault current, the proposed method can be applied to studies on the economical evaluation problems of the newly discussed technologies in solving the fault current problem, such as the optimization of bus separation, the superconducting fault current limiter, and the usage of medium voltage direct current.