A framework to study QV-constraint exchange points in the maximum loadability analysis

Electrical Power and Energy Systems 64 (2015) 347–355

Contents lists available at ScienceDirect

Electrical Power and Energy Systems journal homepage: www.elsevier.com/locate/ijepes

A framework to study QV-constraint exchange points in the maximum loadability analysis R. Salgado ⇑, J. Takashiba Department of Electrical Engineering, Federal University of Santa Catarina, Florianópolis, SC, Brazil

a r t i c l e

i n f o

Article history: Received 17 September 2013 Received in revised form 13 June 2014 Accepted 16 June 2014

Keywords: QV-constraint exchange Limit induced bifurcation Static optimisation Critical loadability Rectangular coordinates

a b s t r a c t This paper presents a framework to analyse critical levels of the power system loadability from the point of view of reactive power generation. Three main aspects are addressed: (1) the preliminary identification of the critical points of reactive power generation in terms of the power system loadability through a predictor scheme, (2) the solution of an optimisation problem to precisely determine these points in the corrector step, and (3) the classification of these critical points in terms of bifurcation type by using by-products of the optimisation process. The conventional power flow equations are expressed in rectangular coordinates, since this type of formulation allows both an accurate estimation of the power network variables in the predictor step and the computation of fast and robust solutions in the correction step. Numerical results obtained with two power networks are shown to illustrate the performance of the proposed methodology. Ó 2014 Elsevier Ltd. All rights reserved.

Introduction Current practices of the power system operation planning require the identification of critical load levels from the point of view of voltage stability. A detailed theoretical framework to model this problem in terms of bifurcation points of the Power Flow (PF) equations is presented in [1], where two types of bifurcation points have deserved special attention: the Saddle-Node Bifurcation (SNB) and the Limit Induced Bifurcation (LIB). The former is characterised by both a singular Jacobian matrix of the algebraic equations representing the power system and the absence of real power flow solution beyond SNB points. The latter is related to changes in the stability condition, as a consequence of a particular variable to reach its limit. It has been shown elsewhere in [2,3] that the set of operational limits which restrains the power flow solutions can be seen as a surface between the feasible and infeasible operating conditions, named solution boundary curve ou transfer limit surface. The limits of reactive power generation, which represent one type of these boundary surfaces, are those with most significative influence on voltage stability [4]. Ref. [5] shows that if the reactive power generation limit is reached at high load levels, the power system is prone to become voltage unstable. This sudden change in the stability resulting from a device to

⇑ Corresponding author. E-mail addresses: [email protected] (R. Salgado), julianatakashiba@ labspot.ufsc.br (J. Takashiba). http://dx.doi.org/10.1016/j.ijepes.2014.06.050 0142-0615/Ó 2014 Elsevier Ltd. All rights reserved.

reaching its limit is referred to as Limit Induced Bifurcation [1]. From the point of view of the power flow problem, the status of the bus whose limit has been hit is modified (conversion PV–PQ), which requires an immediate change in the set of nonlinear equations that represent the power network. These points are considered critical from the point of view of voltage profile and are referred to as QV constraint exchange (or switching) points [6]. Methods developed specifically to identify such points have been proposed in [4,6]. These approaches are based on a predictor–corrector scheme, similar to that of the conventional Continuation Method. A similar approach is shown in [7], with the power flow equations expressed in rectangular coordinates, which provides better results in terms of the predicted solution. Recently, [8] proposed the computation of sensitivity relationships, to classify QV-constraint exchange points. Ref. [9] formulates an optimisation problem in which the reactive power generation limits are represented by inequality constraints on both the field voltage magnitude and the automatic voltage regulator. In [10,11] the limits of capability of the generating units are modelled as complementarity constraints. The first approach proposes a voltage stability constrained optimal power flow, in which the inequality constraints set is enhanced to include the limits on several variables. The second approach provides a detailed analysis of the structure of the boundary loadability surface, based on both the formulation and the results of a constrained optimisation problem. The present paper proposes a framework based on static optimisation to deal with QV-constraint exchange points. The

348

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

concepts of prediction–correction are combined with a strategy to treat the reactive power generation limits. The objective of this combination is to identify these points, to find a corresponding power flow solution and to classify them in terms of the type of bifurcation for voltage stability purposes. The similarities of the proposed approach with respect to the other works found in the literature [4,6] are: (1) the use of a predictor strategy based on sensitivity relationships; (2) the use of a corrector scheme based on both an augmented power flow formulation and the application of Newton’s method; (3) the use of figures of merit based on sensitivity relationships to classify the QV-constraint exchange points. The main contributions of the present work are: (1) the formulation of an optimisation problem with equality constraints expressed in rectangular coordinates, whose properties pointed out in [12,13] ensure an accurate determination of the QV-constraint exchange point and (2) the computation of the figures of merit proposed in [3,8,14] with the aid of the Lagrange multipliers. Numerical results obtained with two test-systems illustrate the main features of the proposed methodology. The remainder of this paper is organised as follows. The basic theory about the LIB points found in power system operation and some strategies to classify these points are stated in Section ‘Preliminaries’. The analytical model proposed herein is presented in Section ‘Proposed approach’. The numerical results obtained through the application of the proposed approach are shown in Section ‘Numerical results’. Section ‘Conclusions’ summarizes the conclusions of the proposed approach.

