Abstract: This study has been initiated to improve the time taken by the CPF method and eliminate the convergence problem for the tested system by enforcing the General Minimal Residual (GMRES) method at the initial point at the start up, herewith referred as CPF-GMRES solver. Analyses have shown that this new algorithm not only reduced the computation time but also eliminated the convergence problem for the tested system.
INTRODUCTION
The catastrophic effect of voltage collapse, which causes total collapse to the operation of the system, raises concern to the power system utility. Studies on voltage collapse phenomenon provide ways to prevent this event from occurring and static analysis has proved to be the best approach.
Although voltage stability is a dynamic phenomenon by nature, the use of steady state analysis method is permitted in many cases in which load flow equations are used to represent the system conditions. In comparison to dynamic study, investigation on long-term voltage stability determined by steady state studies usually offers optimistic results (IEEE/PES, 2001). However, assumption shall be made where all dynamics are died out and all controllers have done their duty. Accurate dynamic simulation is also needed for post mortem analysis and the coordination of protection and control. Methods that recently available take the maximum load flow point of the power grid as the critical point of voltage collapse, corresponding to the state that the Jacobian matrix of the power flow becomes singular from the mathematical point of view (Ayasun et al., 2004).
Voltage stability problems can be analyzed based on Differential Algebraic Equation (DAE) where mathematical models are simplified according to different purposes. It is well known that when parameters are subjected to variations, the equilibria of the DAE power-system model may exhibit three local bifurcations, namely Saddle Node (SN), Hopf and singularity induced (SI) bifurcations. The SN and Hopf bifurcations, which are observed in the Ordinary Differential Equation (ODE) models of power systems, have been extensively studied and they are linked to the voltage collapse and oscillatory instabilities, respectively. On the other hand, the SI bifurcation is due to the singularity of the algebraic equations of the DAE model under some parameter variations. With an improved version based on the decomposition of parameter, dependent polynomials of the SI bifurcation can also be found. It has been shown that the SI bifurcation occurred when system equilibria encountered the singularity manifold and referred to a stability change owing to one of the eigenvalue of the reduced Jacobian matrix associated with the equilibrium diverging to infinity (Chiang et al., 1995).
In early 1990′s, researchers have proposed a new method to replace the conventional Newton-Raphson method for calculation of the power flow, which is called the Continuation Power Flow. They believed that the conventional method has a convergence problem when it reached the critical loading point due to the singularity problem (Momoh and El-Hawary, 2000). They have proven that the Jacobian matrix of the power flow analysis became singular at the voltage stability limit. The use of continuation power flow analysis can solve this problem by reformulating the power flow equations so that they remain well-conditioned at all possible loading conditions (Iba et al., 1991; Ajjarapu and Christy, 1992). This solution allows the system to be stable at stable and unstable equilibrium points.
Continuation Power Flow (CPF) method has been used to compute the system-loading margin, the given load and generation growth pattern and identifying the critical buses for the system. The continuation algorithm used is taken from a well-documented class of techniques used to find a path of equilibrium solutions of a set of nonlinear equations. From its conception, the purpose of the continuation power flow was to find a continuum of power flow solutions for a given load change scenario. The general principle behind the continuation power flow is rather simple. It employs a predictor corrector scheme to find a solution path of a set of power flow equations that have been reformulated to include a load parameter. It starts from a known solution and uses a tangent predictor to estimate a subsequent solution corresponding to a different value of the load parameter. This estimate is then corrected using the same Newton-Raphson technique employed by a conventional power flow (Flueck and Chiang, 1996, 1998).
However, after some modifications, the effectiveness and efficiency of these methods are doubtful due to the fact that the larger the continuation step, the heavier the computational effort is needed (Canizares and Alvarado, 1993). This results the delay of the analysis in the system and causes the voltage collapse during the analysis period. In some cases, the continuation method will create convergence problem when facing the sharp turning point (Canizares and Alvarado, 1993). This shows that the conventional continuation power flow algorithm could not totally eliminate singularity problem for a certain system especially a larger system.
Constraint Reactive Implicit Decoupling (CRIC) method has been tested for reactive power and voltage (QV) curves studies. The combination of the CRIC method and stopping criteria have provided an accurate result as well as significantly reduced the required computational time (Flueck and Chiang, 1996, 1998). An efficient geometric parameterization technique has been proposed by Enio et al. (2007), which introduced an additional parameter to the existing continuation power flow equation. This new method is believed could avoid the Jacobian singularity and is claimed to preserve the conventional Newton′s method and improved the existing method.
This study elucidates information on a new algorithm for continuation power flow with the use of Generalized Minimized Residual (GMRES). Several systems have been chosen to be evaluated using the new algorithm namely IEEE 14, IEEE 118, IEEE 300 and Malaysia′s 275 kV systems. Results have shown that the computational time is significantly reduced and the convergence problem has also been eliminated, thus increasing the system robustness.
GENERALIZED MINIMIZED RESIDUAL (GMRES) METHOD
GMRES is a member of the family of Krylov subspace iterative methods, which produces a sequence xk of approximations to the solution x = A-1b of linear system (Momoh and El-Hawary, 2000). In general, the Krylov subspace iterates are described by:
(1) |
where, x0 is the initial estimate of the solution to Ax = b and Kk (r0, A) is the k-th Krylov subspace:
(2) |
In particular, GMRES creates a sequence xk that minimizes the norm of the residual at step k over the k-th Krylov subspace as follows:
(3) |
At step k, GMRES applies the Arnoldi process to a set of k orthonormal basis vectors for the kth Krylov subspace to generate the next basis vector. When the norm of the newly created basis vector is sufficiently small, GMRES solves the following (k+1) x k least squares problem:
(4) |
where Hk is a (k+1) x k upper Hessenberg matrix of full rank k and
Preconditioning
The convergence rate of iterative methods is relying on spectral properties
of the coefficient matrix. A preconditioner is a matrix that can transform
the linear system into another solution, which has the same solution but
has more positive spectral properties (Barrett, 1994). As an example if
a matrix M in some way, the transformed system is:
(5) |
has the same solution as the original system Ax = b, but the spectral properties of its coefficient matrix M-1A may be more favorable. The use of the preconditioner in iterative methods will cause the solution more complicated. Therefore, it is suggested to apply during the initial setup and per iteration, there is a trade-off between the cost of constructing and applying the preconditioner and the gain in convergence speed. The initial computational step has to be amortized over the iterations or over repeated use of the same preconditioner in multi linear system. The application of preconditioners gives an amount of work proportional to the number of variables (Barrett, 1994).
A preconditioner should carefully approximate the original linear system coefficient matrix so that the iterative solver iterates for the preconditioned system converge faster than the iterative solver iterates for the original system. Also, the preconditioners should be relatively easy to use in a linear solve, since the matrix vector products now contain a linear solve step. For a well-behaved matrix A (e.g., a matrix with tightly clustered eigenvalues) iterative solver will converge quickly to a solution x that satisfies the given user defined residual tolerance. However, for an ill-behaved matrix A (e.g., a matrix with widely scattered eigenvalues) iterative solver will converge slowly to a solution x, if at all. Hence, preconditioning is essential for practical applications.
Types of Preconditioning
As stated in Barrett (1994), there are four types of preconditioning.
They are:
Point Diagonal of J
The point diagonal of J is the easiest of the four to create. The
preconditioner is a diagonal matrix so the linear solve steps is quick.
However, the point diagonal preconditioner does not represent the original
matrix well, thereby only slightly decreasing the amount of work necessary
to find a solution to the linear system.
Incomplete Factorization of J
The incomplete LU (ILU) factorization is a column oriented LU factorization
(in the complete LU factorization, this is sometimes called the kij
algorithm (SI) which discarded all fill-in elements. In fact, fill-in
is not even computed, so that both extra storage and extra floating-point
operations are avoided. The incomplete factorization of J is an attempt
to factorize the matrix J without computing any elements that would cause
fill-in. The resulting factors L and U have the same sparsity pattern
as J. This approach has two advantages. First the number of nonzero of
the factors is known beforehand, so allocation of memory is straightforward.
Second, since floating point operations are confined to the nonzero of
the original matrix, the quantity of work necessary to compute the preconditioner
is also known beforehand and lit does not depend on the ordering. This
is not the case with a complete LU factorization, where the necessary
storage grows with the amount of fill-in and so does the work.
Block diagonal of J (Pθ and Q -V blocks)
The block diagonal of J is also easy to create. The Pθ and Q-V
sub-blocks of the Jacobian are used to precondition the entire original
Jacobian. In this case, the off-diagonal blocks are ignored.
Fast Decoupled Approximation of J
An extension of the block diagonal idea is the Fast Decoupled approach
to power flow. While the Fast Decoupled Jacobian approximation can be
used on its own as a power flow solver, it can also be employed as a preconditioner
in the Newton-GMRES solution algorithm. There are situations in which
the Fast Decoupled approach may not be robust enough to solve the power
flow equations efficiently; hence the Newton approach is required.
In this study, the block diagonal of J preconditioner has been selected as it was found to be the right preconditioner for this research.
CONTINUATION POWER FLOW
The Continuation Power Flow method is a powerful tool to detect the voltage collapse point because of its ability to generate the whole PV curve without having the singularity problem. The Jacobian matrix from power flow (Kundur, 1993):
(6) |
becomes singular at the voltage stability limit. Conventional power flow algorithms are subjected to the convergence problems at operating condition near the stability limit. To overcome this problem, continuation power flow algorithm has been proposed. This algorithm will solve the problem by reformulating the power flow equations and by ensuring the system remains in well-conditioned at all possible loading conditions. In other words, the algorithm gives the solution of the power flow problem for stable and unstable equilibrium points (Kundur, 1993). This continuation power flow uses an iterative process involving predictor and corrector steps.
The use of the predictor is to find an approximation for the next solution. Usually tangent, first-order polynomial, or zero order polynomial predictor is employed. The parameterization is a mathematical way of identifying each solution on the solution curve. Parameterization augments the system of power flow equations. The corrector is usually an application of Newton method to the augmented system of equations. The step length control can be done by optimal fixed step length or by adaptive step length control.
The Predictor
In this part as stated in (Cutsem, 2007) the predicted solution is predicted
in order to find the right solution. It started with a simple power system
equation. The power flow equations are the function of the bus voltages
and bus angles and the bus injections:
(7) |
Replace the power flow equations with the functions of λ.
(8) |
| and |
(9) |
Each set of partial derivatives are evaluated at a particular set of operating conditions. If the power flow equations are linear with the three sets of variables in the region between the old solution and the (close) new one, the equation:
(10) |
| Thus, |
(11) |
Unknown λ is added to the power flow problem without adding a corresponding equation, i.e., in G(θ,V,λ) = 0. There are N equations and N+1 variables, so that in equation, the matrix [Gθ GV, Gλ], has N rows (the number of equation being differentiated) and N+1 columns (the number of variables for which each equation is differentiated).
(12) |
Defining the step size as σ,
(13) |
The tangent vector is given as,
(14) |
This vector provides the direction to move in order to find a new solution (i+1, p) from the old one (i). Suppose that the k-th parameter in the tangent vector is set to be ±1.0. Then,
(15) |
The Corrector
The predicted solution done in the predictor method is corrected for
getting the correct solution (Cutsem, 2007). There are two techniques
available and widely been used by researcher today namely the perpendicular
intersection and the parameterization techniques.
Parameterization
Fixing the continuation parameter and then solving the power flow
equations perform the corrector step.
(16) |
yk(i+1) is the continuation parameter where it corresponds to the k-th element in the tangent vector. The parameter η is the value to which yk is set, which is the value found in the predictor step solution. By using the modified Newton Raphson method, the solution of using the parameterization technique can be found. If there is no convergence, the step size (σ) is cut by half of the previous value and procedure is continuously running from the beginning of the predictor method.
Perpendicular Intersection
The intersection between the power flow equations (the PV curve) and
a plane that is perpendicular to the tangent vector can be found using
the Eq. (17).
(17) |
The next step is to used the Newton-Raphson method for solving the intersection in Eq. 17 (requires only 1-3 iterations since it has a good starting point) and again, if there is no convergence, the step size (σ) is cut by half and the procedure will keep continuing from the beginning of the predictor method.
The difference between the two techniques can be seen from their equations, which lead to different types of approaches. Parameterization technique is done by fixing one of its continuation parameters, either the voltage or the loading parameters whilst the perpendicular intersection technique applies the principal of finding the perpendicular intersection in applied mathematics. The intersection between power flow equation (PV curve) and a plane that is perpendicular to the tangent vector will then be determined. In addition, this technique requires no additional equation in the Jacobian matrix and involves with product to the tangent vector (Cutsem, 2007). The critical point can be achieved from the equation of the tangent vector itself. The dλ in the tangent vector becomes negative, because the loading has reached a maximum point and begun to decrease. The continuation parameter is selected from among λ and the state variables in y. This will be the parameter that has the largest element in the tangent vector.
CPF-GMRES SOLVER
Basic algorithm of GMRES was programmed in MATLAB and tested on several large sparse matrices for evaluation. Similar procedure was also applied to the CPF method. Both methods were then combined to form the CPF-GMRES method or referred as a CPF-GMRES solver.
Figure 1 shows the algorithm for implementing the CPF method. Initial point was determined by solving the base case load flow. Once the step size σ, was chosen, the predicted solution was then computed. From the predicted solution, the continuation parameter, ek was chosen and the corrected solution was calculated using the modified Newton Raphson method in which additional equation is inserted in the Jacobian matrix where η is defined as an appropriate value for the kth element of y.
The GMRES method is designed to solve nonsymmetrical linear systems. The most popular form of GMRES is based on a modified Gram-Schmidt orthonormalization procedure. GMRES methods for power flow has been proposed by several researchers and they have proven that power flow with GMRES method are capable to overcome the shortcomings of previous methods such as Gauss-Seidal, Newton and Fast Decoupled method (Adam, 1996; Flueck and Chiang, 1996, 1998). Equation for power flow is solved iteratively using matrix-vector multiplications without factorization (Adam, 1996; Flueck and Chiang, 1996, 1998).
| Fig. 1: | Algorithm for CPF method |
RESULTS AND DISCUSSION
Four systems have been tested to prove the reliability of the new algorithm namely, IEEE 14 bus system, the IEEE 118 bus system, the IEEE 300 bus system and 275kV Malaysia′s system. The results were then compared with the standard CPF method.
Table 1 shows that the performance of the CPF-GMRES method is surpassed compared to the standard CPF. This is attributed to the improvement made on the algorithm, which eliminates the matrix factorization and contribution from spares matrix-vector multiplications in GMRES as mentioned previous.
Figure 2 (a, b) and Fig. 3 (a, b) show the resulting curves at bus 10 for IEEE 14 bus system. The CPF-GMRES method produces similar result as standard CPF method. Both methods give similar voltage collapse point. This illustrate that the CPF-GMRES method is able to generate full π-V curve as the standard CPF does with the convergence and robustness properties of the CPF-GMRES also have been proven accordingly.
| Table 1: | Result comparing the CPF method and the CPF-GMRES method (IEEE 14) |
| Fig. 2: | (a) Voltage at bus 10 (IEEE 14) using standard CPF method and (b) Voltage at bus 10 (IEEE 14) using CPF-GMRES method |
| Fig. 3: | (a) Voltage at bus 23 (IEEE 118) using standard CPF method and (b) Voltage at bus 23 (IEEE 118) using CPF-GMRES method |
Table 2 shows the comparison made on CPF-GMRES solver and standard CPF method for IEEE 300 bus system and 275kV Malaysia′s system. The sharp turning point and step cutting technique used in standard CPF (Barrett, 1994) creates problem for the solver to converge. It may produce the right answer; however, if the step is cut below the convergence tolerance of the Newton solver, the program will take this incorrect value as a solution. The predictor, when applied to a sharp nose curve close to bifurcation, yields a large step. This creates convergence problems for the corrector part of the method, since there is no crossing with the bifurcation branch of the equilibria. In contrast, the CPF-GMRES solver yields a better performance in completing the whole solution curve with the total time of 14.1505 and 10.0029 sec for IEEE 300 bus system and 275 kV Malaysia system respectively. The solution curves for IEEE 300 bus system is shown in Fig. 4.
The results, as expected, show a very good agreement as obtained using the IEEE data systems. The solution curve in Fig. 5 clearly shows that the standard CPF faced a convergence problem during the analytical analysis whilst with the new solver; the solution curve can be drawn successfully, as shown in Fig. 6. This curve also shows that the new CPF-GMRES solver able to trace the upper and lower equilibrium point as well as voltage collapse point.
| Table 2: | Result comparing the standard CPF and the CPF solvers (IEEE 300) and 275 kV Malaysia system |
| Fig. 4: | Voltage at bus 283 (IEEE 300) using CPF-GMRES |
| Fig. 5: | Voltage at bus 64 ((TNB 275kV) using CPF method |
| Fig. 6: | Voltage at bus 64 ((TNB 275kV) using CPF-GMRES method |
CONCLUSION
The CPF method is a powerful algorithm to overcome the singularity problem in Newton Raphson method. CPF method offers the continuation of the power flow solution and has the ability to detect the bifurcation point or known as the manifold point. However, the drawback of the CPF method is the length of the time taken to get to the solution. Therefore the CPF-GMRES had been proposed to minimize the time taken to get to the solution. Results have proven that the method is reliable in term of the accuracy and speed if compared to the standard CPF method. In addition, the results show a very good agreement as obtained using the IEEE bus data systems.
This study shows that for a large system, the conventional CPF face a convergence problem in getting the solution. This was evidence in the studies done on the IEEE 300 bus system and 275kV Malaysia system. The new method has been applied to the two systems to prove that this method is capable in solving the convergence problem and it has been found to be successful in eliminating the convergence problem.
ACKNOWLEDGMENTS
J. Jasni acknowledges a scholarship from Ministry of Science, Technology and Innovation (MOSTI), Malaysia. The authors are grateful to Tenaga Nasional Berhad, Malaysia for providing the data for this research.