Sensitivity based Thevenin index for voltage stability assessment considering N-1 contingency

Abstract

Systems and method are disclosed for monitor long-term voltage stability in a power system by estimating post-contingency operating point(s) considering PV-PQ transitions based on the current operating point; predicting a sensitivity based Thevenin index (STI) for each contingency using the estimated operating condition; and determining a voltage stability assessment (VSA) considering N−1 contingency.

Claims

1. A method for determining long-term voltage stability in a power system, comprising: estimating post-contingency operating point(s) of the power system considering PV-PQ transitions based on the current operating point; predicting a sensitivity based Thevenin index (STI) for each contingency using the estimated operating condition; and determining a voltage stability assessment (VSA) considering N−1 contingency, wherein the current operating states {tilde over (θ)},{tilde over (V)}, further comprising determining post-contingency operating states as: $\left[\begin{array}{c}{\theta }_c\\ V_c\end{array}\right]={\left[\begin{array}{c}\tilde{\theta }\\ \tilde{V}\end{array}\right]\left[\begin{array}{cc}{f}_{\theta }& f_V\\ g_{\theta }& g_V\end{array}\right]}^{-1}\Delta F_{\mathrm{inj}}$ where ({tilde over (θ)}, {tilde over (V)}) is an operating point under base case; (θ.sub.c,V.sub.c) is an operating point following a contingency; Δ F.sub.inj is a vector containing the active and reactive bus injection changes which are numerically equal to the active and reactive line flow on an outage branch; f.sub.θ and f.sub.V represents the partial derivative of active power injection expression with respect to voltage angle θ and voltage magnitude V; g.sub.θ and g.sub.V denotes the partial derivative of reactive power injection expression with respect to voltage angle θ and voltage magnitude V; and operating a power generator adapted to be connected to a power grid and controlled by a processor.

2. The method of claim 1, further comprising: obtaining current operating states; estimating a post-contingency operating point θ.sub.c,V.sub.c for every considered contingency, where θ.sub.c and V.sub.c represent the predicted bus voltage angle and bus voltage magnitude after a contingency; determining corresponding sensitivities based on the presence of PV-PQ transition following the contingency; and determining the STI for each considered contingency.

3. The method of claim 1, comprising estimating an operating point after a foreseen transmission contingency.

4. The method of claim 1, comprising construction of a Jacobian matrix based on a system topology with a contingency line ij in an outage.

5. The method of claim 1, comprising determining corresponding sensitivities as one predictor step of continuation of power flow (CPF) when dλ is a continuation parameter with: $\left[\begin{array}{ccc}{f}_{\theta }& f_V& f_{\lambda }^P\\ g_{\theta }& g_V& g_{\lambda }^Q\\ 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]$ where dλ represents the load increment; dθ and dV denote the change of voltage angle and voltage magnitude respectively; f.sup.p.sub.λ,g.sup.p.sub.λ indicate the amount of active and reactive power injection change as a function of dλ.

6. The method of claim 5, comprising obtaining sensitivities when there are generators reaching their reactive power limit after the contingency as: $\left[\begin{array}{cccc}{f}_{\theta }& f_V& f_{V_i}& f_{\lambda }^P\\ g_{\theta }& g_V& g_{V_i}& g_{\lambda }^Q\\ g_{\theta }^i& g_V^i& g_{V_i}^i& g_{\lambda }^{Q,i}\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ {\mathrm{dV}}_i\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 1\end{array}\right]$ where bus i is the bus which transforms from PV bus to PQ bus; f.sub.V.sub.i is the partial derivative of active power injection expression at other PV and PQ buses with respect to bus voltage magnitude at bus i, i.e., V.sub.i; g.sub.V.sub.i is the partial derivative of reactive power injection expression at other PQ buses with respect to bus V.sub.i; g.sub.θ.sup.i, g.sub.V.sup.i, g.sub.V.sub.i.sup.i, represent the partial derivative of reactive power injection expressions at bus i with respect to θ,V and V.sub.i; g.sub.λ.sup.Q,i is the partial derivative of reactive power injection expression at bus i with respect to the load scaling factor A.

7. The method of claim 1, comprising determining STI for each contingency and load bus in parallel.

8. The method of claim 1, comprising determining a plurality of contingencies in parallel.