In order to check the eventual existence of singularities in the point of PV–PQ conversion (represented by the change in l from Eq. 0 to 1),  (1) is expanded in a first order Taylor series, at the point xb ; zrk ; qb , in the direction ðDx; Dzk ; DqÞ, that is,

Jgx Dx þ Jgz Dzk þ Jgq Dq ¼ 0  t   ð1  lÞJzx þ lJhx Dx þ ð1  lÞJzz þ lJhz Dzk h i þ ð1  lÞJzq þ lJhq Dq ¼ 0

ð2Þ

where

@ga ðx; zk ; qÞ @g ðx; zk ; qÞ @g ðx; zk ; qÞ Jgz ¼ a Jgq ¼ a @x @zk @q @zk @zk @zk Jzx ¼ Jzq ¼ Jzz ¼ @x @zk @q @hak ðx; zk ; qÞ @hak ðx; zk ; qÞ @hak ðx; zk ; qÞ Jhz ¼ Jhx ¼ Jhq ¼ @x @zk @q Jgx ¼

are respectively the Jacobian matrices of ga ðx; zk ; qÞ; zk and   hak ðx; zk ; qÞ with respect to x; zk e q, computed at xb ; zrk ; qb . From the first Eq. (2),

h i Dx ¼ J1 g x Jg z Dzk þ Jg q Dq

ð3Þ

and since Jzx ¼ Jzq ¼ 0 e Jzz ¼ 1, the second Eq. (2) is re-written as





lJthx Dx þ ð1  lÞ þ lJhz Dzk þ lJhq Dq ¼ 0

ð4Þ

The combination of Eqs. (3) and (4) yields

h i   lJthx J1 g x Jg z Dzk þ Jg q Dq þ ð1  lÞ þ lJhz Dzk þ lJhq Dq ¼ 0

Preliminaries

ð5Þ

Characterisation of the limit-induced bifurcation points

such that the division of all terms of Eq. (5) by Dq and the rearrangement of the resulting expression provides,

The limit-induced bifurcation point may occur when a power system variable reaches its limit. There are two possibilities concerning the QV-constraint exchange point. In the first, the PV–PQ bus transition is smooth, such that the loadability margin is reduced but the equilibrium point remains voltage stable. These points are named dynamic limit induced bifurcation point [15,16]. In the second scenario, the PV–PQ bus transition provokes an abrupt change in the power system stability, such that the operating equilibrium point immediately becomes voltage unstable. At these points, an additional load increase leads the power system to the voltage collapse, as shown in [5]. These points are alternatively referred to as static limit induced bifurcation points [15], switching loadability limits [11] and saddle limit induced bifurcation points [14]. Refs. [1,15] show the Transversality Conditions to be met at these points. Here, we use the indexes proposed in [3,8] to classify the LIB points. The basic theoretical rationale concerning these indexes is summarised as follows. Consider that the reactive power generation limit of the kth bus has been reached. The set of non-linear equations to be solved in order to find power flow solution at the QV-constraint switching point is expressed as,

h i l Jthx J1 g x Jg q  Jhq Dzk h i ¼ Dq 1 þ l J  Jt J1 J  1 hz hx g x g z

ga ðx; zk ; qÞ ¼ 0   h i lim ð1  lÞ zk  zref þ l hak ðx; zk ; qÞ  hak ¼ 0 k

ð1Þ

where ga ðx; zk ; qÞ ¼ 0 are the power flow equations; l is a parameter whose values (zero and unity) indicate which constraint is currently active (voltage magnitude adjusted in the reference value or reactive power generation set at the limit); zk and zref k represent respectively   the current voltage magnitude (V gk ) and its lim ref reference value V gk ; and hak ðx; zk ; qÞ and hak are the current reactive power generation (Q gk ) and the corresponding limit (Q lim g ).

and since at the critical loadability point

h

1 þ lcrit Jhz 

Jthx J1 g x Jg z

i

ð6Þ Dzk Dq

! 1, then,

1 ¼0

ð7Þ

and therefore,

lcrit ¼

1 1þ

Jthx J1 g x Jg z

ð8Þ

 Jhz

The second component of the denominator of Eq. (8) represents the sensitivity relationship of the reactive power generation with respect to the voltage magnitude of bus k, such that Eq. (8) is rewritten as,

lcrit ¼

1

1 DQ 

ð9Þ

gk

DV g k DQ g

The analysis of Eq. (9) indicates that if DV g k < 0, then k 0 < lcrit < 1, that is, the conversion of the PV bus into PQ bus reaches an unstable QV-constraint switching point. Alternatively, a test based on sensitivity relationships can be applied to check if the power flow solutions at constraint exchange points correspond to limit-induced bifurcation points. According to [3], if the conditions,

Dq > 0 and DV k

Dq >0 DQ g k

ð10Þ

are simultaneously satisfied, then the constraint exchange point is a static (or saddle-node) limit-induced bifurcation point.

349

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

Analytical formulation The saddle-node bifurcation points of the PF equations can be determined by solving a set of nonlinear equations representing the optimality conditions of the optimisation problem [17] expressed as, Maximize q

