Inverse Stability Problem and Applications to Renewables Integration

TL;DR

This paper introduces the inverse stability problem in power systems, focusing on characterizing the set of operating conditions a grid converges to under renewable fluctuations, with applications in stability assessment and control.

Contribution

It proposes a novel inverse stability framework using quadratic energy function approximations, addressing a rarely studied problem in power system stability analysis.

Findings

Estimate of inverse stability region using quadratic energy functions

Application to robust stability assessment with renewables

Application to stability-constrained optimal power flow

Abstract

In modern power systems, the operating point, at which the demand and supply are balanced, may take different values due to changes in loads and renewable generation levels. Understanding the dynamics of stressed power systems with a range of operating points would be essential to assuring their reliable operation, and possibly allow higher integration of renewable resources. This letter introduces a non-traditional way to think about the stability assessment problem of power systems. Instead of estimating the set of initial states leading to a given operating condition, we characterize the set of operating conditions that a power grid converges to from a given initial state under changes in power injections and lines. We term this problem as "inverse stability", a problem which is rarely addressed in the control and systems literature, and hence, poorly understood. Exploiting quadratic approximations of the system's energy function, we introduce an estimate of the inverse stability region. Also, we briefly describe three important applications of the inverse stability notion: (i) robust stability assessment of power systems w.r.t. different renewable generation levels, (ii) stability-constrained optimal power flow (sOPF), and (iii) stability-guaranteed corrective action design.

Tables

Table 1. TABLE I: Bus voltages, mechanical inputs, and static loads.

Node	V (p.u.)	$P_k$ (p.u.)
1	1.0284	3.6466
2	1.0085	4.5735
3	0.9522	3.8173
4	1.0627	-3.4771
5	1.0707	-3.5798
6	1.0749	-3.3112
7	1.0490	-0.5639
8	1.0579	-0.5000
9	1.0521	-0.6054

Equations

$$m_k\ddot{{\delta }_k}+d_k\dot{{\delta }_k}+\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}\sin ({\delta }_k-{\delta }_j)=$$

$$d_k\dot{{\delta }_k}+\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}\sin ({\delta }_k-{\delta }_j)=$$

$$\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}\sin {\delta }_{kj}=$$

$$\displaystyle E{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}(\delta,\delta^{*})}=\sum_{k\in\mathcal{G}}\frac{m_{k}\dot{\delta}_{k}^{2}}{2}+\sum_{\{k,j\}\in\mathcal{E}}\int_{\delta_{kj}^{*}}^{\delta_{kj}}a_{kj}(\sin\xi-\sin\delta_{kj}^{*})d\xi$$

$$\displaystyle\dot{E}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}(\delta,\delta^{*})}$$

$$=\sum _{k\in \mathcal{G}}{\dot{\delta }}_k(P_k-d_k\dot{{\delta }_k}-\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}\sin ({\delta }_k-{\delta }_j))$$

$$+\sum _{k\in \mathcal{L}}{\dot{\delta }}_k(P_k-d_k\dot{{\delta }_k}-\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}\sin ({\delta }_k-{\delta }_j))$$

$$+\sum _{\{k,j\}\in \mathcal{E}}{a}_{kj}(\sin {\delta }_{kj}-\sin {\delta }_{kj}^{\ast })({\dot{\delta }}_k-\dot{{\delta }_j})$$

$$=-\sum _{k\in \mathcal{N}}{d}_k(\dot{{\delta }_k}{)}^2\le 0,$$

$$g(\xi -{\delta }_{kj}^{\ast }{)}^2\le (\xi -{\delta }_{kj}^{\ast })(\sin \xi -\sin {\delta }_{kj}^{\ast })\le (\xi -{\delta }_{kj}^{\ast }{)}^2.$$

$$E(\delta ,{\delta }^{\ast })$$

$$\begin{array}{rl}D(\delta ,{\delta }^{\ast })& =g\sum _{\{k,j\}\in \mathcal{E}}{\underline{a}}_{kj}\frac{({\delta }_{kj}-{\delta }_{kj}^{\ast }{)}^2}{2},\\ F(\delta ,{\delta }^{\ast })& =\sum _{k\in \mathcal{G}}\frac{m_k{\dot{\delta }}_k^2}{2}+\sum _{\{k,j\}\in \mathcal{E}}{\bar{a}}_{kj}\frac{({\delta }_{kj}-{\delta }_{kj}^{\ast }{)}^2}{2}.\end{array}$$