9. A system to control a power system, comprising: a power generator; a power grid coupled to the power generator; and a processor to control the power grid with code for: estimating post-contingency operating point(s) of the power system considering PV-PQ transitions based on the current operating point; predicting a sensitivity based Thevenin index (STI) for each contingency using the estimated operating condition; determining a voltage stability assessment (VSA) considering N−1 contingency; and determining corresponding sensitivities as one predictor when dλ is a continuation parameter with: $\left[\begin{array}{ccc}{f}_{\theta }& f_V& f_{\lambda }^P\\ g_{\theta }& g_V& g_{\lambda }^Q\\ 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]$ where dλ represents the load increment, dθ and dV denote the change of voltage angle and voltage magnitude respectively, f.sup.p.sub.λ,g.sup.p.sub.λ indicate the amount of active and reactive power injection changes as a function of dλ; and where f.sub.θ and f.sub.V represents the partial derivative of active power injection expression with respect to voltage angle θ and voltage magnitude V−g.sub.θ and g.sub.V denotes the partial derivative of reactive power injection expression with respect to voltage angle θ and voltage magnitude V.

10. The system of claim 4, further comprising code for: obtaining current operating states; estimating a post-contingency operating point θ.sub.c,V.sub.c considering the presence of PV-PQ transition for every considered contingency; determining corresponding sensitivities following the contingency; and determining the STI for each considered contingency.

11. The system of claim 9, comprising code for obtaining sensitivities when there are PV-PQ transitions as: $\left[\begin{array}{cccc}{f}_{\theta }& f_V& f_{V_i}& f_{\lambda }^P\\ g_{\theta }& g_V& g_{V_i}& g_{\lambda }^Q\\ g_{\theta }^i& g_V^i& g_{V_i}^i& g_{\lambda }^{Q,i}\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ {\mathrm{dV}}_i\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 1\end{array}\right]$ where bus i is the bus which transforms from PV bus to PQ bus; f.sub.V.sub.i is the partial derivative of active power injection expression at other PV and PQ buses with respect to bus voltage magnitude at bus i, i.e., V.sub.i; g.sub.V.sub.i is the partial derivative of reactive power injection expression at other PQ buses with respect to bus V.sub.i; g.sub.θ.sup.i,g.sub.V.sup.i,g.sub.V.sub.i.sup.i represent the partial derivative of reactive power injection expressions at bus i with respect to θ,V and V.sub.i; g.sub.λ.sup.Q,i is the partial derivative of reactive power injection expression at bus i with respect to the load scaling factor λ.

12. The system of claim 9, comprising code for determining STI for each contingency and load bus in parallel.

13. The system of claim 9, comprising code for determining a plurality of contingencies in parallel.

Description

BRIEF DESCRIPTIONS OF FIGURES

(1) The subject matter that is regarded as the invention may be best understood by reference to the following description taken in conjunction with the accompanying drawing figures in which:

(2) FIG. 1 shows an exemplary two-bus Thevenin equivalent system.

(3) FIG. 2 shows an exemplary System transformation.

(4) FIG. 3 shows an exemplary Flowchart of the proposed VSA approach.

(5) FIG. 4 shows a sample VSI for different contingencies at bus 14 under normal load level.

(6) FIG. 5 shows an exemplary VSI for different contingencies at bus 14 under peak load level.

(7) FIG. 6 shows one embodiment computer that host the application.

DETAILED DESCRIPTION

(8) Methods to monitor the long-term voltage stability in power system are detailed next. The system provides a framework for monitoring the voltage stability margin considering N−1 transmission contingencies using sensitivity based Thevenin index (STI).

(9) Other Nomenclatures used herein include:

(10) i,j Index of buses.

(11) ΔV,Δθ Change of voltage magnitude and angle in two subsequent PMU measurements.

(12) Q.sup.min,Q.sup.max Minimum and maximum reactive power limits for a generator.

(13) V.sub.i,θ.sub.i Voltage magnitude and angle at bus i.

(14) f(θ,V) Active power injection expressions.

(15) g(θ,V) Reactive power injection expressions.

(16) The system works with power systems where buses may be of 3 types.