8   0 2 > > > Pgi  P di þ q DPdi  Pi ðe;f Þ ¼ 0 ðPV and PQ busesÞ > <   subjectto ga ðx; qÞ ¼ 0 ) Q g  Q 0d þ q2 DQ d  Q i ðe;f Þ ¼ 0 ðPQ busesÞ i i i > > > > ref 2 : V i  e2i  fi2 ¼ 0 ðPV busesÞ M ha ðx; qÞ P 0 ) Q m g i 6 Q g i 6 Q g i ðPV busesÞ

ð11Þ

where Pgi and Q gi are the active and reactive power generation of the ith bus; P0di and Q 0di refer to the active and reactive power demand in the base case, DPdi and DQ di represent the pre-specified direction of increase of the power load, V ref is the reference voltage i magnitude of the ith PV bus, and Pi ðe; f Þ and Q i ðe; f Þ are the nodal active and reactive power injections expressed as functions of the real and imaginary components of complex voltages. The optimisation variables are x, representing the real (ei ) and imaginary (fi ) components of the complex bus voltage, and q is the load parameter; the subscripts g and d denote generation and load, respectively, the superscripts ref ; m and M denote respectively the reference value, and the minimum and maximum limits. The quadratic parameterisation of the load change shown in Eq. (11) is equivalent to imposing a non-negativity constraint in the load variation, ensuring that the critical loadability is greater than (or at least equal to) the base case loadability. The constraints of Eq. (11) can be expressed in compact form as

ga ðx; qÞ ¼ y0 þ q2 r  g0 ðxÞ ( m yq0 þ q2 rq0  h0 ðxÞ  h ha ðx; qÞ ¼ M yq0 þ q2 rq0  h0 ðxÞ  h

ð12Þ

1 t 1 x T0 x h0 ðxÞ ¼ xt Tq x 2 2

ð13Þ

where T0 and Tq are tri-dimensional arrays, which depend uniquely on the transmission lines parameters, and the other terms have been previously defined. The optimality conditions for the problem stated by Eq. (11) are expressed by 





Ja ðxa Þt ka

rx £a ðx ; q ; k Þ ¼ 0 ¼   rq £a ðx ; q ; k Þ ¼ 0 ¼  1 þ 2qa rta ka rk £a ðx ; q ; k Þ ¼ 0 ¼ ga ðxa ; qa Þ

ð14Þ

subjectto gb ðx; qÞ ¼ 0 )

8 Pg  ðP 0di þ q2 DPdi Þ  P i ðe;f Þ ¼ 0 ðPV and PQ busesÞ > > > i > > < Q g  ðQ 0d þ q2 DQ d Þ  Q i ðe;f Þ ¼ 0 ðPQ busesÞ i i i 2 > > V ref  e2i  fi2 ¼ 0 ðPV busesÞ > i > > : lim Q gk  ðQ 0dk þ q2 DQ dk Þ  Q k ðe; f Þ ¼ 0 ðk  th PV busÞ

M hb ðx; qÞ P 0 ) Q m g i 6 Q g i 6 Q g i P 0 ðPVbusesÞ; ði – kÞ

ð16Þ

where the number of equality constraints is ð2nb  1Þ and equals the number of optimisation variables. The optimality conditions for this problem are

rx £b ðx ; q ; k Þ ¼ 0 ¼ Jb ðxb Þt kb   rq £b ðx ; q ; k Þ ¼ 0 ¼  1 þ 2qb rtb kb rk £b ðx ; q ; kÞ Þ ¼ 0 ¼ gb ðxb ; qb Þ;

ð17Þ

increased by a set of conditions similar to that stated by Eq. (15). Note that the Jacobian matrix Jb ð:Þ is rectangular, with dimension 2nb  1  2nb  2, and that there is an additional Lagrange multiplier corresponding to the extra equality constraint (reactive power balance of the kth bus).

Eqs. (14) and (17) have similar structures, such that the Taylor series expansion of these equations at the point zte ¼ ½xte qe kte  , in t the direction d ¼ ½Dxt Dq Dkt , to the second-order term is generically expressed by

Jðxe þ DxÞt ðke þ DkÞ ¼ Jðxe Þt ke þ Jðke ÞDx þ Jðxe Þt Dk þ JðDxe Þt Dk    1 þ 2ðqe þ DqÞrt ðke þ DkÞ   ¼  1 þ 2qe rt ke  2Dqrt ke  2qe rt Dk  2Dqrt Dk  gðxe þ Dx; qe þ DqÞ  ¼  gðxe ; qe Þ þ 2qrDq þ Dq2 r  1 1  xt T0 Dx  Dxt T0 Dx ; 2 2

ð15Þ

Fðze þ dÞ ¼ Fðze Þ þ Wðze Þd þ g1 ðdÞ 3 Jðxe Þt ke 6 7 Fðze Þ ¼ 4 ð1 þ 2qe rt ke Þ 5

ð19Þ

2

ð20Þ

gðxe ; qe Þ

kai ¼ 0 for hai ðxa ; qa Þ > 0 where Ja ðxa Þ is the Jacobian matrix of the equality and the active inequality constraints; ga ðx; qÞ includes every hai ðx; qÞ that eventually has become active during the optimisation process; ra is the rate of load variation of equations ga ðx; qÞ, and ka is the vector of Lagrange multipliers of the equality constraints and the active

ð18Þ

where gð:; :Þ and Jð:Þ represent alternatively the equality constraints and the corresponding Jacobian matrix of the problems stated by Eqs. (11) or (16). Eq. (18) is written in compact form as,

where

and additionally