$$D(\delta ,{\delta }^{\ast })\le E(\delta ,{\delta }^{\ast })\le F(\delta ,{\delta }^{\ast }),\forall \delta \in \mathcal{P},{\delta }^{\ast }\in \Lambda .$$

$$\mathcal{B}({\delta }_0)=\{\delta :F({\delta }_0,\delta )\le R({\delta }_0)/4\}.$$

$$\Vert L^{†}p{\Vert }_{\mathcal{E},\infty }\le \sin \lambda ,$$

$$\min c(P)$$

$$s.t.$$

$$Q(V,\delta )=Q$$

$$\underline{V}\le V\le \bar{V}$$

$$\underline{S}\le |S(V,\delta )|\le \bar{S}$$

$$\min c(P)$$

$$s.t.$$

$$\underline{V}\le V\le \bar{V}$$

$$\underline{S}\le |S(V,\delta )|\le \bar{S}$$

$$\delta \in \mathcal{A}({\delta }_0).$$

$$E(M,{\delta }^{\ast })+E({\delta }_0,{\delta }^{\ast })\ge D(M,{\delta }^{\ast })+D({\delta }_0,{\delta }^{\ast })$$

$$=g\sum _{\{k,j\}\in \mathcal{E}}{\underline{a}}_{kj}\frac{({\delta }_{M_{kj}}-{\delta }_{kj}^{\ast }{)}^2+({\delta }_{0_{kj}}-{\delta }_{kj}^{\ast }{)}^2}{2}$$

$$\ge g\sum _{\{k,j\}\in \mathcal{E}}{\underline{a}}_{kj}\frac{({\delta }_{M_{kj}}-{\delta }_{0_{kj}}{)}^2}{4}=\frac{D({\delta }_0,M)}{2}.$$

$$E(M,{\delta }^{\ast })>R({\delta }_0)/4>E({\delta }_0,{\delta }^{\ast }).$$

Full text

I Introduction

Renewable generations, e.g., wind and solar, are increasingly installed into electric power grids to reduce $\ce{CO_{2}}$ emission from the electricity generation sector. Yet, their natural intermittency presents a major challenge to the delivery of consistent power that is necessary for today’s grid operation, in which generation must instantly meet load. Also, the inherently low inertia of renewable generators limits the grid’s controllability and makes it easy for the grid to lose its stability. The existing power grids and management tools were not designed to deal with these new challenges. Therefore, new stability assessment and control design tools are needed to adapt to the changes in architecture and dynamic behavior expected in the future power grids.

Transient stability assessment of power system certifies that the system state converges to a stable operating condition after the system experiences large disturbances. Traditionally, this task is handled by using either time domain simulation (e.g., [1]), or by utilizing the energy method (e.g., [2, 3]) and the Lyapunov function method (e.g., [4]) to estimate the stability region of a given equilibrium point (EP), i.e., the set of initial states from which the system state will converge to that EP. In modern renewable power grids, the operating point may take different values under the real-time clearing of electricity markets, intermittent renewable generations, changing loads, and external disturbances. Dealing with the situation when the EP can change over a wide range makes the transient stability assessment even more technically difficult and computationally cumbersome.

In this letter, rather than considering the classical stability assessment problem, we formulate the inverse stability assessment problem. This problem concerns with estimating the region around a given initial state $\delta_{0}$, called “inverse stability region” $\mathcal{A}(\delta_{0})$, so that whenever the power injections or power lines change and lead to an EP in $\mathcal{A}(\delta_{0})$, the system state will converge from $\delta_{0}$ to that EP. Indeed, the convergence from $\delta_{0}$ to an EP is guaranteed when the system’s energy function is bounded under some threshold [2, 3]. In [5], we observed that if the EP is in the interior of the set $\mathcal{P}$ characterized by phasor angular differences smaller than $\pi/2,$ then the nonlinear power flows can be strictly bounded by linear functions of angular differences. Exploiting this observation, we show that the energy function of power system can be approximated by quadratic functions of the EP and the system state, and from which we obtain an estimate of the inverse stability region.