(17) 1. PQ bus—the real power |P| and reactive power |Q| are specified. It is also known as Load Bus. Generally, in a PQ bus, the generated real and reactive power will be assumed to be zero. However, power will be flowing out, thus, the injected real power and reactive power will be both negative. 2. PV bus—the real power |P| and the voltage magnitude |V| are specified. It is also known as Generator Bus. The real power and voltage are specified for buses that are generators. These buses have a constant power generation, controlled through a prime mover, and a constant bus voltage magnitude. 3. Slack bus—to balance the active and reactive power in the system. It is also known as the Reference Bus or the Swing Bus. The slack bus will serve as an angular reference for all other buses in the system, which is set to 0°. The voltage magnitude is also assumed to be 1 p.u. at the slack bus.
The slack bus provides or absorbs active and reactive power to and from the transmission system to provide for losses, since these variables are unknown until the final solution is established. The slack bus is the only bus for which the system reference phase angle is defined. From this, the various angular differences can be calculated in the power flow equations. If a slack bus is not specified, then a generator bus with maximum real power |P| acts as the slack bus. A given scheme can involve more than one slack bus. The most common formulation of the load flow problem specifies all input variables (PQ at loads, PV at generators) as deterministic values. Each set of specified values corresponds to one system state, which depends on a set of system conditions. When those conditions are uncertain, numerous scenarios must be analyzed.

(18) A sensitivity based Thevenin index is discussed next In FIG. 1, a two-bus system with a load bus i and the rest of the system represented by a Thevenin equivalence is provided. To estimate the Thevenin equivalence, i.e., Ēt.sub.h, Z.sub.th, at least two subsequent measurements are required with the assumption that the Thevenin equivalence remains the same during the time interval. Suppose we have two measurement sets at load bus i, i.e., (V.sub.i.sup.(1),Ī.sub.i.sup.(1),S.sub.i.sup.L.sup.(1)) and (V.sub.i.sup.(2),I.sub.i.sup.(2),S.sub.i.sup.L.sup.(1) (1+Δλ)), the following two equations hold:
Ē.sub.th=V.sub.i.sup.(1)+Z.sub.th.Math.Ī.sub.i.sup.(1)  (1)
Ē.sub.th=V.sub.i.sup.(2)+Z.sub.th.Math.Ī.sub.i.sup.(2)  (2)

(19) Thus, Z.sub.th can be calculated as:

(20) $\begin{array}{cc}{\stackrel{\_}{Z}}_{\mathrm{th}}=-\frac{{\stackrel{\_}{V}}_i^{\left(2\right)}-{\stackrel{\_}{V}}_i^{\left(1\right)}}{{\stackrel{\_}{I}}_i^{\left(2\right)}-{\stackrel{\_}{I}}_i^{\left(1\right)}}=-\frac{\Delta {\stackrel{\_}{V}}_i}{\Delta {\stackrel{\_}{I}}_i}& \left(3\right)\end{array}$

(21) Based on the voltage stability index, the LTI can be expressed as the ratio between Z.sub.th and the load impedance Z.sub.L:

(22) $\begin{array}{cc}\mathrm{LTI}=.Math.\frac{{\stackrel{\_}{Z}}_{\mathrm{th}}}{{\stackrel{\_}{Z}}_L}.Math.=.Math.\frac{\Delta {\stackrel{\_}{V}}_i}{{\stackrel{\_}{V}}_i^{\left(1\right)}}.Math.\frac{{\stackrel{\_}{I}}_i^{\left(1\right)}}{\Delta {\stackrel{\_}{I}}_i}.Math.& \left(4\right)\end{array}$

(23) By leveraging the polar form of the complex value and let V.sub.i.sup.(1)=Ve.sup.jθ,V.sub.i.sup.(2)=(V+ΔV)e.sup.j(θ+Δθ), we can write ΔV.sub.i into:

(24) $\begin{array}{cc}\begin{array}{c}\Delta {\stackrel{\_}{V}}_i=\left(V+\Delta V\right)e^{j\left(\theta +\Delta \theta \right)}-{\mathrm{Ve}}^{j\theta }\\ ={\stackrel{\_}{V}}_i^{\left(1\right)}\left(\left(1+\Delta V/V\right)e^{j\Delta \theta }-1\right)\end{array}& \left(5\right)\end{array}$

(25) In addition, ΔĪ, can be described by