ha ðxa ; qa Þ P 0 kai > 0 for hai ðxa ; qa Þ ¼ 0

Maximise q

Solution method

where the components of vectors y0 þ q2 r and yq0 þ q2 rq0 are related to the bus power injections in the base case (y0 and yq0 ) and the rate of load variation (r and rq0 ). Vectors g0 ðxÞ and h0 ðxÞ represent the power injections and the squared voltage magnitudes expressed in terms of the real and imaginary components of the complex voltages. This allows the power injections to be written as a quadratic form [12] given by,

g0 ðxÞ ¼

inequality constraints. Since ka – 0, this matrix is singular at the saddle-node bifurcation point of the power flow equations. At the maximum loadability point the Lagrange multiplier corresponding to either the reactive power generation limit or the squared voltage magnitude of each PV bus will be zero. In order to find the load level and the corresponding power flow solution at each QV-constraint switching point, the set of equality constraints to be satisfied must include the all equality constraints of Eq. (11) augmented by the reactive power balance equation of the kth bus, that is, the one subject to the QV-constraint exchange. Here, the corrected solution at each QV-constraint switching point is found by solving the following optimisation problem:

2

G 6 Wðze Þ ¼ 4 0 Jðxe Þ

0 2rt ke 2qe r

with G ¼ Jðke ÞDx, and

3 Jðxe Þt 7 2qe rt 5 0

ð21Þ

350

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

2 6 g1 ðdÞ ¼ 4

3

JðDxÞt Dk

Wðze Þd2 ¼ g1 ðd1 Þ

7 5

2Dqrt Dk

ð22Þ

Dxt T0 Dx  Dq2 r

Features of the optimal solutions The relationship between the solutions of the optimisation problems stated by Eqs. (11) and (16) depends on the nature of the QV-constraint switching point. Suppose that the power flow solutions (xa and xb ) as well as the loadability factors (qa and qb ) that satisfy the optimality conditions stated by Eqs. (14) and (17)   are equal. Then, there is a sub-matrix Jb xb , of order ð2nb  2Þ  ð2nb  2Þ, such that the Jacobian matrix Jb xb can be partitioned as,

" ¼

#

Jb ðxb Þ

Ja ðxa Þt

Jbk ðxb Þ

" # i kb  t kbk

 ¼ 0  1 þ 2qb rb t

"

r bk

 kb  kbk

#!

¼0 ð24Þ

where kb  and rb ¼ ra are column vectors of order ð2nb  2Þ; and kbk and r bk represent respectively the Lagrange multiplier of the reactive power balance equation in the optimisation model of Eq. (16) and the element of vector rb corresponding to bus k. Eq. (24) is satisfied for kb  ¼ ka and kbj ¼ 0, which confirms that the QV-constraint exchange point is also a saddle (or static) limit induced bifurcation point. In summary, once the solution     xb ; qb ; kb is available, if Eq. (24) can be solved for kb  – 0 and kbk ¼ 0, then the QV-constraint exchange point represents simultaneously both a saddle-node bifurcation point and a static limit induced bifurcation point. Proposed approach The main stages of the proposed procedure to calculate the QVconstraint exchange points are described as follows. Prediction stage For prediction purposes, the vector d of Eq. (20) can be approximately evaluated in two steps, by solving two linear systems with the same coefficient matrix. The first represents a linear approximation of Fðze Þ, that is,

Wðze Þd1 ¼ Fðze Þ

ð27Þ

where d ¼ d1 þ ad2 and a is a step factor calculated to minimise the squared norm of Fðze þ dÞ, as proposed in [18]. Step factor computation The solution obtained in the correction stage corresponds to the QV-constraint exchange point of the kth PV bus. The inequality constraints on the reactive power generation of the other PV buses must be also satisfied. In order to find out which bus reaches the reactive power generation limit, during the prediction step we use the complete Taylor series expansion of hðx; qÞ (see Eqs. (12) and (13)) at the point ðxe ; qe Þ, in the direction ðaDx; aDqÞ, that is,

hðxe þ aDx; qe þ aDqÞ ¼ hðxe ; qe Þ   þ a Jq ðxe ÞDx þ 2qe Dqrq0 þ a2 ðh0 ðDxÞ þ Dq2 rq0 Þ

ð23Þ

Jbk ðxb Þ

        where Jb xb ¼ Ja xa , and Jbk xb is the extra row of Jb xb corresponding to the reactive power balance of the kth bus. From the point of view of saddle-node bifurcation, the set of equality constraints of Eq. (11) includes alternatively the squared voltage equation or the reactive power balance equation of the kth bus. Thus, in order to satisfy the conditions stated by Eq. (17) one of the Lagrange multiplies corresponding to these equations will be zero. Here, it is assumed that there is only one singularity at the maximum loadability point and then rankðJb ðx ÞÞ ¼ 2nb  3. Now, consider the first two Eq. (17) and the partition stated by Eq. (23), that is,

h

such that the predicted solution of Eq. (20) is given by

zp ¼ ze þ d

The solution of the problems stated by Eqs. (11) and (16) through Newton’s method is based on the first order approximation of Eq. (19). A strategy to use the second order term expressed by Eq. (22) is presented in [18].

Jb ðxb Þ

ð26Þ

ð25Þ

The second linear system takes advantage of the second-order term of Eq. (19). Thus, estimates d1 can be improved by adding the increments obtained from the solution of the linear system,