The remarkable advantage of the inverse stability certificate is making it possible to exploit the change in EP to achieve useful dynamical properties. We will briefly discuss three applications of this certificate, which are of importance to the integration of large-scale renewable resources:

• Robust stability assessment: For a typical power system composed of several components and integrated with different levels of renewable generations, there are many contingencies that need to be reassessed on a regular basis. Most of these contingencies correspond to failures of relatively small and insignificant components, so the post-fault dynamics is probably transiently stable. Therefore, most of the computational effort is spent on the analysis of non-critical scenarios. This computational burden could be greatly alleviated by a robust transient stability assessment toolbox that could certify the system’s stability w.r.t. a broad range of uncertainties. In this letter, we show that the inverse stability certificate can be employed to assess the transient stability of power systems for various levels of power injections.
• Stability-constrained OPF: Under large disturbances, a power system with an operating condition derived by solving the conventional OPF problem may not survive. It is therefore desirable to design operating conditions so that the system can withstand large disturbances. This can be carried out by incorporating the transient stability constraint into OPF together with the normal voltage and thermal constraints. Though this problem was discussed in the literature (e.g., [6]), there is no way to precisely formulate and solve the stability-constrained OPF problem because transient stability is a dynamic concept and differential equations are involved in the stability constraint. Fortunately, the inverse stability certificate allows for a natural incorporation of the stability constraint into the OPF problem as a static constraint of placing the EP in a given set.
• Stability-guaranteed corrective actions: Traditional protection strategies focus on the safety of individual components, and the level of coordination among component protection systems is far from perfect. Also, they do not take full advantage of the new flexible and fast electronics resources available in modern power systems, and largely rely on customer-harmful actions like load shedding. These considerations motivated us to coordinate widespread flexible electronics resources as a system-level customer-friendly corrective action with guaranteed stability [7]. This letter presents a unconventional control way in which we relocate the operating point, by appropriately redispatching power injections, to attract the emergency state and stabilize the power systems under emergency situations.

II Inverse stability problem of power systems

In this letter, we utilize the structure-preserving model to describe the power system dynamics [8]. This model naturally incorporates the dynamics of the generators’ rotor angle and the response of load power output to frequency deviation. Mathematically, the grid is described by an undirected graph $\mathcal{A}(\mathcal{N},\mathcal{E}),$ where $\mathcal{N}=\{1,2,\dots,|\mathcal{N}|\}$ is the set of buses and $\mathcal{E}\subseteq\mathcal{N}\times\mathcal{N}$ is the set of transmission lines $\{k,j\},k,j\in\mathcal{N}$. Here, $|A|$ denotes the number of elements of set $A$. The sets of generator buses and load buses are denoted by $\mathcal{G}$ and $\mathcal{L}$. We assume that the grid is lossless with constant voltage magnitudes $V_{k},k\in\mathcal{N},$ and the reactive powers are ignored. Then, the grid’s dynamics is described by [8]:

[TABLE]

where equation (1a) applies at the dynamics of generator buses and equation (1b) applies at the dynamics of load buses. Here $a_{kj}=V_{k}V_{j}B_{kj},$ where $B_{kj}$ is the (normalized) susceptance of the transmission line $\{k,j\}$ connecting the $k^{th}$ bus and $j^{th}$ bus. $\mathcal{N}_{k}$ is the set of neighboring buses of the $k^{th}$ bus (see [9] for more details). Let $\delta(t)=[\delta_{1}(t)\;...\;\delta_{|\mathcal{N}|}(t)\;\dot{\delta}_{1}(t)\;...\;\dot{\delta}_{|\mathcal{N}|}(t)]^{\top}$ be the state of the system (1) at time $t$ (for simplicity, we will denote the system state by $\delta$). Note that Eqs. (1) are invariant under any uniform shift of the angles $\delta_{k}\to\delta_{k}+c.$ However, the state $\delta$ can be unambiguously characterized by the angle differences $\delta_{kj}=\delta_{k}-\delta_{j}$ and the frequencies $\dot{\delta}_{k}.$