(26) $\begin{array}{cc}\begin{array}{c}\Delta {\stackrel{\_}{I}}_i-{\stackrel{\_}{I}}_i^{\left(2\right)}-{\stackrel{\_}{I}}_i^{\left(1\right)}=\frac{{\stackrel{\_}{S}}_i^{L^{{\left(1\right)}^{*}}}\left(1+\Delta \lambda \right)}{{\mathrm{Ve}}^{-j\left(\theta +\Delta \theta \right)}}-\frac{{\stackrel{\_}{S}}_i^{L^{{\left(1\right)}^{*}}}}{{\mathrm{Ve}}^{\left(-j\theta \right)}}\\ ={\stackrel{\_}{I}}_i^{\left(1\right)}\left(\frac{1+\mathrm{\Delta \lambda }}{\left(1+\Delta V/V\right)e^{-j\theta }}-1\right)\end{array}& \left(6\right)\end{array}$

(27) With (5) and (6), the ratio between Zt.sub.h and Z.sub.L can be derived as:

(28) $\begin{array}{cc}\frac{{\stackrel{\_}{Z}}_{\mathrm{th}}}{{\stackrel{\_}{Z}}_L}=-\frac{\left(1+\Delta V/V\right).Math.\left(\left(1+\Delta V/V\right)-e^{-j\Delta \theta }\right)}{1+\mathrm{\Delta \lambda }-\left(1+\Delta V/V\right)e^{-j\Delta \theta }}& \left(7\right)\end{array}$

(29) Since Δθ is small, e.sup.−jΔθ≈1−jΔθ, the following equation can be obtained:

(30) $\begin{array}{cc}\frac{{\stackrel{\_}{Z}}_{\mathrm{th}}}{{\stackrel{\_}{Z}}_L}=-\left(1+\frac{\Delta V}{V}\right).Math.\frac{\left(\left(1+\Delta V/V\right)-\left(1-j\Delta \theta \right)\right)}{1+\mathrm{\Delta \lambda }-\left(1+\Delta V/V\right)\left(1-j\Delta \theta \right)}& \left(8\right)\end{array}$

(31) In (8), ΔV is small compared to

(32) $V,\frac{\Delta V}{V}$
can be neglected in the first term. Similarly, the product of ΔV and Δθ is also small and can be dropped from the denominator in the second term. The simplified version of (8) can be described as:

(33) $\begin{array}{cc}\frac{{\stackrel{\_}{Z}}_{\mathrm{th}}}{{\stackrel{\_}{Z}}_L}=-\frac{\left(\Delta V/V+j\Delta \theta \right)}{\Delta \lambda -\Delta V/V+j\Delta \theta }=-\frac{\left(\frac{\Delta V}{\Delta \lambda }.Math.\frac{1}{V}+j\frac{\Delta \theta }{\Delta \lambda }\right)}{1-\frac{\Delta V}{\Delta \lambda }.Math.\frac{1}{V}+j\frac{\Delta \theta }{\Delta \lambda }}& \left(9\right)\end{array}$

(34) Thus, the LTI can be formulated as a function of Δλ:

(35) $\begin{array}{cc}\mathrm{LTI}\left(\Delta \lambda \right)=.Math.\frac{{\stackrel{\_}{Z}}_{\mathrm{th}}}{{\stackrel{\_}{Z}}_L}.Math.=\sqrt{\frac{{\left(\frac{\Delta V}{\Delta \lambda }.Math.\frac{1}{V}\right)}^2+{\left(\frac{\Delta \theta }{\Delta \lambda }\right)}^2}{{\left(1-\frac{\Delta V}{\Delta \lambda }.Math.\frac{1}{V}\right)}^2+{\left(\frac{\Delta \theta }{\Delta \lambda }\right)}^2}}& \left(10\right)\end{array}$

(36) The assumption that the Thevenin equivalence remains unchanged during the time interval is only valid when the load increment Δλ is very small. The idea LTI, i.e., STI, can be achieved by using the limit of (10) with Δλ.fwdarw.0:

(37) 0$\begin{array}{cc}\mathrm{STI}=\underset{\Delta \lambda .fwdarw.0}{\lim }\mathrm{LTI}\left(\Delta \lambda \right)=\sqrt{\frac{{\left(\frac{dV}{d\lambda }.Math.\frac{1}{V}\right)}^2+{\left(\frac{d\theta }{d\lambda }\right)}^2}{{\left(1-\frac{dV}{d\lambda }.Math.\frac{1}{V}\right)}^2+{\left(\frac{d\theta }{d\lambda }\right)}^2}}& \left(11\right)\end{array}$