ð28Þ

where Jq ðxe Þ is the Jacobian matrix of equations h0 ðxÞ and a is a scalar to be used as a step factor. The combination of Eqs. (12) and (28) provides

a þ ab þ a2 c 6 0

ð29Þ

where

a ¼ hðxe ; qe Þ  h

lim

b ¼ Jq ðxe ÞDx þ 2qe Dqrq0

ð30Þ

c ¼ h0 ðDxÞ þ Dq2 rq0 For each generated reactive power device there is a quadratic equation of the form

ai þ abi þ a2 ci 6 0 which is a component of Eq. (29). The sequence of quadratic equations represented by Eq. (29) is solved and the minimum positive value of a, such that 0 < a 6 1, is used to adjust the increments Dei ; Dfi ; Dki ; Dqi . This procedure allows to determine which reactive power source is the first to hit its limit and prevents the violation of any other limit of reactive power generation in the prediction step. Note that, due to the use of the second-order term of the Taylor series expansion, this strategy provides the best estimate of the power system variables (in terms of accuracy) resulting from a load increase. Correction stage The exact solution of the optimisation problem stated by Eq. (16) is obtained in the correction stage. For this purpose, the reactive power balance equation of the bus selected in the stage of step size calculation is included in the set of equality constraint corresponding to the optimality conditions stated by Eq. (17). The predicted solution is used as the initial estimate in the correction stage, as this facilitates the convergence of the iterative process. Classification of the QV-constraint switching point The solution of Eq. (17) provides the sensitivities of the squared load parameter q2 with respect to the active power injections, reactive power injections and quadratic voltage magnitudes of the PV buses as a by-product of the optimisation process. These sensitivities are numerically represented by Lagrange multipliers kpi ; kqi and kv i , that is,

351

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

kpi ¼

Dq DP i

kqi ¼

Dq DQ i

kv i ¼

Dq DV 2i

ð31Þ

and, since

Dq ¼ DV i

Dq DV 2i

!

DV 2i DV i

! ¼ 2V i kv i

such that the sensitivity relationship between the reactive power generation and the voltage magnitude of the bus subject to QV-constraint exchange is given by,

DQ i 2V i kv i ¼ DV i kQ i

ð32Þ

Then, after solving the optimisation problem stated by Eq. (16), the computation of the sensitivities as given by Eq. (32) as well as the classification of the QV-constraint switching point based on Eqs. (9) and (10) is straightforward. Algorithm The procedure to apply the proposed methodology can be summarised in the following sequence of steps: 1. determine the power flow solution (for the base case or for the maximum loadability level); 2. estimate the predicted solution by solving the linear systems of Eqs. (25) and (26); 3. determine the smallest step size for which a PV bus reaches its reactive power generation limit, through Eqs. (29) and (30). If there is no such bus, the process is finished; otherwise, proceed to the next step; 4. solve Eq. (17) to obtain the corrected solution; 5. classify the LIB point based on Eqs. (9) or (10) and return to step 2. Due to the use of the accurate predicted solution (computed in step (2) of the algorithm) and the way the optimisation problem is (power flow) constrained, only a few (usually 2–5) iterations are necessary to solve the optimisation problem of step (2). Numerical results With the aim of illustrating the performance of the proposed methodology, main results of the computational simulation performed with two test systems are presented in this section. 5-Bus test system The one-line diagram of the 5-bus test system is shown in Fig. 1. The series resistance and the shunt admittance of the transmission lines are neglected. The series reactance of each transmission line is equal to 0.1 pu. The rate of load increase at bus 4 is 0.5 pu(MW) 3 1

and 0.15 pu(Mvar) whereas the load of bus 5 is constant. The reactive power generation limits of buses 1, 2 and 3 are respectively 110, 92.5 and 90 Mvar. The power flow solutions for the base case (q ¼ 0) and for the maximum loadability level (q ¼ 3:808) are shown in Tables 1 and 2. Between these two solutions there are two QV-constraint exchange points. The reactive power generation of buses 3 and 2 (in this order) have reached the limit when the load parameter is 2.068 and 3.808, respectively. Table 3 summarizes the main characteristics of these QV-constraint exchange points. Columns 2–4 show respectively the reference voltage magnitude for the base case, the voltage magnitude at the maximum loadability point and the reactive power generation. Column 5 shows the load parameter corresponding to each QV-constraint exchange point. Columns 6–8 show the Lagrange multipliers corresponding to the active power injections, the reactive power injections and the squared voltage magnitudes, respectively. The figures of merit used to classify the limit induced bifurcation points, in accordance with Eqs. (9) and (10) are shown in columns 7–10. In order to show additional details of these QV-constraint exchange points, Figs. 2 and 3 show respectively the variations of both the voltage magnitude and the reactive power generation in terms of the load parameter. In the base case, Q g2 ¼ 0:87 Mvar and Q g3 ¼ 20:58 Mvar. As the power load increases, the upper limit QM g 3 ¼ 90 Mvar is reached for q ¼ 2; 0688, with the voltage magnitude of bus 3 set at its reference value. According to the figures of merit presented in the columns 7, 8, 9 and 10 of Table 3, this QVconstraint exchange corresponds to a dynamic limit induced bifurcation. With the additional load increase, the voltage magnitude of bus 3 decreases from 1.20 pu (base case) to 1.124 pu, until the maximum loadability and the reactive power generation limit of bus 2 (Q M are simultaneously reached, for g 2 ¼ 92:5 Mvar) q ¼ 3:808, with V 2 ¼ 1:200 pu. In this case, figures of merit presented in the columns 8, 9 and 10 of Table 3 indicate that this point is related to a static limit induced bifurcation. The analysis of Fig. 2 shows that the solutions of upper and lower voltage magnitude level coalesce in a critical point of loadability, beyond which no real solutions for power flow equations exist, confirming the nature of the LIB point. Subsequently, the reactive power generation limit of bus 2 was increased from 92.5 Mvar to 110 Mvar. The power flow solution at the maximum loadability level (q ¼ 3:890) and the main characteristics of the QV-constraint exchange points are shown in Tables 4 and 5, respectively. In this case, the sign of the figures of merit shown in columns 7–10 of Table 5 indicate the dynamic nature of both QV-constraint switchings. This is confirmed by the variations of the voltage magnitude and reactive power generation with respect to the load parameter shown in Figs. 4 and 5, where it can be seen that there is no PV bus with simultaneous reactive power generation at the limit and corresponding voltage magnitude at the reference value. In addition, note that the increase of 17.5 Mvar in the reactive power limit of bus 2 results in an increase in the voltage magnitude V 2 from 1.200 pu to 1.226 pu, confirming the location of the QV-constraint switching point of the previous case in the lower portion of the PV-curve.