Normally, a power grid operates at an operating condition of the pre-fault dynamics. Under the fault, the system evolves according to the fault-on dynamics. After some time period, the fault is cleared or self-clears, and the system is at the so-called fault-cleared state $\delta_{0}$ (the fault-cleared state is usually estimated by simulating the fault-on dynamics, and hence, is assumed to be known). Then, the power system experiences the so-called post-fault dynamics. The transient stability assessment certifies whether the post-fault state converges from $\delta_{0}$ to a stable EP $\delta^{*}.$ Mathematically, the operating condition $\delta^{*}=[\delta_{1}^{*}\;...\;\delta_{|\mathcal{N}|}^{*}\;0\;...\;0]^{\top}$ is a solution of the power-flow like equations:

[TABLE]

With renewable generations or under power redispatching, the power injections $P_{k}$ take different values. Also, the couplings $a_{kj}$ can be changed by using the FACTS devices. Assume $\underline{a}_{kj}\leq a_{kj}\leq\bar{a}_{kj}.$ In those situations, the resulting EP $\delta^{*}$ also takes different values. Therefore, we want to characterize the region of EPs so that the post-fault state always converges from a given initial state $\delta_{0}$ to the EP whenever the EP is in this region. Though the EP can take different values, it is assumed to be fixed in each transient stability assessment because the power injections and couplings can be assumed to be unchanged in the very fast time scale of transient dynamics (i.e., 1 to 10 seconds). We consider the following problem:

• •

Inverse (Asymptotic) Stability Problem: Consider a given initial state $\delta_{0}.$ Assume that power injections and the line susceptances can take different values. Estimate the region of stable EPs so that the state of the system (1) always converges from $\delta_{0}$ to the EP in this region.

This problem will be addressed with the inverse stability certificate to be presented in the next section.

III Energy function and inverse stability certificate

III-A Stability assessment by using energy function

Before introducing the inverse stability certificate addressing the inverse stability problem in the previous section, we present a normal stability certificate for system with the fixed power injections and line parameters. For the power system described by Eqs. (1), consider the energy function:

[TABLE]

Then, along every trajectory of (1), we have

[TABLE]

in which the last equation is obtained from (2). Hence, $E(\delta,\delta^{*})$ is always decreasing along every trajectory of (1).

Consider the set $\mathcal{P}$ defined by $|\delta_{kj}|\leq\pi/2,\{k,j\}\in\mathcal{E},$ and the set $\Phi=\{\delta\in\mathcal{P}:E(\delta,\delta^{*})<E_{\min}(\delta^{*})\},$ where $E_{\min}(\delta^{*})=\min_{\delta\in\partial\mathcal{P}}E(\delta,\delta^{*})$ and $\partial\mathcal{P}$ is the boundary of $\mathcal{P}.$ $\Phi$ is invariant w.r.t. (1), and bounded as the state $\delta$ is characterized by the angle differences and the frequencies. Though $\Phi$ is not closed, the decrease of $E(\delta,\delta^{*})$ inside $\Phi$ assures the limit set to be inside $\Phi.$ As such, we can apply the LaSalle’s Invariance Principle and use a proof similar to that of Theorem 1 in [5] to show that, if $\delta_{0}$ is inside $\Phi$ then the system state will only evolve inside this set and eventually converge to $\delta^{*}.$ So, to check if the system state converges from $\delta_{0}\in\mathcal{P}$ to $\delta^{*},$ we only need to check if $E(\delta_{0},\delta^{*})<E_{\min}(\delta^{*}).$

III-B Inverse stability certificate

For a given initial state $\delta_{0}\in\mathcal{P}$, we will construct a region surrounding it so that whenever the operating condition $\delta^{*}$ is in this region then $E(\delta_{0},\delta^{*})<E_{\min}(\delta^{*}).$ Hence, the grid state will converge from $\delta_{0}$ to $\delta^{*}$ according to the stability certificate in Section III-A. Indeed, we establish quadratic bounds of the energy function for every $\delta^{*}$ in the set $\Lambda$ defined by inequalities $|\delta_{kj}|\leq\lambda<\pi/2,\forall\{k,j\}\in\mathcal{E}.$ Let $g=\dfrac{1-\sin\lambda}{\pi/2-\lambda}>0.$ In [5], we observed that for $\delta^{*}\in\Lambda,\xi\in\mathcal{P},$