(38) To obtain the STI, two sensitivity terms dV/dλ and dθ/dλ should be calculated. The procedure can be interpreted as one predictor step of CPF when dλ is selected to be the continuation parameter, which is provided in equation (12).

(39) $\begin{array}{cc}\left[\begin{array}{ccc}{f}_{\theta }& f_V& f_{\lambda }^P\\ g_{\theta }& g_V& g_{\lambda }^Q\\ 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]& \left(12\right)\end{array}$

(40) In (12), f.sub.θ,f.sub.V,g.sub.θ,g.sub.V represent the partial derivative of the active and reactive power injection expressions with respect to bus voltage magnitudes and angles, i.e., the standard Jacobian matrix in power flow analysis at an operating point. The term f.sup.P.sub.λ,g.sup.P.sub.λ indicate the amount of active and reactive power injection changes as a function of dλ.

(41) Equation (12) should be evaluated at an operating point to obtain the sensitivities. However, the foreseen N−1 contingencies are hypothetic and there is no information regarding the operating point. To address this issue, the post-contingency states are estimated based on current operating point. Then equation (12) is evaluated at the estimated operating point to achieve the sensitivities and further the STI for each contingency.

(42) The preferred embodiment addresses the voltage stability assessment (VSA) considering contingency. The approach leverages the sensitivity based Thevenin index (STI) which involves evaluating the Jacobian matrix at current operating condition. Since the contingency case is hypothetic, there is no information regarding the operating condition after a foreseen contingency. The preferred embodiment first estimates the post-contingency operating point as well as possible PV-PQ transitions based on the current operating point. Then the STI for each contingency can be predicted using the estimated operating condition. Numerical results based on IEEE 14-bus system demonstrate the accuracy of the preferred embodiment in predicting the voltage stability margin under contingency. Moreover, the on-line implementation of the preferred embodiment is promising since it only involves solving several linear equations.

(43) A. Estimation of Post-Contingency Operating States

(44) The system estimates the operating point after a foreseen transmission contingency. The core idea is provided in FIG. 2.

(45) As shown in FIG. 2, system A represents the power network in base operating condition, where S.sub.i.sup.inj and S.sub.ij indicate the net power injection at bus i and line flow from bus i to bus j, respectively. Similar description applies to S.sub.j.sup.inj and S.sub.ji. The difference between system A and system B lies in two parts: 1) system B does not have branch ij; 2) the net power injections at bus i and j of system B are changed to be S.sub.i.sup.inj−S.sub.ij and S.sub.j.sup.inj−S.sub.ji, respectively. System A and B are equivalent and have the same power flow results. The post-contingency system, i.e., branch ij in outage, is provided as system C. Compared to system B, the power injections at bus i and j are changed by the amount of S.sub.ij and S.sub.ji. Thus, the post-contingency operating states can be calculated as:

(46) $\begin{array}{cc}\left[\begin{array}{c}{\theta }_c\\ V_c\end{array}\right]={\left[\begin{array}{c}\tilde{\theta }\\ \tilde{V}\end{array}\right]\left[\begin{array}{cc}{f}_{\theta }& f_V\\ g_{\theta }& g_V\end{array}\right]}^{-1}\Delta F_{\mathrm{inj}}& \left(13\right)\end{array}$
where ({acute over (θ)},{tilde over (V)}) is the operating point under base case; (θ.sub.c,V.sub.c) is the operating point following a contingency; ΔF.sub.inj is the vector containing the active and reactive bus injection changes which are numerically equal to the active and reactive line flow on the outage branch. In one embodiment, the Jacobian matrix in (13) is constructed based on the topology of system B, i.e., the contingency line ij is in outage

(47) When a severe transmission contingency occurs, some generators may reach their reactive power limits and the corresponding bus type will change from PV to PQ. To both accurately estimate the post-contingency states and the voltage stability margin, it is paramount to predict the possible PV-PQ transitions. The piecewise linear sensitivity method is used to predict the Q limits violation. First, a K factor which reflects the severity of reactive power violation at a PV bus is defined as