2 Table 1 5-Bus test system: power flow solution – base case.

5

4

Fig. 1. 5-Bus test system: one-line diagram.

Bus

V (pu)

d ð Þ

1 2 3 4 5

1.200 1.200 1.200 1.173 1.185

0.00 2.98 0.99 2.65 2.03

Total

P G (MW) 5.00 50.00 50.00 0.00 0.00 105.0

Q G (MVAr) 0.22 0.87 20.58 0.00 0.00 21.67

P Dbase (MW) 30.00 0.00 0.00 15.00 60.00 105.0

Q Dbase (Mvar) 0.00 0.00 0.00 15.00 0.00 15.00

352

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

Table 2 5-Bus test system: power flow solution – static limit induced bifurcation. Bus

V (pu)

d ð Þ

P G (MW)

Q G (MVAr)

P Dbase (MW)

Q Dbase (Mvar)

1 2 3 4 5

1.200 1.200 1.124 0.806 0.931

0.00 4.92 7.04 37.60 21.73

196.51 50.00 50.00 0.00 0.00

101.75 92.50 90.00 0.00 0.00

31.15 0.00 0.00 205.37 60.00

0.00 0.00 0.00 72.11 0.00

296.51

284.25

296.51

72.11

Total

Table 3 5-Bus test system: QV-constraint exchange points – static limit induced bifurcation. Bus

V ref (pu)

V crt (pu)

QM g (MVAr)

q (pu)

kpi (pu)

kqi (pu)

kv i (pu)

@Q i @V i

3 2

1.200 1.200

1.124 1.200

90.00 92.50

2.068 3.808

0.113 0.028

2.081 0.409

20.053 0.202

23.125 1.181

(pu)

1.1

Bus

V (pu)

d ð Þ

PG (MW)

QG (MVAr)

P Dbase (MW)

Q Dbase (Mvar)

1 2 3 4 5

1.200 1.226 1.137 0.822 0.946

0.00 5.08 7.14 37.26 21.64

200.09 50.00 50.00 0.00 0.00

86.07 110.00 90.00 0.00 0.00

30.58 0.00 0.00 209.50 60.00

0.00 0.00 0.00 73.35 0.00

300.09

286.07

300.09

73.35

1

voltage magnitude (pu)

0.045 0.458

Table 4 5-Bus test system: power flow solution – saddle node bifurcation.

1.2

0.9 0.8 0.7

Total

0.6 0.5 0.4

bus 1 bus 2 bus 3 bus 4 bus 5

0.3 0.2

1

1.5

2

2.5

3

3.5

4

loading factor (pu) Fig. 2. 5-Bus test system: V  q – static limit induced bifurcation.

1

reactive power generation (pu)

lcrit (pu)

0.8

0.6

0.4

0.2

0 bus 1 bus 2 bus 3

−0.2 1

1.5

2

2.5

3

3.5

4

loading factor (pu) Fig. 3. 5-Bus test system: Q g  q – static limit induced bifurcation.

IEEE 118-buses test system This system has a total demand of 3668.00 MW and 1438.00 MVAr. A rate of load change equal to 1.0 % of base load