[TABLE]

Hence, for all $\delta^{*}\in\Lambda,\delta\in\mathcal{P},$ we have

[TABLE]

Define the following functions

[TABLE]

Using (5) and (6), we can bound the energy function as

[TABLE]

For a given initial state $\delta_{0}$ inside the set $\mathcal{P},$ we calculate the “distance” from this initial state to the boundary of the set $\mathcal{P}:R(\delta_{0})=\min_{\delta\in\partial\mathcal{P}}D(\delta_{0},\delta).$ Let $\mathcal{B}(\delta_{0})$ be the neighborhood of $\delta_{0}$ defined by

[TABLE]

The following is our main result regarding inverse stability of power system, as illustrated in Fig. 1. **

Theorem 1

Consider a given initial state $\delta_{0}$ inside the set $\mathcal{P}.$ Assume that the EP of the system takes different values in the set $\mathcal{A}(\delta_{0})=\Lambda\cap\mathcal{B}(\delta_{0})$, where the set $\mathcal{B}(\delta_{0})$ is defined as in (9). Then, the system state always converges from the given initial state $\delta_{0}$ to the EP.

Proof: See Appendix VI. $\qquad\square$

Remark 1

In this paper, we limit the grid to be described by the simplified model (1) which captures the dynamics of the generators’ rotor angle and the response of load power output to frequency deviation. More realistic models should take into account voltage variations, reactive powers, dynamics of the rotor flux, and controllers (e.g., droop controls and power system stabilizers). It should be noted that the results in this paper is extendable to more realistic models. The reason is that all the key results in this paper rely on the analysis of the system’s energy function, in which we combine the energy function-based transient stability analysis with the quadratic bounds of the energy function. In high-order models, the energy function is more complicated, yet it still can be bounded by quadratic functions. Hence, we can combine the energy function-based transient stability analysis of the more realistic models in [2, 3, 4] with the quadratic bounds of the energy function to extend the results in this paper to these higher-order models. This approach may also extend to other higher-order models in the port-Hamiltonian formulation (e.g., [10, 11]) by applying the appropriate approximation on the Lyapunov functions established for these models.

IV Applications of Inverse Stability Certificate

IV-A Robust stability assessment

The robust transient stability problem that we consider involves situations where there is uncertainty in power injections $P_{k},$ e.g., due to intermittent renewable generations. Formally, for a given fault-cleared state $\delta_{0},$ we need to certify the transient stability of the post-fault dynamics described by (1) with respect to fluctuations of the power injections, which consequently lead to different values of the post-fault EP $\delta^{*}$ as a solution of the power flow equations (2). Therefore, we consider the following robust stability problem [5]:

• •

Robust stability assessment: Given a fault-cleared state $\delta_{0},$ certify the transient stability of (1) w.r.t. a set of stable EPs $\delta^{*}$ resulted from different levels of power injections.

Utilizing the inverse stability certificate, we can assure robust stability of renewable power systems whenever the resulting EP is inside the set $\mathcal{A}(\delta_{0})=\Lambda\cap\mathcal{B}(\delta_{0}).$ To check that the EP is in the set $\Lambda,$ we can apply the criterion from [12], which states that the EP will be in the set $\Lambda$ if the power injections $p=[P_{1},...,P_{|\mathcal{N}|}]^{T}$ satisfy

[TABLE]

where $L^{\dagger}$ is the pseudoinverse of the network Laplacian matrix and the norm $\|x\|_{\mathcal{E},\infty}$ is defined by $\|x\|_{\mathcal{E},\infty}=\max_{\{i,j\}\in\mathcal{E}}|x(i)-x(j)|.$ On the other hand, some similar sufficient condition could be developed so that we can verify that the EP is in the set $\mathcal{B}(\delta_{0})$ by checking the power injections. This will help us certify robust stability of the system by only checking the power injections.