(48) $\begin{array}{cc}K=\{\begin{array}{c}\frac{Q^{\max }-{\tilde{Q}}^{\left(l\right)}}{Q^{\left(l+1\right)}-{\tilde{Q}}^{\left(l\right)}},Q^{\left(l+1\right)}\ge Q^{\max }\\ 1,Q^{\min }\le Q^{\left(l+1\right)}\le Q^{\max }\\ \frac{{\tilde{Q}}^{\left(l\right)}-Q^{\min }}{{\tilde{Q}}^{\left(l\right)}-Q^{\left(l+1\right)}},Q^{l+1}\le Q^{\min }\end{array}& \left(14\right)\end{array}$

(49) The complete process is provided as below:

(50) 1) Set l=0 and obtain the power flow results, {tilde over (θ)}.sup.(l), {tilde over (V)}.sup.(l), {tilde over (Q)}.sup.(l), in base operating condition, i.e., before the contingency.

(51) 2) Apply ΔF.sub.inj.sup.(l) to system l and leverage equation (13) to predict the system l+1 states θ.sup.(l+1), V.sup.(l+1), Q.sup.(l+1).

(52) 3) Calculate K factor for all the PV buses. If the values of all K are equal to 1 or there is no PV bus, go to step 6); otherwise, go to step 4).

(53) 4) Find the PV bus with the lowest K, say bus i with K.sub.i. Apply a portion of ΔF.sub.inj.sup.(l), i.e., K.sub.iΔF.sub.inj.sup.(l) to system l to predict the intermediate system l+1 states {tilde over (θ)}.sup.(l+1),{tilde over (V)}.sup.(l+1), {tilde over (Q)}.sup.(l+1).

(54) 5) Change bus i to PQ bus by fixing its reactive power generation to its limit. Let Δ F.sub.inj.sup.(l+1)=(1−K.sub.i)ΔF.sub.inj.sup.(l), set l=l+1 and go to step 2).

(55) 6) Output the predicted final states θ.sup.(l+1),V.sup.(l+1), Q.sup.(l+1).

(56) B. Implementation of the VSA Embodiment

(57) With the post-contingency states, the STI under each transmission contingency can be predicted. However, equation (12) should be modified to calculate the corresponding sensitivities when there are generators reaching their limits. As an example, suppose bus i is transformed from PV bus to PQ bus. Equation (12) is modified as follows

(58) $\begin{array}{cc}\left[\begin{array}{cccc}{f}_{\theta }& f_V& f_{V_i}& f_{\lambda }^P\\ g_{\theta }& g_V& g_{V_i}& g_{\lambda }^Q\\ g_{\theta }^i& g_V^i& g_{V_i}^i& g_{\lambda }^{Q,i}\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}d\theta \\ \mathrm{dV}\\ {\mathrm{dV}}_i\\ d\lambda \end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\\ 1\end{array}\right]& \left(15\right)\end{array}$

(59) Since V.sub.i becomes a variable, one column and row should be added to the previous Jacobian. For the added row, g.sub.θ.sup.i,g.sub.V.sup.i,g.sub.V.sub.i.sup.i, represent the partial derivative of reactive power injection expressions at bus i with respect to θ, V and V.sub.i. In one embodiment, if bus i does not have any load, g.sub.λ.sup.Q,i is equal to zero because the reactive limit of the generator at that bus is reached; otherwise, the value of g.sub.λ.sup.Q,i depends on how will the reactive load varies as dλ changes. Equation (15) can be easily extended for any number of PV-PQ bus transitions.

(60) The flowchart of the preferred embodiment is depicted in FIG. 3. The procedures are described as below (note that the number shown in the figure is clearly defined in the procedures):

(61) 1) Obtain the current operating states {tilde over (θ)},{tilde over (V)} [110].

(62) 2) Estimate the post-contingency operating point θ.sub.c,V.sub.c for every considered contingency with the method presented above [112].

(63) 3) If there is no PV-PQ transition following the contingency [114], calculate the corresponding sensitivities using equation (12) [116]; otherwise, obtain the sensitivities using equation (15) [118].

(64) 4) Calculate the STI for each considered contingency with equation (11) [120,122].

(65) This system determines the voltage stability assessment (VSA) considering N−1 contingency using sensitivity based Thevenin index (STI) which involves evaluating the Jacobian matrix at current operating condition. Since the N−1 contingency case is hypothetic, there is no information regarding the operating condition after a foreseen contingency. The preferred embodiment first estimates the post-contingency operating point as well as possible PV-PQ transitions based on the current operating point. Then the STI for each contingency can be predicted using the estimated operating condition. Numerical results based on IEEE 14-bus system demonstrate the accuracy of the preferred embodiment in predicting the voltage stability margin under contingency. Moreover, the on-line implementation of the preferred embodiment is promising since it only involves solving several linear equations.