(with constant power factor) was assumed for each bus. At the point of maximum loadability, 32 buses (from the total of 54 voltage control buses) have reached the maximum limit of generation of reactive power. Except for the last one, these QV-constraint switching points are not typically associated with the loss of voltage stability [14]. For this reason, the main features of the QV-constraint exchange points relative only to the five last buses that have reached the upper limit of reactive power generation are shown in Table 6. The critical loadability power flow solution is obtained for q ¼ 60:763, corresponding to a power demand of 5896.80 MW and 2311.78 Mvar. Column 5 shows the value of the load parameter for which the limit of reactive power generation was reached and columns 7–10 show the figures of merit used to classify the QV-constraint exchange points. Figs. 6 and 7 show respectively the variations of the voltage magnitude and reactive power generation of the last 5 buses that have reached the limit of generation of reactive power as functions of the load parameter. The voltage magnitude of these buses is constant until the limit of reactive power generation is reached. At this point, there is a QV-constraint exchange. As the load increases, the bus voltage magnitude, which is no longer a voltage control variable, decreases, until the point that the maximum loadability is reached. In the present case, the last QV-constraint switching point, corresponding to the PV bus 61, is reached for q ¼ 60:763. No further increase of the power demand is allowed, since the point of maximum loadability has been also reached for q ¼ 60:763. The sensitivity relationships shown in columns 7–10 and the parameter lcrit of Table 6 indicate that this point corresponds to a static limit induced bifurcation. In order to check the consistency of classification inferred with the results of the table shown in the previous case, the reactive power generation limit of bus 61 was increased from 200 Mvar to 290 Mvar. Under these conditions, the maximum loadability is reached for q ¼ 61:796, that is, a total demand of 5934.66 MW and 2326.62 MVAr can be supplied, satisfying the reactive power generation limit. The same 32 buses of previous case reached the

353

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355 Table 5 5-Bus test system: QV-constraint exchange points – dynamic limit induced bifurcation. Bus

V ref (pu)

V crt (pu)

QM g (MVAr)

q (pu)

kpi (pu)

kqi (pu)

kv i (pu)

@Q i @V i

lcrit (pu)

3 2

1.200 1.200

1.137 1.226

90.00 110.00

2.068 3.839

0.113 0.028

2.081 0.054

20.053 1.654

23.125 74.103

0.045 0.013

1.06

1.2

1.04

1.1

1.02

voltage magnitude (pu)

voltage magnitude (pu)

1 0.9 0.8 0.7 0.6 0.5 0.4

1.5

2

2.5

3

3.5

0.96 0.94 0.92

bus 66 bus 61 bus 46 bus 59 bus 85

0.88

0.2 1

0.98

0.9

bus 1 bus 2 bus 3 bus 4 bus 5

0.3

1

0.86

4

1

1.1

1.2

loading factor (pu) Fig. 4. 5-Bus test system: V  q – saddle node bifurcation.

1.4

1.5

1.6

Fig. 6. IEEE118-Bus test system: dynamic limit induced bifurcation – variation q  V.

2

1

reactive power generation (pu)

1.3

loading factor (pu)

bus 66 bus 61 bus 46 bus 59 bus 85

1.5

reactive power generation (pu)

0.8

0.6

0.4

0.2

0 bus 1 bus 2 bus 3

1

0.5

0

−0.5

−0.2 1

1.5

2

2.5

3

3.5

4

−1

loading factor (pu)

1

1.1

1.2

1.3

1.4

1.5

1.6

loading factor (pu)

Fig. 5. 5-Bus test system: Q g  q – saddle node bifurcation.

upper limit of reactive power generation. The characteristic of the QV-constraint exchange points of for that last 5 buses which reached the reactive power generation limit are shown in Table 7. Except for the figures of merit of bus 61, the results for the other 31 buses do not change. These figures indicate that the last

Fig. 7. IEEE118-bus test system: dynamic limit induced bifurcation – variation q  Qg.

QV-constraint exchange is located at a dynamic limit induced bifurcation point. This is confirmed by the behaviour of the voltage magnitude of bus 61. The increase in the limit of reactive power

Table 6 IEEE118-bus test system: QV-constraint exchange points – static limit induced bifurcation. Bus

V ref (pu)

V crt (pu)

QM g (MVAr)

q (pu)

kpi (pu)

kqi (pu)

kv i (pu)

@Q i @V i

66 72 99 54 61

1.050 0.980 1.010 0.955 0.988

0.993 0.946 0.997 0.941 0.988

200.00 100.00 90.00 300.00 200.00

51.741 58.828 59.535 60.140 60.763

1.852 7.094 6.063 4.095 1.478

7.284 9.242 10.261 3.629 1.210

134.049 42.781 60.650 28.173 0.217

36.366 9.147 11.679 15.340 0.354

(pu)

lcrit (pu) 0.026 0.123 0.091 0.072 0.738

354

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355

Table 7 IEEE 118-bus test system: QV-constraint exchange points – dynamic limit induced bifurcation. Bus

V ref (pu)

V crt (pu)

QM g (MVAr)

q (pu)

kpi (pu)

kqi (pu)

kv i (pu)

@Q i @V i

66 72 99 54 61

1.050 0.980 1.010 0.955 0.988

1.008 0.939 0.990 0.963 1.029

200.00 100.00 90.00 300.00 290.00

51.741 58.828 59.535 60.140 61.161

1.852 7.094 6.063 4.095 0.945

7.283 9.246 10.261 3.629 0.342

134.049 42.781 60.650 28.173 15.283

36.366 9.147 11.679 15.340 88.131

(pu)

lcrit (pu) 0.026 0.123 0.091 0.072 0.011

@q same of the previous case, the sensitivity relationship @Q becomes g negative, as a result of the increase in the reactive power generation limit of bus 61. These results characterise QV-constraint exchange as a dynamic limit induced bifurcation point.

1.06 1.04

voltage magnitude (pu)

1.02

Conclusions

1

0.98 0.96 0.94 0.92 0.9

bus 66 bus 61 bus 46 bus 59 bus 85

