LOAD-FLOW SOLUTIONS FOR ILL-CONDITIONED POWER SYSTEMS BY A NEWTON-LIKE METHOD

3648 LOAD-FLOW SOLUTIONS FOR ILL-CONDITIONED POWER SYSTEMS BY A NEWTON-LIKE METHOD G. S.C. Tripathy Dept. of Electrical Indian Institute of New Delhi - 110016 Durga Prasad Engineering Technology India ABSTRACr In this paper mathematician K.M. Brown's method is used to solve load-flow problems. The method is particularly effective for solving of ill-conditioned nonlinear algebraic equations. It is a variation of Newton's method incorporating Gaussian elimination in such a way that the most recent infonnation is always used at each step of the algorithm; similar to what is done in the Gauss-Seidel process. The iteration converges locally and the convergence is quadratic in nature. A general discussion of ill-conditioning of a system of algebraic equations is given, and 'it is also show by the fixed-point formulation that the.-proposed method falls in the general category of sucessive approximation methods. Digital computer solutions by the proposed method are given for cases for which the standard load-flow methods failed to converge, namely 11-, 13and 43-bus ill-conditioned test systems. A comparison of this metlhod with the stanidard load-flow methods is also presented for the well-conditioned AEP 30-' and 57-bus systemns. planning. admittance matrix has gained-widespread popularity because of its quadratic convergence characteristics. it has However, limitations in small-core- computer applications where the weakly convergent Gauss-Seidel method is generally resuitble Extnsie meory methd isgenrall more suitable. Extensive meiry in ticalculations orered vatedthe xplotatio of parsty wih elimination and skilful progranuning in the Newt,on n~th od f3]. Recently, advantage has been taken of loose physical interaction between MW and MVAr flows by mathematically decoupling the MW-e and MVAr-V-calculations [4]. Load-flow calculations are performed in system An early approach was the Gauss-Seidel iterative method [1] using the nodal-admittance-method and this was further improved by using the nodal-impedance matrix method [2,9].- Newton's method [3] using a nodal raequiremt lare pow system Despite the substantial progress there are still some difficulties with som'eof the above methods [5,6]. Features which cause instability and divergence in load-flow solution methods are: ' . of (i) position .: the reference slack bus (ii) existence of negative line reactance (iii) certain types of radial systems (iv) high ratio of long-to-short line reactance for lines terminating on the same bus Power Engineering Society for presentation at the IEEE PES 1982 Winter Meeting, New York, Nlew York, January 31- 82 WM 021-4 A paper recomended and approved by the3 IEEE Power System Engineering Committee of the IEEE February 5, 1982. Man^uscript s-ubmitted January 21, 1981; made available for printing Ootober 23, 1981. - . - . x=vectorfukwnaibles IEEE Transactions on Power Apparatus and Systems, Vol. PAS-101, No. 10 October 1982 G.S. Hope O.P. Malik Dept. of Electrical Engineering University of Calgary Calgary, Alta. T2N lN4, Canada (v) choice of acceleration factors. INTRODUJCIION Networks which have the above features are described by a system of ill-conditioned nonlinear algeIn other words, a small change in braic equations. parameter produces a large'change in the solution. In this paper a new algorithm using K.M. Brown's method [7,8] is presented. It is useful when formlation of the load-flow problem results in ill-conditioned equation and the nodal admttance matrix. This is a quadratically convergent Newton-like method based upon Galussian elimination. In Brown's method, each equation is expanded in the approximate Taylor's series; however, the most recent information available is immediately used in the construction of the next fuction This procedure is similar to the procedure argment in the Gauss-Seidel process for the solution of nonlinear sets of equations. This contrasts sharply with Newton's method in which all equations are treated simIt has been proven that Brown's method ultaneously [7] is not equivalent to Newton's method. An optional feature of Brown's method is that the partial derivatives of the equations are replaced by their first difference quotient approximations similar to the discrete is shown by the fixed-point falls Newton's [13] that It the proposed load-flow methodformulation method. into the general category of successive approximation methods such as Gauss-Seidel proposed Brown's load-fl tal comuter results of the and Newton's method. Digi- algorithm are given for 11-, 13- and 43-bus ill-condipower 1,4] naey Gus-id,Nwotionedmtos systems [6,12], forwhichstandardloadflw ~~~~~~~~~~~and shw proble they systems [9]. mor standard load-flow methods [1,3,4,5] test systems, for 11-, 13- and 43-bus ill-conditioned is presented -as well as for1 AEP 14-, 30- and 57-bus well-conditioned test LIST OF SYMBOLS k = condition number [J] = Jacobian matrix A comaisoneo poor cn rg problems, they Brown's load-flow algorithm with the the proposed show poor convergence. A,co, arlson of Although method [5] work Denmead method [2] verge. Fletcher-Powell the Brameller andfor ill-conditione'd Raphson and fast decoupled algorithms failed to con- =total ~~~~~~~~~N number of ukown variables x° = initial assumption to the vector of knowns f= ith fuction of system of equations f gj = ith function of system of equations , x* = solution vector F = iteration function of Brown's method fi = af /ax. partial derivative of function fi w.r.t. x. 3 gi = ag./ax. partial w. r.t. x . D deriv5ative of fuction gi T = transpose of matrix (superscript) n = iteration cout (superscript) 1982 IEEE ~~~~~0018-9510/82/1000-3648$00.75 ©B 3649 P. = ith bus real power mismatch 1. = ith bus reactive power mismatch Q. = ith bus reactive power mismatch Pi(cal) = ith bus calculated real power QQ(cal) = ith bus calculated reactive power P i(sch) = ith bus scheduled real power Q i(sch) = ith bus scheduled reactive power P. = ith bus net real power 1 Qi = ith bus net reactive power V .,V = voltage magnitude at 1 bus i and j i e., O.= voltage phase angles at bus i and j Y..= G. +jB. = (i,j)th element of bus admittance matrix [H] = transverse upper triangular matrix formulated in forward part of Brown's method E = specified tolerance = machine tolerance Appendix K.M. Brown's method is described in detail in 1. The application of this method to a system of nonlinear equations is also demonstrated in Appendix 1. Brown's method can be formalized by writing terms of the iteration function F = (fl, N) used, beginning with a starting guess x°, to form successively: stable. f**. x = P(xn) n = 0, 1,p2 .( 2, The iteration fuction, F, in this method is given by n xn+1 =F( ) 0 3 F.(x1, i (F n * x x1) = Xi XN) Xi i-l - , J /gN-i+l,x. (gN-i+l~~~~~~~~~~~~~~~,x. x1 ' i 1, 2, .,N (4) (5) where I = 0 whenever m > n j=m and gi is defined as follows: g= f1(x1, x2, ..., A= eigenvalue. xN) XN-i+l LN-i+2 * LN) N**, *L** LN, i ILL-CONDITIONED SYSTEMS AND EIGENVAL[E ANALYSIS A computer formulation of a problem- is defined to be ill-conditioned if computed values are very -sensitive to small changes in input value. A matrix may have some eigenvalues which are very sensitive to small changes in its elements while others are comparatively insensitive. It is convenient to have some number called a 'Condition Number' which,defines the condition of a matrix with respect to the computing problem. Ideally, it should give. some overall assessment of the rate of change of the solution with respect to changes in coefficients. It may be recalled that the size of the condition number k([J]) of a Jacobian matrix [J] is defined to be IIJII IIJ-11 . It'gives a good indication of the sensitivity of [JI- to small perturbation in [JM. It is also called the spectral condition nunber because of its dependence on the spectral norm 1 u 1*. Thle condition number of a symmetric positive definite matrix [JTJ], whose eigenvalues are all real and positive can be computed by, max (1) min The condition nunber of [J] is the square root of k([JTJ]). For a well-conditioned system, the value of k is 1. A very high value of k indicates that the sysIf k exceeds 0a, where a tem is ill-conditioned. equals decimal precision of the digital computer used, it is not possible to obtain a solution [11]. Solutions to eqn. (1) are. given in Table I for the AEP 30-bus , and 43-bus ill-condiwell-conditioned and 11-, tioned test systems. K.M. BROWN'S METHOD .~ -, -. -- Consider the following real continuously differentiable system of N nonlinear equations in N unknowns, x. (i = 1, 2, ..., N). In vector notation f. (x) = 0 ; i-= 1, 2, .., N (2) * * ,-1 K.M. Brown proposed a. local method [7,8] which handles the functions of eqn. (2) one at a time so that information obtained from working with f1 can be. incorporated. when working with f2, etc. A successive substitution scheme is used rather thanl the. simTultaneous treatment of the f. which is the characteristic of Newton's method. It is a quadratically convergent local technique (in the vicinity of a root). It is fast and gi =fi(Xl, X2 N .. LN) (6) fN(Xl1 L2 L3 . 1, * The linear terms L N-l1 are themselves functlons oi x. and are obtained recursively by successive substitution in the system: L = x. 1 1 y n ) (L J Xi) ) (7) -n i = N, N-1, ..., 2 /g N4+1,X 1. k([JTJ]) Ix I with' g = f1 and L1 = xl, so that gN is just a function of the single variable xl. Now expanding, linearizing and solving for x we obtain x n n (8) ) xl x gN,x The xkhl li6nce, is ys vec'tor x* H ve 1 i renam'ed as-.L, a the t sy. and o tem of eqn (7) is back-solved to get improved approxiis mations to the other compoents of x*. Hre taken as the value obtained. for L. when back-s6lving I eqn., (7). The "successive substitu%ion!' nature of the algorithm allows the mst recent information available to. be used in the construction of next function argument, similar to the Gauss-Seidel process for linear and nonlinear systems of equations. mation to the first point xl thus obtained is used as the next approxi- component xl of the solution M M-tri M pe a reetatim . For the sake of -definiteness, it may be stated that the variables are eliminated in the.order xN, xN, T,Using thechain rule for differentiax2. t ~~~~~~tion o expand each derivative g. , gives the following matrix representation for tlieXforward part of the method: -n+ n - x ) = - [g] [H] (x (9) where the muatrix [H] = (h..) is given by 3650 TABLE I Maximum and Minimum Eigenvalues and Condition Number Type of System No. of Buses 30 Maximum Eigenvalue Eigenvalue Minimum Condition Number (k) * Re ]rs emar WVell-conditioned | syteml-conditione system Ill-conditioned system .1087xl13 .1222xlO3 .2905x10 2 .2322x10 0 .1126x100 X13 .4681x10 .1086xl04 .2014x10 . 2560x105 deratelywll~~~~~~~~~~~~~~~conditio'ned (fair Ill-conditioned Moeael (bad k) el k) 11 13 43 .1442xlO-1 14 .2426xl04 .9476xl0_1 Ill-conditioned (bad k) Ill-conditioned (bad k) * Ideal value of condition number: k = 1 flj h.. = i+l (-l) f Ifij lj f 2j 1,N-i f 1N 1,N-i+2 ff 2N 2,N-i+2 f f ~~~~~~j=1,..,N i = 2, J lInitialize | all the variables- ' f | i, N N partial pivoting ~~~~~~~~~~~~~~a Create an array which permits e f'ect -Statieainc-_t> S and f1N-+2 1-l1,N-i+2 fi i-1,N 1 Lerchange rows D | out having to physically in- witi columns pivoting. When the condition number of the matrix [H] << the condition number of the Jacobian [J] convergence occurs. Thus Brown's method gives convergence in cases where Newton's method fails. , . .. (10) where the argunents f.j are progressive arguments generated successively anA [g] is a vector-valued function given by eqn. (10). It is observed that hi; = 0 for j > N-i+l, that is, the matrix [H] is transverse upper triangular. The matrix [H] can be obtained by transforming the Jacobian matrix [J] into transverse upper triangular form using Gaussian elimination with partial Jwhere K - runction numberl K+ ve = KMIN . K-I Solue the first KMIN rows of a triangularized linear system for improved values ofe x in terms of' previous one*s I vet Computational Efficiency of Brown's Method A count of the number of function values of the fi i needed per iteration of Brown's method, is given in Appendix 1. The first step requires N+l evaluations of fl, the second step N evaluations of f2, the third N-l evaluations of f3, etc., so that the total is N+1 es of the kth function and find the largest absolute alue or thp a up approximate partial derviati- absolut apro.artial der vprox-partial deO e l~~argest YS YES NO i=2 i2 = 12~ N+3)(1 1 2 Newton's. method requires N partial derivative evaluations and N function component evaluations per iterative step. Thus there is a corresponding savings in storage locations required: from (N2+N) for Newton's method to +3N) locations for Brown's method. However, it must be stressed that here the savings in function values applies only to saving function values of the fi of the original system (2), because Brown's method adds a number of other fuctions to be evaluated namely the linear functions, L}{. The evaluation of the Lkr o included in the cout above, 1(N2+N L ar no '-N+3) The equnceof omptatonsis how ina smpl flow-chart for Brown's method in Fig. 1 FIXED-POINT FORMUATION Set up coefricients for kth ro or' triangular linear system used to back solve for the f'irst B Try a different i in tial .approximation xck substitute to obtain YES__j -I(N Fig. 1. Simnplified flow-chart for Brow's method The fixed-point formulation of Ref. 13 is used to show that Brow's method also falls into the general category of successive approximation methods like the Gass-Seidel method and Newton-Raphson method, differing only slightly in their iteration fuction. ~~~~~~~solution of a vector set of equations f. (x) = 0 The iS sought by means of a suitable recursive Tormla [14] put into the form: 3651 xn+l n(x ) = xn = p(xn) [Rn] f(xn) (12) (13) where the mapping 1p(xn) is given by By a successive approximation procedure, solution to the fixed-point formulation is given by n = 0 1 2,... (14) xn+l = xn [R'] h(xn, xn+l) ' ' ' ' - ' where [Rn] is a square matrix and h(.,.) is a continuous vector-valued fumction. - same. The stopping criterion is when all the function values are less than the prespecified tolerance. 5. print. Calculate all bus powers and line flows, and The Gauss-Seidel method uses a simple diagonal matrix [R] having scalar acceleration factors for the real and imaginary parts of bus voltages. For Newton's method matrix [R] corresponds to the inverse of the Jacobian matrix, whereas in Brown's method it corresponds to the inverse of Brown's Jacobian matrix [H] given in eqn. (9). LOAD-FLOW BY BROWN'S METHOD Nodes are numbered so that at each step of the Gaussian elimination the next node to be eliminated is the one having the least number of nonzero elements. This is generally preferred for Newton-Raphson method [3,16]. The same strategy is used,for Brown's method also. Henice, the above procedure can be easily incorporated into the existing Newton-Raphson load-flow programs. c) Results and CIarison: Line diagrams of 11-, 13-, and 3-us ill-conditioned test systems are shown in Figs. 2, 3 and 4 respectively. The in system data of these ill-conditioned systems are given Appendix II. Load and generation data are shown only for 13-bus system, whereas, for 11-, and 43-bus systems only injection information is available. The system in Fig. 3 is difficult to solve [6] because of the two series capacitors and the position of the slack-generator. The systems in Figs. 2 and 4 are also difficult to solve [4,12] because of low X/R ratios and some negative line reactances. This is also clear from Table I which gives the eigenvalue analysis and condition number of the Jacobian matrix I [J]show the load-flow equations. The results of in Table that the condition number of the illi oned s t a conditi much highe in systems Jacobian matrix is much higher in conditioned to the well-conditioned systems. comparison Slack-node 1 2 s a) Problem formulation: The load-flow problem involves the solution of a system of nonlinear algebraic equations, that is f. (x) = 0 * i = 1, 2, ..., N (2) For the solution to be unique, onee equation must be oloig specified for each unknown variable, h following equations define more clearly the load-flow problem on the basis of the bus admittance matrix. Corresponding to fi in eqn. (2), the load-flow mismatch equations equal to the number of unknown voltage phase-angles and magnitudes are given by Forcifiedsolution unknwni eqaThe Mismatch: Real Power APi. =P. S I i(Sch) icl i(Cal)=0 (15) (16) ( 7 Reactive Power Mismatch: AQi = n Qi(Sch) Qi(Cal) + B.. The equations for the calculated real and reactive bus powers at bus i are given by P. = 1j= I V.V.(G.. cos(6. - 0.) 1 3 1 13 1 3 sin(e. O.)) 1 J ( 1 Qi = n * * (17) 11 l V.V.(G.. sin(0. - 0.) - Bi; cos(0. - 0.)) (18) Fig. 2. 11 Bus ill-conditioned power system P a 500 The total number of nonlinear and (16) that equations (15) Brown's method is equal to: (n - 1) equations and (n - 1 - number of power mismatch equations. algebraic load-flow are to be solved by real power mismatch PV buses) reactive 0a ? generotor 650 Pao Slack Pa O 0.? Paso b) Solution algorithm: Equations (15) and (16) are to be solved alternately as in [3], for each bus by K.M. Brown's method as follows: 1. Read network line data, generation and load data, and form the Y-bus matrix elements according to the optimal ordering described below. 2. Formulate the system of equations for real and reactive bus power mismatch equations (16 and 17) by using optimal ordering. 3. Initialize the unknown vector xn with alternate voltage phase-angles and voltage miagnitudes, i.e., voltage magnitudes to 1.0 p.u. and angles to 0°. 4. Solve the bus power mismatch equations for the new estimates xl+± by Brown's method and update the Piaso _2 .2 P'O Q a30' Fig. 3. Pa 500 sa? 13 Bus ill-conditioned power system 3652 Slack-.node A comparison of the proposed Brown's load-flow method with the standard methods [1,3,4,5] in terms of convergence characteristics is given in Table II. The proposed Brown's method converges in fewer iterations for both ill- and well-conditioned systems. 12 22 23 28 16D1 (!) (9 37 38 ( 31 21) 30 g 33) ( 36 A comparison for the well-conditioned systems between Brown's load-flow method and the standard Newton's method in terms of ratio of computing time per iteration is given in Table III. Since the NewtonRaphson method does not work for any of the 11-, 13-, and 43-bus ill-conditioned test systems, comparison results are not given for this method. Digital computer @results were obtained by using non-optimized Fortran compilation on ICL 2960 computer with George 2+ operating systems. Even though Brown's method takes 15% more time per iteration than Newton's method, for large systems the total number of iterations is less. This makes Brown's load-flow method comparable to Newton's - \9method in total computing time. -The storage requirement of Brown's method is slightly more (not exceeding 10%) than that of Newton's method. However, this disadvantage is offset by the superiority of Brown's method for ill-conditioned systems. CONCLUSIONS 3 13 18 Fig. 4. 43 Bus ill-conditioned power system The proposed Brown's load-flow method has been applied to the ill-conditioned test systems and AEP 14-, 30-, and 57-bus well-conditioned test systems. The results are shown in Fig. 5 for the convergence rate of the largest absolute mismatch for ill-conditioned as well as for well-conditioned test systems. It can be observed that the convergence rate of the absolute maximum mismatch is quite rapid. The proposed method converges close to the solution on the first iteration and requires few iterations to obtain the final solution. 11 bus system 1- 101 The main objective of this paper is to show that the K.M. Brown technique is suitable for solution of ill-conditioned power systems equations. The technique has quadratic convergence characteristics. The procedure adopted in this method is relatively simple and can be easily incorporated in the existing Newton-Raphson algorithm. The optimal ordering strategy and sparsity programming technique with compact storage scheme can be easily implemented. Also transformers with no-load tap changers and phase shifters can be taken into account while developing production type load-flow programs. It is not necessary to include transformer resistance in the algorithm to help the load-flow to converge to a solution. Brown's load-flow algorithm converges to the solution of ill-conditioned system in a few iterations, whereas the standard loadflow methods either show poor convergence or diverge. -3 bus system bs\1I. bus system 13 bus system I ll - cond; t ioned ACKNOWLEDGEMNTS The authors would like to acknowledge Dr. S. Iwamoto of the Department of Electrical Engineering, Tokai University, Japan for providing data for 11- and 43-bus ill-conditioned systems. REFERENCES t X \ *-- 57 bus system 30 bus syslem Well-conditioned X0 o° [1] A.F. Glimn, and G.W. Stagg, "Automatic calculation of load-flows", AIEE Trans., 1957, PAS-76, pp. 817-825. [2] A. Brameller, and J.K. Demnead, "Some improved methods for digital network analysis", Proc. IEE, 1962, 109, A, pp. 109-116. by Newton's I method"v, l0.2 w A. iJ3L '<>-: \\ X 2 ITERATIONS 3 4 \ \ \[3] W.F. Tinney, and C.E. Hart, "Power flow solution decoupled load-flow" B. Stott, and 0. Alsac, IEEE Trans., 1974, PAS-93, pp. 859-869. [5] A.M. Sasson, "tNonlinear programming solutions for load-flow, minimum loss and economic dispatch problems", IEEE Trans., 1969, PAS-88, pp. 399-409. [4] 1449-1460. ~~~~~~pp. EEE Trans., 1967, PAS-86, 10.4l Fig. 5. 5 "'Fast Gonvergence rate of Brown' s load-flow method for well-conditioned and ill-conditioned systems TABLE II Convergence Characteristics of Different Methods (Number of iterations) System Type Illconditioned systems No. of Buses 11 13 43 14 30 57 GaussSeidel NewtonRaphson Decoupled divergent divergent 22 4 4 Fast Fletcher (Steps) 47 divergent divergent divergent divergent divergent divergent 24 33 59 divergent 3 3 4 divergent 32 68 135 Wlellconditioned systems TABLE III 41 1966). tems of equations", Comparison of Computing Time per Iteration One iteration time in terms ratios (Newton's method as 1 unit) System System Type No. of buses obuss 14 30 57 ~~~~~of [15] R. Fletcher, and M.J.D. Powell, "A rapidly convergent descent method for minimization", Computer Journal, Vol. 6, pp. 163-168, 1963. [16] Newton's method Brown's method Wellconditioned systems 1.0 1.0 1.0 1.17 1.12 1.16 W.F. Tinney and J.W. Walker, "Direct solutions of sparse network equations by optimally ordered triangular factorization", Proc. IEE, Vol. 55, pp. 1801-1809, Nov. 1967. [17] S. Iwamoto, Personal communication with the auAPPENDIX I [6] K. Zollenkopf, "Load-flow calculation using loss minimization techniques", Proc. IEE, 1968, 115(1), pp. 121-127. K.M. Brown's Method [7] K.M. Brown, "A quadratically convergent Newtonlike method based upon Gaussian elimination", SIAM J. Numer. Anal., 1969, pp. 560-569. A step by step synthesis of Brown's method [7] of solving the system of nonlinear algebraic equations is as follows: Let the system be in vector notation as [8] K.M. Brown, "Computer oriented algorithms for solving systems of simultaneous nonlinear algebraic equations", in G.D. BRYNE and C.A. HALL, (eds.), "Nunerical solutions of systems of nonlinear algebraic equations", (Academic Press, 1973), pp. 281-348. [9] L.L. Freris, and A.M. Sasson, "Investigation of the load-flow problem", Proc. IEE, 1968, 115, (10) pp. 1459-1470. f(x) Step 1 0 [10] G.W. Stagg, and A.E. El-Abiad, "Computer methods in power system a;nalysis", (McGraw-Hill, 1968). Let x' denote an approximation to the solution x* of (9). Eipand the first function in the Taylor series expansion about point _,. Alternatively f, can be expanded in an approximate Taylor series with the actual partial derivatives replaced by first difference quotient approximations. In the later case the. user does not have to provide the partial derivative expressions in his program. Retaining only first order terms gives a first order approximation 1 [11] J.H. Wlilkinson, "The algebraic eigenvalue probled'l, (Clarendon P-ress, Oxford, 1965). [12] S. Iwamoto, and Y. Tamura, "A load-flow calculation method for ill-conditioned power systems"', presented at the IEEE PES Summer Meeting, Vanicouver, Canada, 15-20th July, 1979, paper A79 441~~~~~~~~f (x) e f1 (x ) af. + f (xn) +fx (x%(x -xE) l N- f.(x'n + -n ) f.(x)=0*i1,. N(9) Brown's Method 5 4 5 3653 2 3 3 (Academic Press, New York, ;i -1 2, (x - xn) + f 2( n) rx 2 - x) (20) _ i = , ,.., [13] J. Meisel, and R.D. Barnard,. "Application of fixed-point techniques to load-flow studies", IEEE Trans., 1970, PAS-89, pp. 136-140. [14] A.M. Ostrowski, ) ) h~ . .. (21) "Solutions of equations and sys- where ej denotes the j th unit vector and the scalar hn is normally chosen such that hn = 0(1 If(xn) |). With this choice it can be proven 1181 that the discrete 3654 Newton's method has second order convergence. This is the same rate of convergence as the ordinary Newton's method. the Ws are obtained by back-substitution in the N-1 rowed Itriangular linear system which now has the form If xn is close enough to x*, fl(x ) Owe can equate (70) to zero and solve for the ffrst variable. This is the variable, say XN, whose corresponding partial derivative, f (xD, is largest in absolute value. This gives lxN N-1 L= - g - I J=1 /gf Nl J, i )(Lg - xj, /g-NN+l,x.2 5 XN . j=1 (fXNflx j n -X. N )(Xj X f/f N (22) with g- f1 and L1 = xl so that gn is just a function of the single variable xl. Now expanding, linearizing and solving for xl, we obtain The constants f'1 /f'n1 [7] has shown that under the usual hypotheses for Newton's method there is always at least one non-zero partial derivative, and, of course, the corresponding approximate partial (21) is non-vanishing. Thus the solution procedure (22) is well defined. By choosing the approximate partial derivative of largest absolute value as the division, a partial pivoting effect is achieved similar to what is often done when using the Gaussian elimination process for solving linear systems. This enhances the numerical stability of the method. We observe from (22) that xN is a linear function of the N-1 variables xl, x2, ..., XN . For purX poses of clarity the right-hand side of is renamed as L.(x1, x2, *.., xNl) lxi lxN are saved for later use. Brown X, i1 = -'n/a glg, (26) 1 2) Thepoint x thus obtained is used as the next approxiT mation, x'1+ to the first coponent, xl, of the solution vector x*. Renane xl as L1 and back solve the Li system (25) to get an improved approximation of' comonents of x*. Here x4+1 is taken as the value obtained for L. wheii back solving (25). Choice of h in the Computer Implementation In order to guarantee quadratic convergence, the following strategy is used for choosing hM, by which to increment x' when working with the i finction fi in eqn. (21): n (221 Step 2 Define a new function g2 of the N-1 variables xl, .., xN_1 which is related to the second function f2 Of (19) as follows hn7 where = max Cn. 1) 5 x 10 -6+2) 2(X **..., XN-1) f2(X1i, ..., XN-1 Ln(xi, 1' ' { fln n fN-1m 2^***^|gn| ) n ; *001 x xn |} ..(27) Now expand g2 in an approximate Taylor series about the point (x, ..., N_). Linearize (ignore higher order terms) and solve for that variable, say xN_, whose corresponding partial derivative, x is largest in magnitude: N-1 N-2 = ) xN- j1~N-11 22X/2xj - ) 92/2x j=1 ... **(23) expansion where S is the machine tolerance, i.e. the nunber of significant digits carried by the machine. EXAMPLE: Application of K.M. Brown's method Let us consider a nonlinear system with two variables: f (X)=X -2x + 1 = ° X_l is a linear function of the remaining N-2 variables. Remember the right-hand side of (24) is denoted by Li ...`, xN2). Again this forms the ratio, is stored 1, ..., N-2, and 92x 2N-1 for future use in any camputer inplementation of this algorithm. Each step of the algorithm adds one more linear expression to a linear system. During the (k+l)st step of the algorithm, it is necessary to evaluate gk+l, i.e fk1 for various arguments. The values of the last k components of the argument of fk+1 are obtained by back-substitution in the linear system LN, LN-1, .. ., 'N-k+1 which has been built up. These arguments are required to determine the quantities gktl and ,n+l , j = 1, ..., N-k, needed for the elimination of (24) (xl, N2NI/gteration gn/gn 922N-1Ieato x; 0 f (x) = x + 2x - 3 = O 2( 1 2 which has the solution at [1,1] T. For x° = [0,0]T as the initial approximation, the application procedure is demonstrated here using exact partial derivatives. ~~~~~~~~~~~~ Expand f1(x) i the Taylor series about the poit retain only first-order terms and thus obtain the Iinear approximation:+ x + x 0 f1(X) f1(X + f1X (Xi X0 + f1X (X2 2 where teRSo XP (2) ... the (k+l)st variable, say xN..k, by the basic process of expansion, linearization and solution of the resulting expression. The process results, for each k, in the (k+l)st variable, say ._,i expressed as a linear combination, L-k of th~ekremaining N-k+l variables, Step N n ov o h Equ.ate te oeq. eoad ov o h .2)t variable x2, whose corresponding partial derivative is largest in absolute value: xl = 0. 5 1-' lx q.(9 ozr x2 At this stage gN EN f(x1 L2, L3, ... LN) where Substituting this value into function f2 (x), we have 3655 Thus2xl 2.0 = 2 Substitute this value in function define as g2(x): f2(x) in (28) and 1.95 re- Expanding the function g2(x) into the Taylor series (31) 9 2(x) = g2( + g2x (x2 - x2) g 2 g2(x) = 2.5x2 - =j Iteration 211 By repeating the samne procedure with new estimates 1 2.0 f ) 4.0 4.0 -2.0lx2 f1(x(X 4.0 ; f = 4.0 ; f1 where f1 = 2.5 2 Solving for the variable x2 from (31) gives the new estimate after the second iteration g(x ) = -1.0 2 gl Substituting in 22 f x2 - X fx o2 en(30) x( = 1;nd equatingh and equating the RHS of eqn. (30) to ze-ro gives By repeating this same procedure the problem converged in five iterations to xi = 1.0 and x2 = 1.0. x= 0.5x2 + 0.75 f (x) = f(x1) + f f1x f1(x 1 x1 =1.20 (x1-- xl + 2 2x20)9 APPENDIX II System Data for Ill-conditioned Systems A2.2 Description of the 13 bus ill-conditioned system The line and bus data has been taken from [6]. TABLE A2.3 A2.1 Description of the 11 bus ill-conditioned system Line The line and bus data has been taken from [17] and is available only in the form of Y-bus matrix elements and net bus powers. TABLE A2.1 Y-bus matrix elements of 11 bus system data for 13 bus ill-conditioned system To node 1 2 3 2 1 6 2 Branch number 1 2 3 4 6 Gi km k -m 1 2 2 2 3 3 4 4 4 4 1 Bkm k - m 5 6 7 7 8 8 8 11 9 10 10 11 Gkm B km node 0 0 4 3 5 5 5 7 From Resistance (p.u.) 0.0040 0.0040 0.0040 0.0074 0.0481 Reactance (p.u.) 0.0850 Suscep(pu) 0 0 tance 1 -14.939 0.0 14.148 2 0.0 12.051 -33.089 2 6.494 0.0 3 13.197 4 -12.051 2.581 -10.282 3 3.789 5 4 7-2.581 -0.592 0.786 4 12.642 -74.081 2.177 5 0.0 56.689 6 0.0 0.786 -0.592 7 ~ Operating~ 2.581 0.0 3.226 -2.213T 2.893 -0.138 -0.851 0.283 0.104 1.346 -0.374 ~ ~ ~ 0.283~ ~ 5 6 7 8 8 9 10 11 9 10 11 11 -5.889 -55.556 -4.304 ~ ~~ -5.468 1.379 1.163 -2.785 -1.042 -6.110 3.742 -2.785 2.959870 0.,0090 0.0121 0.0105 0 9 8 9 10 9 10 10 MVA 11 11 Base 7-01000 6 12 11 12 12 0075 0.,0086 0 0.2020 -0.1500 0.1665 .16 10 0.0947 0.0947 0.1430 0.4590 0.1080 0.2330 0 0.436 0.246 0.016 0.712 0.620 0 0.508 048 0 ~1codiio systm4-8 0bus6 1107 TABLE A2.2 Operating condition of 11 bus system__________________ TABLE A2.4 Transformer data for 13 bus system Net reactive - 0 . 068~~Bas powe-r Q (p.u.) Bus No. Voltage magnitude6 V (p.u.) 1.024* Phasee angle (deg.) 0.0* Net real. P - 0. 090 power 003 Branch number 100 W 1 2 From node 0 0 4 TO node 1 2 (p.u.) Tap setting + 5%90 +10% +10% 1 3 4 3 6 0.0 -0.128 -0.1028 -0.090 -0.080 -0 .062 -0.068 0.0 *Slack-bus input data 3656 TABLE A2. 5 Bus data for 13 bus system TABLE A2. 7 Operating condition of 43 b-us system Lea Load Bus Bus Assumed Busbus V (p.u.) 1. 0* 1.0 '1.0 1.0 1.0* voltage Assumed phasue- Generation agle 1650 560 0 0 0 0 0 1 0 0 0 2 0 0 0 0 0 3 0 0 0 4 30 50 500 5 1.037* 0 0 0 0 6 1.063 0 0 0 7 1.100* 0 0 500 0.943* 8 0 0 9 0 1.100* 50 30 0 0 10 1.0 50 32 0 0 1.0 11 0 0 0 0 12 1.0 *Input data A2.3 Description of the 43 bus ill-conditioned system The line and bus data has been taken from [17] and is available only in the form of Y-bus matrix elements and net bus powers. TABLE A2.6 Y-bus matrix elements of 43 bus system :_____________________________________________ _ . B kG k -m Gk Bk e (deg.) 0. 0* 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 pG 4) (4VAr)) QG (M) PL (WAr) QL Voltage No. 1 2 3 4 5 6 7 8 9 * magnitude V (p.u.) 1. 136* Phasee (deg.) 0. 0* angle Net real P power - Net (p.u.) power (P.U.) - reactive 1 2 2 3 4 5 5 6 8 8 9 10 10 12 13 14 15 15 16 18 19 21 21 22 23 25 25 27 28 1 2 6 3 4 5 8 12 8 23 10 10 17 12 18 14 15 20 16 18 22 21 29 26 29 25 27 39 29 30 30 30 30 38 29 31 31 32 32 33 38 -30.609 1 0.0 481.288 -1545.194 2 -34.368 108.124 2 -5.714 3 0.0 61.331 -69.160 4 277.195 -916.892 5 20.513 6 0.0 10.638 7 0.0 452.840 -482.861 8 -163.902 167.191 9 12.342 9 -12.045 12.045 -20.855 10 5.714 11 0.0 -10.000 13 0.0 0.0 6.015 13 -15.015 14 0.0 340.398 -916.783 15 15.791 15 0.0 -8.576 17 0.0 -5.714 19 0.0a -164.292 272.805 20 104.312 -143.609 21 -104.312 133.623 22 9.023 23 0.0 -157.677 161.760 24 87.150 -106.814 25 65.824 26 -56.100 0.0 -8.572 28 -202.775 256.136 29 3.766 29 0.0 125.789 -524.464 30 0.0 4.131 30 2 30.609 0.0 5 -277.195 873.583 15 -169.726 534.322 4 6.015 0.062.874 13 -61.331 7 21.277 0.0 34.368 -118.699 6 7 -20.000 0.0 9 -288.938 295.777 9 300.983 -317.044 16 8.796 0.0 2.857 0.0 11 -2.857 0.0 11 92.381 -100.709 13 31.640 25 -31.050 15.400 0.0 43 8.649 19 0.0 28 -170.673 357.003 -5.714 0.0 17 19 164.292 -280.783 .km~ -ki -23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 10 11 12 13 14 15 16 17 18 19 20 21 22 0.0 -0.160 0.0 -0.530 0.0 -1.600 0.0 0.0 0.0 -0.800 0.0 -0.800 0.0 0.0 -0.240 0.0 -0.880 0.0 0.0 0.0 -0.640 0.0 -0.800 -0.320 0.0 0.0 0.0. 1.160 2.900 0.285 0.0 0.580 -0.005 0.0 -1.440 0.0 0.0 -0.800 -2.240 0.0 0.0 0.0 -0.120 0.0 -0.400 0.0 -1.200 0.0 0.0 0.0 -0.600 0.0 -0.600 0.0O -0 .480 0.0 -0.180 0.0 -0.660 0.0 0.0 0.0 -0.480 0.0 -0.600 -0.240 0.0 0.0 0.0 0.520 2.570 0.300 0.0 0.560 0.030 0.0 -1.020 0.0 0.0 -0.300 -1.680 0.0 0.0 -0.640 24 22 23 24 27 26 28 29 37 32 20 34 38 35 38 36 38 0.0 0.0 0.0 0.0 0.0 0.0 -13.038 31 37 -30.769 33 33 3.320 34- 34 6.852 35 35 6.180 36 36 2.703 37 37 9.023 0.0 -8.572 0.0 373.447 -612.837 318.089 -372.311 7.895 0.0 30.769 0.0 40 -125.789 485.547 -15.002 9.267 0.0 164.292 -282.281 321.579 -328.810 -8.572 0.0 0.0 *Slack-bus input data Dr. S.C. Tripathy, Professor of Electrical Engineering at I.I.T. Delhi, India, is presently a visiting professor at The University of Calgary. His research interests are in power system analysis and control. DelTli, with interests in r. G. Dlurga Prasad is a research scholar at I.I.T. power system analysis. 0.0 0.0 0.0 0.0 0.0 0.0 38 39 40 41 43 38 41 40 41 43 -22.398 39 0.0 0.0 15.015 39 125.789 -508.837 40 0.0 -15.015 42 309.806 -408. 029 39 512.581 -663.260 43 -309. 806 392.255 0.0 21.622 42 -20.000 42 0.0 -6.180 -2.703 -21.348 13.038 -3.320 -7.365 'Dr. O.P. Malik, SMIEEE, Professor of Electrical Engineering afid Associate Dean - Academic, The University of Calgary is interested in the real-time digital control of electrical machines and power systems. Dr. G.S. Hope, SIEBE, Professor of Electrical Engineering, Thie University of Calgary, has research interests in the area of digital systems, anld real-time control and protection of power systems. 3657 Discussions J. Nanda, D. P. Kothari and D. L. Shenoy (Indian Institute of Technology, New Delhi,. India): The authors are to be commended for establishing their stand that K. M. Brown's method is reliable for solving ill-conditioned power systems. Authors' definition and consequently the condition they impose for "condition number" for ill-conditioning of computer formulation of a problem are quite interesting indeed. It is stated in the text of the paper that for a well conditioned system, the value of the condition number K is unity. However, in Table I which shows the maximum and minimum eigen values and the condition number K for well-conditioned 30 bus and ill-conditioned 11.13 and 43 bus systems. The condition number for the 30 bus test system is given as 0.468 x 103, which is quite large. Authors have indicated in the remarks column this number K as fair. Would the authors indicate the ranges in which the condition number could be said to be good (well-conditioned), fair (moderately wellconditioned) and bad (ill-conditioned systems)? Table I 11 Bus System Solution as obtained by Newton-Raphson (rectangular co-ordinate version) method. MW i 2 3 4 5 As regards 11 bus system results. The discussers have the following observations to make: 1) 11 bus test system as reported by the authors pertains to an Australian Distribution System originally reported by Mr. B. Stott. 2) For the same 11 bus system the discussers made several tests which revealed that this system has converged by Newton-Raphson (rectangular co-ordinate version) with a flat voltage start (1 +jo) for all load buses in 6 iterations within an accuracy of 0.01 MW/MVAR, while B. Stott's FDLF failed (1]. The results-for the 11 bus system (obtained by Newton-Raphson method) are given in Table-1. Could the authors please clarify why they could not get convergence? REFERENCE [1] P. G. Murthy, D. L. Shenoy, J. Nanda and D. P. Kothari, "Performance of Typical Power Flow Algorithms with reference to Indian Power Systems", presented at Second Symposium on Power Piant Dynamics and Control, Bharat Heavy Electricals Ltd., Hyderabad, India, February, 1979. Manuscript received February 18, 1982. Load in MVAR - Voltage obtained Mag Angle in in p.u. degree 1.024 1.0552 1.0444 1.0306 1.0329 1.0501 0.8310 0.9131 1.1264 0.8266 1.0186 7 8 9 10 11 6 15.8 0.0 12.8 0.0 16.5 9.0 0.0 0.0 2.6 0.0 0.0 6.2 0.0 8.0 6.8 0.0 0.0 0.9 0.0 5.7 0.000(slack bus) 1.215 4.608 2.830 4.830 2.916 12.066 14.857 15.911 21.033 23.996 S. C. Tripathy, G. D. Prasad, O. P. Malik, and G. S. Hope: The authors are thankful for a stimulating discussion on the present paper. The condition number k is a pointer to the extent of separation of system eigenvalues and ideally its value is one if the minimum and maximum values are -equal. It gives a qualitative assessment of the illconditioning effect, so it is improper to define the range of its values quantitatively. Regarding the second question on the solution of the 1 1-bus system, the results of the angle in degrees given in Table 1 by the discussers are in doubt because of their positive signs. Manuscript received April 15, 1982.