(66) In one embodiment, the calculation of STI for each contingency and load bus is independent so the computational speed can be enhanced if parallelization technique is leveraged.

(67) Exemplary case studies are presented next. The IEEE-14 bus system provided in MATPOWER software package is selected to test the performance of the proposed VSA approach. It includes five generators, 14 buses and 20 transmission branches. Since the data only gives one load pattern, we assume that load profile as normal load level, whose active and reactive loads are 259 MW and 73.5 MVAr, respectively. It is also assumed that under peak load level, the active and reactive loads are 1.2 times the values in normal load level. We consider all the N−1 branch outages except the following two branches: line 1-2 is a key transmission element whose loss will cause the divergence of power flow under normal load pattern; line 7-8 only connects a generator bus in one end whose loss will cause islanding system.

(68) First, we implement the preferred embodiment under normal load condition. Table I provides the VSI, i.e., |Z.sub.th|/|Z.sub.L|, at every load bus for two critical contingencies. The STI* denotes the benchmark results, which are calculated based on exact post-contingency power flow solutions. We compare the accuracy of the preferred embodiment with existing solutions where the VSI is based on the voltage magnitude deviation. To accommodate the comparison, we shift the VSI to |Z.sub.th|/|Z.sub.L|. Nevertheless, we leverage the same method as done conventionally to identify the post-contingency TE. The last row in Table I, i.e., σ, denotes the average relative error (averaged among all load buses) of the two approaches with respect to the benchmark results. As can be observed from Table I, under the normal load level, the preferred embodiment provides more accurate results than conventional systems for every load bus. The average errors of the proposed method are only 0.31% and 0.96% for contingency 1-5 and 5-6, which are far less than 5.32% and 5.41% given by the conventional system.

(69) FIG. 4 depicts the VSIs for all the 18 contingencies at bus 14, which is the critical load bus. The VSI at bus 14 for the base case, i.e., without outage, is 0.2771. It can be seen that the majority of the predicted post-contingency VSI are higher than 0.2771, which indicates that the voltage stability margin will decrease if those contingencies occur. In one embodiment, if the VSI exceeds 0.45, alarms will be sent to the control center for situational awareness. Under the normal load level, no alarms will be triggered for bus 14 since the predicted post-contingency VSIs are all below 0.45. Also, the preferred embodiment outperforms conventional systems for better estimation of the post-contingency VSI by using the novel framework for monitoring the voltage stability margin considering N−1 transmission contingencies using STI and solving several linear equations, which makes it promising for on-line implementation.

(70) TABLE-US-00001 TABLE I COMPARISON OF VSI UNDER NORMAL LOAD LEVEL Contingency Line (i-j) 1-5 5-6 Bus STI* Conventional Preferred Embodiment STI* Conventional Preferred Embodiment  4 0.2874 0.2704 0.2865 0.1976 0.1948 0.1969  5 0.2780 0.2591 0.2770 0.1570 0.1551 0.165   9 0.3590 0.3417 0.3579 0.3577 0.3415 0.3544 10 0.3638 0.3461 0.3627 0.3763 0.3565 0.3725 11 0.3695 0.3506 0.3684 0.4085 0.3788 0.4036 12 0.3821 0.3622 0.3809 0.4443 0.4045 0.4381 13 0.3797 0.3601 0.3785 0.4401 0.4024 0.4342 14 0.3784 0.3606 0.3773 0.4092 0.3845 0 4046 σ (%) — 5.32   0.31   — 5.41   0.96  

(71) Next, we test the preferred embodiment for the peak load condition. In Table II, the VSI at every load bus for two critical contingencies are provided. Still, the preferred embodiment gives more accurate results in predicting the post-contingency VSI. For contingency 1-5, the average error σ is 24.82% by using the conventional system. This value decreases to 4.25% by leveraging the preferred embodiment.