0.88 0.86 1

1.1

1.2

1.3

1.4

1.5

1.6

loading factor (pu) Fig. 8. EEE118-bus test system: static limit induced bifurcation – variation q  V.

This study emphasises the advantages of both formulating the power flow equations in rectangular coordinates and modelling the determination of QV-constraint exchange points as an optimisation problem. The quadratic characteristic inherent to former increases the accuracy of the predicted solution, since the availability of the increments in the power flow variables allows to determine precisely the variation of the generated reactive power. Then a scheme of prediction combined to the computation of a step factor has been developed to estimate each QV-constraint exchange point. On the other hand, further than an accurate power flow solution and the corresponding loadability factor obtained by solving an optimisation problem, the Lagrange multipliers, computed as a by-product of the iterative process, can be used to classify the power flow solutions from the point of view of LIB points. The numerical results presented here confirm the features of the proposed methodology. Future work in this subject include the determination of bifurcation points induced by power flow limits in the transmission lines.

3 bus 66 bus 61 bus 46 bus 59 bus 85

reactive power generation (pu)

2.5

Acknowledgment The authors acknowledge financial support from the Brazilian government funding agency Conselho Nacional de Desenvolvimento Científico e Tecnológico, CNPq/Brazil.

2

1.5

References

1

0.5

0

−0.5

−1 1

1.1

1.2

1.3

1.4

1.5

1.6

loading factor (pu) Fig. 9. EEE118-bus test system: static limit induced bifurcation – variation q  Q g .

generation results in a corresponding variation in the voltage magnitude from 0.988 to 1.029 pu, which confirms that the QV-exchange point of the previous case is located at the lower portion of the PV curve corresponding to the bus 61. Figs. 8 and 9 show the variation of the voltage magnitude and reactive power generation of the buses referred previously, in terms of the loadability level. With respect to bus 61, note that despite the QV-constraint switching point remains virtually the

[1] Venkatasubramanian V, Schättler H, Zaborsky J. Dynamics of large constrained nonlinear systems – a taxonomy theory. Proc IEEE 1995;83(11):1530–61. [2] Hiskens IA, Davy R. Exploring the power flow solution space boundary. IEEE Trans Power Syst 2001;16(3):559–65. [3] Kataoka Y, Shinoda Y. Voltage stability limit of electric power systems with generation reactive power constraints considered. IEEE Trans Power Syst 2005;20(2):951–62. [4] Hiskens IA, Chakrabarti BB. Direct calculation of reactive power limit points. Int J Electr Power Energy Syst 1996;18(2):121–9. [5] Dobson I. Voltage collapse precipitated by the immediate change in stability when generator reactive power limits are encountered. IEEE Trans Circ Syst – I Fundam Theory Appl 1992;39(9):559–65. [6] Yorino N, Li HQ, Sasaki H. A predictor/corrector scheme for obtaining q-limits points for power flow studies. IEEE Trans Power Syst 2005;20(1):130–7. [7] Salgado RS. Determination of critical reactive power points in power system loadability studies. In: Proceedings of the IFAC power plant and power system control; 2012. p. 1–6. [8] Echavarren FM, Lobato E, Rouco L. Steady state analysis of the effect of reactive generation limits in voltage stability. Electr Power Syst Res 2009;09:1292–9. [9] Vournas CD, Karystianos M, Maratos NG. Bifurcation points and loadability limits as solutions of constrained optimisation problems. In: IEEE 1995 winter power meeting panel session on challenges to OPF, vol. 3; 2000. p. 1883–8. [10] Rosehart W, Roman C, Schellenberg A. Optimal power flow with complementarity constraints. IEEE Trans Power Syst 2005;20(2):813–22. [11] Karystianos M, Maratos NG, Vournas CD. Maximising power system loadability in the presence of multiple binding complementarity constraints. IEEE Trans Power Syst 2007;54(8):1775–87.

R. Salgado, J. Takashiba / Electrical Power and Energy Systems 64 (2015) 347–355 [12] Iwamoto S, Tamura Y. A fast load method retaining nonlinearity. IEEE Trans Power Ap Syst 1978;97(5):1586–99. [13] Makarov Y, Hill D, Hiskens I. Properties of quadratic equations and their application to power system analysis. Int J Electr Power Energy Syst 2000; 22(3):313–23. [14] Cañizares CA, Mithulananthan M, Berizzi A, Reeve J. On the linear profile of indices for the prediction of saddle-node and limit-induced bifurcation points in power systems. IEEE Trans Circ Syst – Regular Papers 2003;50(12):1588–95. [15] Avalos RJ, Cañizares CA, Milano F, Conejo AJ. Equivalency of continuation and optimization methods to determine saddle node and limit-induced

355

bifurcations in power systems. IEEE Trans Circ Syst – Regular Papers 2009; 56(1):210–32. [16] Venkatasubramanian V, Schättler H, Zaborsky J. Dynamics of large constrained nonlinear systems – a taxonomy theory. Proc IEEE 1995;83(11):1530–61. [17] Cañizares C. Applications of optimization to voltage collapse analysis. In: IEEE 1995 winter power meeting panel session on challenges to OPF; 1998. p. 1–8. [18] Salgado RS, Zeitune AF. A direct method based on tensor calculation to determine maximum loadability power flow solutions. Electr Power Syst Res 2013;103:114–21.