IV-B Stability-constrained OPF

Stability-constrained OPF problem concerns with determining the optimal operating condition with respect to the voltage and thermal constraints, as well as the stability constraint. While the voltage and thermal constraints are well modeled via algebraic equations or inequalities, it is still an open question as to how to include the stability constraint into OPF formulation since stability is a dynamic concept and differential equations are involved [6].

Mathematically, a standard OPF problem is usually stated as follows (refer to [6] for more detailed formulation):

[TABLE]

where $c(P)$ is a quadratic cost function, the decision variables $P$ are typically the generator scheduled electrical power outputs, the equality constraints (12)-(13) stand for the power flow equations, and the inequality constraints (14)-(15) stand for the voltage and thermal limits of branch flows through transmission lines and transformers. Assume that the stability constraint is to make sure that the system state will converge from a given fault-cleared state $\delta_{0}$ to the designed operating condition, and that the reactive power is negligible. With the inverse stability certificate, the stability constraint can be relaxed and formulated as $\delta\in\mathcal{A}(\delta_{0}).$ Basically, the inverse stability certificate transforms the dynamic problem of stability into a static problem of placing the prospective EP into a set. In summary, we obtain a relaxation of the stability-constrained OPF problem as follows:

[TABLE]

Solution of this optimization problem in an optimal operating condition at which the cost function is minimized and the voltage/thermal constraints are respected. Furthermore, the stability constraint is guaranteed by the inverse stability certificate, and the system state is ensured to converge from the fault-cleared $\delta_{0}$ to the operating condition.

IV-C Emergency control design

Another application of the inverse stability certificate, that will be detailed in this section, is designing stability-guaranteed corrective actions that can drive the post-fault dynamics to a desired stability regime. As illustrated in Fig. 2, for a given fault-cleared state $\delta_{0}$, by applying the inverse stability certificate, we can appropriately dispatch the power injections $P_{k}$ to relocate the EP of the system so that the post-fault dynamics can be attracted from the fault-cleared state $\delta_{0}$ through a sequence of EPs $\delta^{*}_{1},...,\delta^{*}_{\ce N}$ to the desired EP $\delta^{*}_{\bf desired}$. In other words, we subsequently redispatch the power injections so that the system state converges from $\delta_{0}$ to $\delta_{1}^{*}$, and then, from $\delta_{1}^{*}$ to $\delta_{2}^{*},$ and finally, from $\delta_{\ce N}^{*}$ to $\delta^{*}_{\bf desired}$. This type of corrective actions reduces the need for prolonged load shedding and full state measurement. Also, this control method is unconventional where the operating point is relocated as desired, instead of being fixed as in the classical controls.

Mathematically, we consider the following problem:

• •

Emergency Control Design: Given a fault-cleared state $\delta_{0}$ and a desired stable EP $\delta^{*}_{\bf desired},$ determine the feasible dispatching of power injections $P_{k}$ to relocate the EPs so that the post-fault dynamics is driven from the fault-cleared state $\delta_{0}$ through the set of designed EPs to the desired EP $\delta^{*}_{\bf desired}$.

To solve this problem, we can design the first EP $\delta^{*}_{1}$ by minimizing $\|L^{{\dagger}}p\|_{\mathcal{E},\infty}$ over all possible power injections. The optimum power injection will result in an EP which is most far away from the stability margin $|\delta_{kj}|=\pi/2,$ and hence, the stability region of the first EP $\delta^{*}_{1}$ probably contains the given fault-cleared state. To design the sequence of EPs, in each step, we carry out the following tasks:

• •

Calculate the distance $R(\delta_{i-1}^{*})$ from $\delta_{i-1}^{*}$ to the boundary of the set $\mathcal{P},$ i.e., $R(\delta_{i-1}^{*})=\min_{\delta\in\partial\mathcal{P}}D(\delta,\delta_{i-1}^{*}).$ Noting that minimization of $D(\delta,\delta_{i-1}^{*})$ over the boundary of the set $\mathcal{P}$ is a convex problem with a quadratic objective function and linear constraints. Hence, we can quickly obtain $R(\delta_{i-1}^{*})$.