(72) TABLE-US-00002 TABLE II COMPARISON OF VSI UNDER PEAK LOAD LEVEL Contingency Line (i-j) 1-5 5-6 Bus STI* Conventional Preferred Embodiment STI* Conventional Preferred Embodiment  4 0.4694  0.3345 0.4474 0.3674  0.2081  0.3260  5 0.4585  0.3111 0.4358 0.3054  0.1642  0.2692  9 0.5705  0.4389 0.5471 0.6055  0.4229  0.5395 10 0.5759  0.4441 0.5524 0.6323  0.3477  0.5642 11 0.5719  0.4377 0.5484 0.6739  0.3784  0.6022 12 0.5795  0.4441 0.5557 0.7154  0.3555  0.6413 13 0.5817  0.4473 0.5581 0.7107  0.3566  0.6370 14 0.5970  0.4666 0.5734 0.6748  0.3669  0.6047 σ (%) — 24.82   4.25   — 48.25   10.82  

(73) FIG. 5 provides the VSIs for different contingencies at bus 14 under peak load level. As can be seen from FIG. 5, for every contingency, the VSI calculated by the preferred embodiment is closer to the benchmark results than the VSI given conventional system. The largest difference between the benchmark results and the predicted VSIs occurs in contingency 5-6, where the benchmark VSI is 0.6748 and an alarm should be signaled. The preferred embodiment calculates the VSI as 0.6047, which still allows the generation of alarm. However, the conventional system outputs 0.3669 for the VSI, which will not signal any alarm. In addition, the VSI for contingency 1-5, 2-3, 2-4, 5-6 and 7-9 are above 0.45. Thus, the corresponding alarms are generated to trigger some preventive controls, which improves the system operator's situational awareness on the criticality of the current operating condition.

(74) The invention may be implemented in hardware, firmware or software, or a combination of the three. Preferably the invention is implemented in a computer program executed on a programmable computer having a processor, a data storage system, volatile and non-volatile memory and/or storage elements, at least one input device and at least one output device.

(75) By way of example, a block diagram of a computer to support the system is discussed next in FIG. 7. The computer preferably includes a processor, random access memory (RAM), a program memory (preferably a writable read-only memory (ROM) such as a flash ROM) and an input/output (I/O) controller coupled by a CPU bus. The computer may optionally include a hard drive controller which is coupled to a hard disk and CPU bus. Hard disk may be used for storing application programs, such as the present invention, and data. Alternatively, application programs may be stored in RAM or ROM. I/O controller is coupled by means of an I/O bus to an I/O interface. I/O interface receives and transmits data in analog or digital form over communication links such as a serial link, local area network, wireless link, and parallel link. Optionally, a display, a keyboard and a pointing device (mouse) may also be connected to I/O bus. Alternatively, separate connections (separate buses) may be used for I/O interface, display, keyboard and pointing device. Programmable processing system may be preprogrammed or it may be programmed (and reprogrammed) by downloading a program from another source (e.g., a floppy disk, CD-ROM, or another computer).

(76) Each computer program is tangibly stored in a machine-readable storage media or device (e.g., program memory or magnetic disk) readable by a general or special purpose programmable computer, for configuring and controlling operation of a computer when the storage media or device is read by the computer to perform the procedures described herein. The inventive system may also be considered to be embodied in a computer-readable storage medium, configured with a computer program, where the storage medium so configured causes a computer to operate in a specific and predefined manner to perform the functions described herein.

(77) The invention has been described herein in considerable detail in order to comply with the patent Statutes and to provide those skilled in the art with the information needed to apply the novel principles and to construct and use such specialized components as are required. However, it is to be understood that the invention can be carried out by specifically different equipment and devices, and that various modifications, both as to the equipment details and operating procedures, can be accomplished without departing from the scope of the invention itself.

(78) In sum, our approach based on STI addresses the VSA considering a series of foreseen transmission N−1 contingencies. The wide area measurements based STI is interpreted as an idea LTI and can be calculated by evaluating the Jacobian matrix at current operating condition. To extend its application in monitoring the voltage stability margin for N−1 contingencies, the post-contingency operating point as well as the PV-PQ transitions are first estimated based on the current operating condition. Then the estimated operating point is used to calculate the STI for each N−1 contingency. Numerical case studies based on the IEEE 14-bus system demonstrate the effectiveness of the preferred embodiment and its benefits in the situational awareness enhancement.

(79) The foregoing has described a measurement-based voltage stability monitoring and control scheme. While specific embodiments of the present invention have been described, it will be apparent to those skilled in the art that various modifications thereto can be made without departing from the spirit and scope of the invention. Accordingly, the foregoing description of the preferred embodiment of the invention and the best mode for practicing the invention are provided for the purpose of illustration only and not for the purpose of limitation.