• •

Determine the set $\mathcal{B}(\delta_{i-1}^{*})$ and the set $\mathcal{A}(\delta_{i-1}^{*}).$

• •

The next EP $\delta_{i}^{*}$ will be chosen as the intersection of the boundary of the set $\mathcal{A}(\delta_{i-1}^{*})$ and the line segment connecting $\delta_{i-1}^{*}$ and $\delta_{\bf desired}^{*}.$

• •

The power injections $P_{k}$ that we have to redispatch will be determined by $P_{k}^{(i)}=\sum_{j\in\mathcal{N}_{k}}a_{kj}\sin\delta_{i_{kj}}^{*}$ for all $k.$ This power dispatch will place the new EP at $\delta^{*}_{i}$ which is in the inverse stability region of the previous EP $\delta_{i-1}^{*}$. Therefore, the controlled post-fault dynamics will converge from $\delta_{i-1}^{*}$ to $\delta_{i}^{*}$.

This procedure strictly reduces the distance from EP to $\delta_{\bf desired}^{*}$ (it can be proved that there exists a constant $d>0$ so that such distance reduces at least $d$ in each step). Hence, after some steps, the EP $\delta_{\ce N}^{*}$ will be sufficiently near the desired EP $\delta_{\bf desired}^{*}$ so that the convergence of the system state to the desired EP $\delta_{\bf desired}^{*}$ will be guaranteed.

To illustrate that this control works well in stabilizing some possibly unstable fault-cleared state $\delta_{0}$, we consider the 3-machine 9-bus system with 3 generator buses and 6 frequency-dependent load buses as in Fig. 3. The susceptances of the transmission lines are as follows: $B_{14}=17.3611p.u.,B_{27}=16.0000p.u.,B_{39}=17.0648p.u.,B_{45}=11.7647p.u.,B_{57}=6.2112p.u.,B_{64}=10.8696p.u.,B_{78}=13.8889p.u.,B_{89}=9.9206p.u.,B_{96}=5.8824p.u.$ The parameters for generators are: $m_{1}=0.1254,m_{2}=0.034,m_{3}=0.016,d_{1}=0.0627,d_{2}=0.017,d_{3}=0.008.$ For simplicity, we take $d_{k}=0.05,k=4\dots,9.$ Assume that the fault trips the line between buses $5$ and $7,$ and make the power injections to fluctuate. When the fault is cleared this line is re-closed. We also assume the fluctuation of the generation (probably due to renewables) and load so that the voltages $V_{k}$ and power injections $P_{k}$ of the post-fault dynamics are given in Tab. I. The stable EP $\delta_{\bf desired}^{*}$ is calculated as $[-0.1629\;0.4416\;0.3623\;-0.3563\;-0.3608\;-0.3651\;0.1680\;0.1362\;0.1371]^{\top}.$ However, the fault-cleared state, with angles $[0.025\;-0.023\;0.041\;0.012\;-2.917\;-0.004\;0.907\;0.021\;0.023]^{\top}$ and generators angular velocity $[-0.016\;-0.021\;0.014]^{\top},$ is outside the set $\mathcal{P}.$ It can be seen from Fig. 4 that the uncontrolled post-fault dynamics is not stable since $|\delta_{45}|$ and $|\delta_{57}|$ quickly evolve from initial values to $2\pi,$ which will activate the protective devices to trip the lines.

Using CVX software [13] to minimize $\|L^{{\dagger}}p\|_{\mathcal{E},\infty},$ we obtain the new power injections at buses 1-6 as follows: $P_{1}=0.5890,P_{2}=0.5930,P_{3}=0.5989,P_{4}=-0.0333,P_{5}=-0.0617,$ and $P_{6}=-0.0165.$ Accordingly, the minimum value of $\|L^{{\dagger}}p\|_{\mathcal{E},\infty}=0.0350<\emph{\emph{sin}}(\pi/89).$ Hence, the first EP obtained from equation (2) will be in the set defined by the inequalities $|\delta_{kj}|\leq\pi/89,\forall\{k,j\}\in\mathcal{E},$ and can be approximated by $\delta^{*}_{1}\approx L^{{\dagger}}p=[0.0581\;0.0042\;0.0070\;0.0271\;0.0042\;0.0070\;-0.0308\;-0.0486\;-0.0281]^{\top}$. The simulation results confirm that the post-fault dynamics is made stable by applying the optimum power injection control, as showed in Fig. 5. Using the above procedure, after one step, we can find that $\delta_{2}^{*}=0.9259\delta^{*}_{\bf desired}+0.0741\delta_{1}^{*}$ is the intersection of the set $\mathcal{A}(\delta_{1}^{*})$ and the line segment connecting $\delta_{1}^{*}$ and $\delta^{*}_{\bf desired}.$ This EP is inside the inverse stability region of $\delta_{1}^{*},$ and hence the system state will converge from $\delta_{1}^{*}$ to $\delta_{2}^{*}$ when we do the power dispatching $P_{k}^{(2)}$ corresponding to $\delta_{2}^{*}$. On the other hand, $\delta_{2}^{*}$ is very near the desired EP $\delta^{*}_{\bf desired}$ and it is easy to check that $\delta^{*}_{\bf desired}$ is in the inverse stability region of $\delta_{2}^{*},$ and thus the system state will converge from $\delta_{2}^{*}$ to the desired EP $\delta^{*}_{\bf desired}$. Such convergence of the controlled post-fault dynamics is confirmed in Figs. 6-7.

V Conclusions

Electric power grids possess rich dynamical behaviours, e.g., nonlinear interaction, prohibition of global stability, and exhibition of significant uncertainties, that challenge the maintenance of their reliable operation and pose interesting questions to control and power communities. This letter characterized a surprising property termed as “inverse stability”, which was rarely investigated and poorly understood (though some related inverse problems were addressed in [14]). This new notion could change the way we think about the stability assessment problem. Instead of estimating the set of initial states leading to a given operating condition, we characterized the set of operating conditions that a power grid converges to from a given initial state under changes in power injections and lines. In addition, we briefly described three applications of the inverse stability certificate: (i) assessing the stability of renewable power systems, (ii) solving the stability-constrained OPF problem, and (iii) designing power dispatching remedial actions to recover the transient stability of power systems. Remarkably, we showed that robust stability due to the fluctuation of renewable generations can be effectively assessed, and that the stability constraint can be incorporated as a static constraint into the conventional OPF. We also illustrated a unconventional control method, in which we appropriately relocate the operating point to attract a given fault-cleared state that originally leads to an unstable dynamics.

VI Appendix: Proof of Theorem 1

For each EP $\delta^{*}\in\mathcal{A}(\delta_{0})$, let $M$ be the point on the boundary of the set $\mathcal{P}$ so that $E(M,\delta^{*})=E_{\min}(\delta^{*})=\min_{\delta\in\partial\mathcal{P}}E(\delta,\delta^{*}),$ as showed in Fig. 1. From (8), we have

[TABLE]

Note that $D(\delta_{0},M)\geq R(\delta_{0})$ as $R(\delta_{0})$ is the distance from $\delta_{0}$ to the boundary of the set $\mathcal{P}.$ This, together with (VI), leads to $E(M,\delta^{*})+E(\delta_{0},\delta^{*})\geq R(\delta_{0})/2.$ From (8) and (9), we have $E(\delta_{0},\delta^{*})\leq F(\delta_{0},\delta^{*})<R(\delta_{0})/4.$ Hence,

[TABLE]

Therefore, for any $\delta^{*}\in\mathcal{A}(\delta_{0}),$ we have $E(\delta_{0},\delta^{*})<E_{\min}(\delta^{*}).$ By applying the stability analysis in Section III-A, we conclude that the initial state $\delta_{0}$ must be inside the stability region of the EP $\delta^{*}$ and the system state will converge from the initial state $\delta_{0}$ to the EP $\delta^{*}.$

VII Acknowledgements

This work was supported by the MIT/Skoltech, Ministry of Education and Science of Russian Federation (Grant no. 14.615.21.0001.), and NSF under Contracts 1508666 and 1550015.