Dynamic stability of electric power grids: Tracking the interplay of the network structure, transmission losses and voltage dynamics
aa r X i v : . [ n li n . AO ] A ug Dynamic stability of electric power grids: Tracking the interplay of the network structure, transmission lossesand voltage dynamics
Leonardo Rydin Gorjão
1, 2, a) and Dirk Witthaut
1, 2, b)1)
Forschungszentrum Jülich, Institute for Energy and Climate Research -Systems Analysis and Technology Evaluation (IEK-STE), 52428 Jülich,Germany Institute for Theoretical Physics, University of Cologne, 50937 Köln,Germany
Dynamic stability is imperative for the operation of the electric power system. This article provides analyticalresults and effective stability criteria focusing on the interplay of network structures and the local dynamics ofelectric machines. The results are based on an extensive linear stability analysis of the third-order model forsynchronous machines, comprising the classical power-swing equations and the voltage dynamics. The articleexplicitly covers the impact of Ohmic losses in the transmission grid, which are often neglected in analyticalstudies. Necessary and sufficient stability conditions are uncovered, and the requirements of different routesto instability are analysed, awarding concrete mathematical criteria applicable to all scales of power grids,from transmission to distribution and microgrids.
The secure supply of electric power is reliant onthe stable, coordinated operation of thousands ofelectric machines connected via the electric powergrid. The machines run in perfect synchrony withfixed voltage magnitudes and adjustable relativephases. The ongoing introduction of renewablepower systems poses several challenges to the sta-bility of the system, as line loads and tempo-ral fluctuations increase gravely. This trend takesplace in both the transmission grid at high volt-ages, as well as in distribution grids and micro-grids at middle and low voltages. This article con-tributes to the understanding of dynamical sta-bility of electric power systems and provide a de-tailed analysis of the third-order model for syn-chronous generators, which includes the transientdynamics of voltage magnitudes. Special empha-sis is laid on the impact of Ohmic losses in thetransmission, which are often neglected in ana-lytical treatments of power system stability. Theanalytical results thus find applicability on all sizescales of power grids, from transmission to iso-lated microgrids, for openly tackling systems withlosses in an analytical mathematical manner. Fur-thermore, the results are independent of the net-work construction and entail explicit criteria forthe connectivity of the power grid and the physi-cal requirements needed to ensure stability in thepresence of resistive terms. a) Electronic mail: [email protected] b) Electronic mail: [email protected]
I. INTRODUCTION
An unwavering operation of the electric power systemsis vital to our daily life, to the smooth operation of theeconomy and industry. An improved understanding of theelectrical system’s internal structure and properties is es-pecially relevant at this time as further renewable ener-gies enter the electric power-grid system . The stabilityof power-grid dynamics, facing the current revolution inboth technological development and market orientationtowards renewable energies, demands a thorough investi-gation of the properties of both the transmission as wellas the distribution grid.Conventional power grids typically amount to heavymass synchronous generators operating at a fixed fre-quency. These power grids are of various scales, spanningcontinents to single islands. More recently, the concept ofmicrogrids has emerged : partially independent powergrids in smaller environments that are coupled to a mainpower grid. These power grids operate at lower voltagesthan conventional transmission grids, which are capableor producing their own power, and work partially inde-pendently to an overarching power grid. Microgrids em-bedded in overarching power grids are still ruled by acommon understanding of fixed nominal frequency, e.g. Hz in Europe, amongst many other stability criteria.Models for electric machines, transmission and distri-bution power grids exist and are routinely used in ex-tensive numerical simulations . Analytical treatmentsmostly invoke significant simplification to keep the prob-lem tractable. Although studies with complex models ex-ist, cf. Schiffer et al and Dörfler et al , results arescarce for extended networks including resistive terms,given the difficulty of tackling dissipative systems math-ematically. Transmission grids, working at very high volt-ages, are able to circumvent this problem by minimisinglosses via controlling the power transmission carefully,but distribution grids—along with microgrids—may sus-tain considerable losses whilst transmitting power .1o surmount the issue of non-inductive power trans-fer analytically, this article puts forth a set of mathe-matical stability criteria for power grids, based on thethird-order model for synchronous generators. The crite-ria can be employed for various scales of power grids—forboth transmission and distribution grids—evidentiatingthe limitations entailed by the existence of resistive termson the operability of power-grid systems. In particular,the article undertakes the task of intertwining results forgraph theory with the characteristics of the power-gridconstruction and their physical properties.The article is structured in the following manner: Sec-tion II introduces the basic dynamical model studied inthis paper. Presented is an analysis of the third-ordermodel, comprising transient voltage dynamics and con-sidering extended grids with complex topology and resis-tive losses. Section III tackles the linear stability analy-sis of the equations of motion, a formal reduction of theproblem to a matrix formulation and develops a mathe-matical apparatus to unveil sufficient and necessary crite-ria for stability in a general sense. Section IV then utilisesthe developed concepts to derive analytical stability con-ditions for lossy systems, presenting criteria for stabilitynot only for the power-angle and the voltage dynamics,but also for a mixed type of instabilities, also comprisinga direct link to graph-theoretical measures. The conclu-sions follow subsequently in Section V. II. MODELLING SCALE-INDEPENDENTNETWORK-BASED POWER GRIDSA. Third-order model for synchronousgenerators
The third-order model for synchronous machines, de-noted as well as a one- or q -axis model, describes thetransient dynamics of coupled synchronous machines .It embodies the power- or rotor-angle δ ( t ) , relative tothe power-grid reference frame, the angular frequency ω ( t ) = ˙ δ ( t ) , in a co-rotating reference frame, and thetransient voltage E q ( t ) , in the q -direction of a co-rotatingframe of reference of each machine in the system. It ex-cludes sub-transient effects, i.e., higher-order effects, andassumes that the transient voltage E d in the d -directionof the co-rotating frame vanishes.Sub-transient effects play a small role, especially in thecase of studying power grids in the vicinity of the steadystate. The truncation of the transient voltage E d in the d -axis is imposed out of necessity to have an analyticallytractable model. Still, the resulting dynamical system israther complex such that analytical results are scarce andmostly restricted to lossless power grids. Hence, the scopeof the analysis here is two-fold: To present the detailsof tackling rotor-angle and voltage stability, whilst notshunning away from complex network topologies, includ-ing Ohmic losses.The equations of motion for one generator are given by , ˙ δ = ω,M ˙ ω = − Dω + P m − P el ,T ˙ E = E f − E + ( X − X ′ ) I, (1)where henceforth E ≡ E q denotes solely the voltage alongthe q -axis, and the dot the differentiation with respect totime. P m denotes the effective mechanical input powerof the machine, E f the internal voltage or field flux, and P el denotes the electrical power out-flow. The parameters M and D are the inertia and damping of the mechanicalmotion and T the relaxation time of the transient voltagedynamics. The voltage dynamics further depend on thedifference of the static X and transient X ′ reactancesalong the d -axis, where X − X ′ > in general, and thecurrent I along d -axis.The active electrical power P el j exchanged with thepower grid, and the current I j at the j -th machine readrespectively P el j = N X ℓ =1 E j E ℓ [ B j,ℓ sin( δ j − δ ℓ )+ G j,ℓ cos( δ j − δ ℓ )] ,I j = N X ℓ =1 E ℓ [ B j,ℓ cos( δ j − δ ℓ ) − G j,ℓ sin( δ j − δ ℓ )] , where the E j and δ j are the transient voltage and therotor angle of the j -th machine, respectively. The param-eters G j,ℓ and B j,ℓ denote the real and imaginary partsof the nodal admittance matrix and embody the networkstructure. Generally, B j,ℓ > and G j,ℓ < for all j = ℓ .This article is especially concerned with the role of Ohmiclosses, which are described by the real parts of G j,ℓ . Allquantities are usually made dimensionless using appro-priate scaled units referred to as the ‘pu system’ or ‘perunit system’ .The equations of motion (1) thus take the form for N coupled synchronous machines as ˙ δ j = ω j ,M j ˙ ω j = P mj − D j ω j + N X ℓ =1 E j E ℓ [ B j,ℓ sin( δ ℓ − δ j )+ G j,ℓ cos( δ ℓ − δ j )] ,T j ˙ E j = E fj − E j +( X j − X ′ j ) N X ℓ =1 E ℓ [ B j,ℓ cos( δ j − δ ℓ ) − G j,ℓ sin( δ j − δ ℓ )] . (2)The generality of the equations allows for graph- ornetwork-nodal reductions to be applied to the system,in a conceptual manner. Although a synchronous ma-chine in a power grid cannot physically be ignored, thepassive elements, e.g. connecting nodes or passive resis-tors, are mere algebraic structures. A power grid of syn-chronous generators and Ohmic loads can be reduced2y eliminating the load nodes, or buses, using a Kronreduction . The reduced system consists of the gen-erator nodes only, linked via an effective network, whichis typically fully connected. While in the original powergrid − ( G j,ℓ + iB j,ℓ ) directly gives the admittance of a line ( j, ℓ ) , the parameters B j,ℓ , G j,ℓ , P mj , and I j now repre-sent effective values characterising the reduced network.Most analytical studies so far neglected, under reason-able assumptions, the line losses of the power-grid struc-ture. The terms ∼ G j,ℓ are assumed negligible in com-parison to the terms ∼ B j,ℓ . Such arguments are rea-sonable for the high-voltage transmission grid, but aremostly unfounded for distribution and microgrids, wherethe resistance and inductance of transmission lines arecomparable . In addition, losses become more consider-able in magnitude when the transmitted power is large.This manuscript puts forth a study of the system in fullform, not discarding the interplay of susceptance andconductance, i.e., fully integrating losses. B. Stable state of operation
The stationary operation of the voltages and power-angles of the machines comprising the power grid is thecornerstone of operability of power grids. Constant volt-ages and perfect phase-locking, i.e., a point in configu-ration space where all E j , ω j and δ j − δ ℓ are constantin time, is the desired state. The latter restriction re-quires that all machines rotate at the same frequency δ j ( t ) = Ω t + δ ◦ j for all j = 1 , . . . , N , leading to the con-ditions ˙ ω j = ˙ E j = 0 , ˙ δ j = Ω , ∀ j = 1 , . . . , N. (3)In dynamical system terms, this is a stable limit cycleof the system. From a physical perspective, all points onthe limit cycle are equivalent as they only differ by aglobal phase δ which is irrelevant for the operation of thepower grid. One can thus choose one of these points as arepresentative of the limit cycle and refer to it as an ‘equi-librium state’. The superscript · ◦ is used to denote thevalues of the rotor-phase angle, frequency, and voltage inthis equilibrium state. Likewise, perturbations along thelimit cycle do not affect the power grid operation and canthus be excluded from the stability analysis. This shallbe made explicit in Definition 1.For the third-order model (2) an equilibrium state ofthe power grid is given by the nonlinear algebraic equa- tions Ω = ω ◦ j , P mj − D j Ω + N X ℓ =1 E ◦ j E ◦ ℓ (cid:2) B j,ℓ sin( δ ◦ ℓ − δ ◦ j )+ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) (cid:3) , E fj − E ◦ j +( X j − X ′ j ) N X ℓ =1 E ◦ ℓ (cid:2) B j,ℓ cos( δ ◦ j − δ ◦ ℓ ) − G j,ℓ sin( δ ◦ j − δ ◦ ℓ ) (cid:3) , noting that many equilibria—stable and unstable—canexist in networks with sufficiently complex topology, al-though this does not preclude performing a linear stabil-ity analysis . III. LINEAR STABILITY ANALYSISA. Fundamental equations and linearisation
One of the cornerstones of dynamical systems study islinear stability analysis . The local stability propertiesof an equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) , i.e., stability with respectto small perturbations around an equilibrium point, canbe obtained by linearising the equations of motion of thesystem (2).To perform a linear stability analysis of (2), one intro-duces the perturbations ξ j , ν j and ǫ j , such that δ j ( t ) = δ ◦ j + ξ j ( t ) , ω j ( t ) = ω ◦ j + ν j ( t ) , E j ( t ) = E ◦ j + ǫ j ( t ) . The perturbations ξ j , ν j , and ǫ j , once inserted intothe equations of motion (2), either lead to a growthor decay of the system’s characteristics over time. Therotor-angle perturbation ξ j , the frequency perturbation ν j , and the voltage perturbation ǫ j can—individually orcollectively—decay to zero or grow indefinitely. The sys-tem, around the equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) , is either stableor unstable, correspondingly. This is also known as ‘ex-ponential stability’ or ‘small signal stability’.Applying the linearisation into (2) whilst simultane-ously gauging onto a rotating frame of reference, withrotation frequency Ω as in (3), and preserving only termslinear in ξ j , ν j , and ǫ j , yields ˙ ξ j = ν j ,M j ˙ ν j = − D j ν j − N X ℓ =1 (Λ j,ℓ +Γ j,ℓ ) ξ ℓ + N X ℓ =1 ( A ℓ,j + C j,ℓ ) ǫ ℓ ,T j ˙ ǫ j = − ǫ j + ( X j − X ′ j ) N X ℓ =1 ( H j,ℓ + K j,ℓ ) ǫ ℓ (4) + ( X j − X ′ j ) N X ℓ =1 ( A j,ℓ + F j,ℓ ) ξ ℓ , Λ , Γ , A , C , F , H , K ∈ R N × N (writ-ten in component form above) are given by Λ j,ℓ = (cid:26) − E ◦ j E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P k = j E ◦ j E ◦ k B j,k cos( δ ◦ k − δ ◦ j ) for j = ℓ Γ j,ℓ = (cid:26) − E ◦ j E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) for j = ℓ P k = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) for j = ℓA j,ℓ = (cid:26) − E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) for j = ℓ P k E ◦ k B j,k sin( δ ◦ k − δ ◦ j ) for j = ℓC j,ℓ = (cid:26) E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P k E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓF j,ℓ = (cid:26) − E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P k = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓH j,ℓ = B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) ,K j,ℓ = − G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) . (5)The diagonal matrices M , D , X , and T (all in R N × N )comprise the elements M j , D j , ( X j − X ′ j ) , and T j for j = 1 , . . . , N , respectively All these diagonal matricesare strictly positive definite. An explicit derivation canbe found in the Supplemental Material A.The linearised system (4) takes a compact matrix for-mulation, where the linearised terms are elegantly com-bined into the Jacobian matrix J ∈ R N × N , by definingthe vectors ξ = ( ξ , . . . , ξ N ) ⊤ , ν = ( ν , . . . , ν N ) ⊤ , and ǫ = ( ǫ , . . . , ǫ N ) ⊤ , each in R N , with the superscript · ⊤ denoting the transpose of a matrix or vector. The lin-earised equations can be written as dd t ξνǫ = J ξνǫ , with J = − M − ( Λ + Γ ) − M − D M − ( A ⊤ + C ) T − X ( A + F ) T − ( X ( H + K ) − ) . (6)The Jacobian J can be brought to a different form thatclearly portrays the interplay between the matrices com-prising the susceptance B j,ℓ and the conductance terms G j,ℓ of the power lines and generators/motors, J = M −
00 0 T − X × − Λ − D A ⊤ A H − X − + − Γ 0 CF K . (7)This decomposition is conspicuously designed to work outthe impact of Ohmic losses. The left matrix in the brack-ets includes all terms ∼ B j,ℓ and embodies several sym-metries, which are further discussed in Section III D. The right matrix composed of the block matrices Γ , C , F , K embodies all the matrices associated with resistive losses ∼ G j,ℓ . The cleavage into two parts will prove usefulhence onward. B. Definition of Stability
An equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearly stable if pertur-bations in the linearised system (4) decay exponentially.In general, this is the case if and only if all eigenvaluesof the Jacobian matrix J have a negative real part .In the present case one has to take into account thatthe dynamical system incorporates a fundamental sym-metry. A shift of all nodal phase angles by a constantvalue Ψ( α ) : δ δ + α, with α ∈ S , S N → S N , does not have any physical effects: all flows, currents andstability properties remain unaffected. A geometric in-terpretation of this symmetry is obtained by viewing thedesired operation of the power grid as a limit cycle. As allpoints along the cycle are equivalent for power grid oper-ation, one can take an arbitrary point as a representativeof the limit cycle and refer to it as ‘the equilibrium’.As a consequence of this symmetry, any perturbationcorresponding to a global phase shift or a shift alongthe limit cycle, respectively, should be excluded from thestability analysis. This allows reducing the analysis tothe perpendicular subspaces of this symmetry, which aredefined as D (3) ⊥ = (cid:8) ( ξ , ν , ǫ ) ∈ S N × R N | ⊤ ξ = 0 (cid:9) , D (2) ⊥ = (cid:8) ( ξ , ǫ ) ∈ S N × R N | ⊤ ξ = 0 (cid:9) , D (1) ⊥ = (cid:8) ξ ∈ S N | ⊤ ξ = 0 (cid:9) . These subspaces are always one dimension smaller thanthe over-branching space. The subscript D ( · ) ⊥ refers to theorthogonality devised here, i.e., these spaces are orthog-onal to the stable limit-cycle manifold.Having defined the spaces of operation, one turns tothe Jacobian matrix (7) to unravel the definition of linearstability. Consider the eigenvalues µ , µ , . . . , µ N ∈ C N of the Jacobian defined via J ξνǫ = µ ξνǫ . There is always is one vanishing eigenvalue µ = 0 cor-responding to the global shift of all nodal phases, as dis-cussed above. One excludes this mode from the definitionof stability and orders the remaining eigenvalues accord-ing to their real parts, without loss of generality, µ = 0 , ℜ ( µ ) ≤ ℜ ( µ ) ≤ · · · ℜ ( µ N ) . ). Definition 1.
The equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearlystable if ℜ ( µ n ) < , for all eigenvalues n = 2 , . . . , N of the Jacobian matrix J defined in (6) . The above symmetry manifests itself equivalently inthe arising Laplacian structure composed by the matrix Λ + Γ , which satisfies N X ℓ =1 Λ j,ℓ + Γ j,ℓ = 0 , ∀ j = 1 , . . . , N. These matrices, in composite, are the Laplacian matrix,sometimes also referred to as Kirchhoff matrices, on theunderlying graph structure . They are ontologicallysingular matrices with at least one zero eigenvalue, asso-ciated with the eigenvector (1 , , . . . , ⊤ . This validatesthe aforementioned claim in Definition 1.Note that the matrix Λ is symmetric, thus correspond-ing to the Laplacian of an ordinary weighted graph. Incontrast, Γ is skew-symmetric such that it must rather beseen as the Laplacian of a signed directed graph. Whilethe properties of ordinary Laplacians are well known ,its directed counterparts constitute an active field of re-search, see e.g. . Amongst their several mathematicalproperties, of special relevance here are the eigenvaluesand associated eigenvectors of a graph’s or network’sLaplacian matrix. Lemma 1.
Let L ∈ R be an arbitrary Laplacian matrix,i.e., a real matrix that satisfies N X ℓ =1 L j,ℓ = 0 ∀ j = 1 , . . . , N, then all eigenvalues λ j are real and can be ordered as λ = 0 , λ ≤ λ ≤ · · · ≤ λ N . The eigenvalue λ is referred to as the Fiedler value oralgebraic connectivity of a network and the associatedeigenvector v F is the Fiedler vector . For a proof of reality, see e.g. . The algebraic con-nectivity λ is a measure of the connectivity of a graph,embodying its topological structure and the connected-ness of a graph. An undirected graph with positive edgeweights has an algebraic connectivity λ greater than ifand only if the graph is connected. This proves to be par-ticularly relevant when studying possible failure in powergrids that can entail the cleavage of the power grid in twoor several disconnected sub-graphs. C. Stability and algebraic decomposition
Similar to Definition 1, the stability of the linearisedsystem (4) can be evaluated by the positiveness or nega-tiveness of the Jacobian matrix J in (6). Having this in mind, one can appropriately construct the system’s Jaco-bian under a hermitian (symmetric) and non-hermitian(non-symmetric) description. Henceforth, the denom-ination hermitian, non-hermitian, and skew-hermitianwill additionally refer to symmetric, non-symmetric, andskew-symmetric, respectively. The generalisation appliesstraightforwardly to complex-valued matrices. Definition 2.
A complex matrix Z ∈ C N × N is positivedefinite on a subspace D ⊂ C N if ℜ ( y † Zy ) > , ∀ y ∈ D , y = 0 , and simply positive definite if D = C N . Conversely fornegative definiteness . Lemma 2.
Consider a complex matrix Z ∈ C N × N . Thefollowing statements are equivalent.1. Z is positive definite2. The hermitian part Z † := 12 ( Z + Z † ) (8) is positive definite.3. All eigenvalues of Z have positive real part.4. All eigenvalues of Z † have positive real part.Conversely for negative definiteness.Proof. (1) ⇔ (2) : A non-hermitian matrix Z can be de-composed into a hermitian and a skew-hermitian partvia Z = 12 ( Z + Z † ) + 12 ( Z − Z † ) , (9)where the first term in parenthesis is the hermitian part,as given in (8), and the second term ( Z − Z † ) is theskew-hermitian part. For every vector y ∈ C , ℜ (cid:2) y † Zy (cid:3) = 12 ℜ h y † ( Z + Z † ) y i + 12 ℜ h y † ( Z − Z † ) y i = ℜ (cid:2) y † Z † y (cid:3) . That is, the conditions for definiteness of the matrices Z and Z † are identical. (1) ⇔ (3) and (2) ⇔ (4) followfrom diagonalising the respective matrices. Negativenessfollows, mutatis mutandis , by the equivalent argument.The stability of the system (2) can now be studiedapplying Lemma 2 to the Jacobian. The hermitian part J † of the Jacobian J becomes the relevant part to bestudied, yet the system permits still another reduction.5 . Stability of a reduced Jacobian The reduced subspace of ( ξ , ǫ ) ∈ D (2) ⊥ suffices to studythe stability of the system, given the Jacobian matrix J .The structure of the equations of motion (1), in regardto the immediate equivalence of ˙ δ j = ω j , and thereforethe same for the linearisation (7), ˙ ξ = ν , immediatelysuggests that the variations of the rotor-angle δ j , givenby ξ j , are linearly dependent, one-to-one, on the variationof the frequency ω j , given by ν j . This implies that onecan simplify the linear stability analysis to the subspaceof angle and voltage variations. Lemma 3.
The equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearly stableif the matrix ΞΞ = (cid:18) − Λ − Γ A ⊤ + CA + F H + K − X − (cid:19) , (10) is negative definite on D (2) ⊥ .Proof. The eigenvalue problem for the Jacobian is J ξνǫ = µ ξνǫ , where the decomposition yields ν = µ ξ , (11a) − ( Λ + Γ ) ξ − Dν + ( A ⊤ + C ) ǫ = µ M ν , (11b) ( A + F ) ξ + (( H + K ) − X − ) ǫ = µ X − T ǫ . (11c)Substituting (11a) into (11b) and separating the termsmultiplied by the eigenvalue, one obtains − ( Λ + Γ ) ξ − µ Dξ + ( A ⊤ + C ) ǫ = µ M ξ , ( A + F ) ξ + (( H + K ) − X − ) ǫ = µ X − T ǫ , which can be written in matrix form as (cid:20)(cid:18) − Λ A ⊤ A H − X − (cid:19) + (cid:18) − Γ CF K (cid:19)(cid:21) (cid:18) ξǫ (cid:19) = (cid:18) µ M + µ D µ X − T (cid:19) (cid:18) ξǫ (cid:19) , where the matrices M , D , X , and T are all both diago-nal and positive definite. The matrix Ξ can also be seenexplicitly. Rewriting the above yields ( Ξ + Q ( µ )) (cid:18) ξǫ (cid:19) = , (12)where Q ( µ ) = (cid:18) − µ M − µ D − µ X − T (cid:19) . The proof follows by reductio ad absurdum : Take Ξ to be negative definite on D (2) ⊥ and assume µ ≥ is a possible eigenvalue. The matrix Q ( µ ) is then negativesemi-definite, thus so must be Ξ + Q ( µ ) , which is impos-sible under (12). A eigenvalue with µ ≥ does not existand the system is linearly stable.Conversely let Ξ be positive definite on D (2) ⊥ . Let β ≥ be an eigenvalue with corresponding eigenvector ( ξ , ǫ ) ⊤ such that Ξ (cid:18) ξǫ (cid:19) = β (cid:18) ξǫ (cid:19) . Each eigenvector must satisfy ǫ = , since equality wouldimply Λ ξ = , which is the case if either ξ = or ξ ∝ .The former is impossible since it results in the vector ( ξ , ǫ ) ⊤ = . The latter is excluded given ( ξ , ǫ ) ⊤ ∈ D (2) ⊥ .Thus, one can evaluate ( ξ , , ǫ ) ⊤ J ( ξ , , ǫ ) with the Jaco-bian given by (7), yielding ξ ǫ ⊤ M −
00 0 T − X × − Λ − D A ⊤ A H − X − + − Γ 0 CF K ξ ǫ = ξ ǫ ⊤ M −
00 0 T − X β ξǫ = ξ ǫ ⊤ β M − ξ β T − Xǫ = β ǫ ⊤ T − Xǫ ≥ , whereby the last inequality in given since X and T areboth positive definite. The Jacobian is not negative def-inite on D (3) ⊥ and the equilibrium is not linearly stable.This concludes the proof.The previous lemma shows that linear stability can beassessed on the basis of the reduced Jacobian matrix Ξ .A further simplification is possible by applying Lemma 2to Ξ , leading to the first main result of this work. Proposition 1.
An equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearlystable if the hermitian part of the reduced Jacobian matrix Ξ † , Ξ † = (cid:18) − Λ A ⊤ A H − X − (cid:19) + (cid:18) − Γ d NN (cid:19) , (13) is negative definite on D (2) ⊥ . The matrices Γ d = ( Γ + Γ † ) and N = ( C + F † ) are given by Γ d j,j = N X k = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) ,N j,j = N X k = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) . Ξ can befound in the Supplemental Material B.The form of the matrix Ξ † reveals the roles of the angleand voltage subspace and their interactions via the sub-matrices A and N . To further understand the role of an-gles and voltages and to derive explicit stability criteria,one can apply a decomposition of the stability criterionemploying the Schur, or Albert, complement . Proposition 2 (Sufficient and necessary stability condi-tions) . I. The equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearly stable if andonly if (a) the matrix Λ + Γ d is positive definite on D (1) ⊥ and (b) the matrix H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ is negative definite, where · + is theMoore–Penrose pseudoinverse.II. The equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) is linearly stable if andonly if (a) the matrix H − X − is negative defi-nite and (b) the matrix Λ + Γ d + ( A + N ) ⊤ ( H − X − ) − ( A + N ) is positive definite on D (1) ⊥ .Proof. Presented is solely the proof for criterion I. Anequivalent procedure works for criterion II and can befound in the Supplemental Material C. The reduced Ja-cobian matrix Ξ † can be decomposed as Ξ † = U ⊤ SU , (14)with S ∈ R N a diagonal block matrix S = (cid:18) − Λ − Γ d H − X − +( A + N )( Λ + Γ d ) + ( A + N ) ⊤ (cid:19) , and a transformation matrix U = (cid:18) − ( Λ + Γ d ) + ( A + N ) ⊤ (cid:19) . An explicit matrix formulation and derivation of thisdecomposition can be found in the Supplemental Ma-terial C.Notice that U is of full rank and maps the vector ( , ) ⊤ onto itself. Hence, U also maps the relevant sub-space D (2) ⊥ onto itself. Assume S is positive definite on D (2) ⊥ , then for every x ∈ D (2) ⊥ , x = , yielding x ⊤ Ξ † x = ( U x ) ⊤ S ( U x ) > . Similarly, assume that Ξ † is positive definite on D (2) ⊥ .Thus for every y ∈ D (2) ⊥ , y = one has y ⊤ Sy = ( U − y ) ⊤ Ξ † ( U − x ) > . The transformation U does not affect the definiteness: Ξ † is positive definite on D (2) ⊥ if S is positive definite on D (2) ⊥ and vice versa. Proposition 1 implies that the equilibrium is linearlystable if Ξ † , or equivalently S , is positive definite on D (2) ⊥ .Since S is block diagonal, definiteness of the entire matrixis equivalent to the definiteness of both blocks and theproposition follows.Given the whole mathematical apparatus developedabove, a short recapitulation is helpful. The third-ordermodel for synchronous machines was introduced and thelinearised dynamics around a given equilibrium were cal-culated and cast into a compact vectorial form. Subse-quently, the symmetries of the problem were discussedand several simplifications were made to show that lin-ear stability can be assessed on the basis of the reducedJacobian Ξ † . This analysis also illustrates the impact oflosses for linear stability, which enter the reduced Jaco-bian only via the two diagonal matrices Γ d and N . Fi-nally, the Schur complements allowed a further decompo-sition of the stability criteria, which elucidates the rolesof voltages and phase angles and allows for the derivationof explicit stability criteria. IV. EXPLICIT STABILITY CRITERIAA. Angle vs. Voltage Stability
The decomposition of the reduced Jacobian in Propo-sition 2 is of fundamental importance to this work, asit evinces the roles of the rotor angle and the voltagedynamics for the stability of the third-order model.Consider first the isolated power-angle dynamics, as-suming that the voltages ǫ remain fixed. Setting ǫ = 0 ,the linearised equations of motions read dd t (cid:18) ξν (cid:19) = (cid:18) − M − ( Λ + Γ ) − M − D (cid:19) (cid:18) ξν (cid:19) . Performing the same simplification as in the previous sec-tion, one finds that the isolated rotor angle dynamics arelinearly stable if and only if the matrix Λ + Γ d is positivedefinite on D (1) ⊥ .Similarly, consider the isolated voltage dynamics byassuming that the rotor angle remains fixed. Setting ν = ξ = 0 , the linearised equations of motions read dd t ǫ = T − X ( H + K − X − ) ǫ . Hence, one finds that the isolated voltage dynamics arelinearly stable if and only if the matrix H − X − is neg-ative definite, by once more studying the hermitian partof the matrix, as Lemma 2 indicates.In conclusion, one finds that the criteria I.(a) and II.(a)in Proposition 2 ensure the stability of the isolated rotorangle or voltage subsystem, respectively. Linear stabilityof the entire system is ensured only if, in addition, thecomplementary criteria I. (b) or II. (b) are satisfied.7o further elucidate the nature of the stability con-ditions, consider the full stability criterion I. in Proposi-tion 2. Assume that criteron I. (a) is satisfied, i.e., Λ + Γ d is positive definite on D (1) ⊥ , and the rotor angle subsys-tem is linearly stable. The complementary criterion I. (b)can then be written as H − X − ≻ ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ , (15)where ≻ denotes positive definiteness (equivalently, ≺ negative definiteness). This condition is far stricter thanthe condition of pure voltage stability, H − X − ≻ .Hence, stability of the two isolated subsystems is notsufficient, instead they must comprise a certain ‘securitymargin’ quantified by the right-hand side of (15) in orderto maintain linear stability.Making use of the angle-voltage decomposition, onecan derive explicit necessary and sufficient stability cri-teria. To this end, first consider the isolated subsystemsand subsequently the composite dynamics of the full sys-tem. Note that the lossless case has been discussed in ,thus here the focus is placed on the impact of Ohmiclosses. B. Voltage stability
Criterion II. (a) in Proposition 2 entails the stability ofthe isolated voltage subsystem. A violation implies theinstability of the voltage dynamics, and as a consequencealso the instability of the entire system, including therotor angle and frequency dynamics.Most remarkably, criterion II.(a) includes only the ma-trices H and X , which are also present in the losslesscase . Ohmic losses in the transmission lines thus affectvoltage stability only indirectly via the position of therespective equilibrium, in particular via the equilibriumrotor angles δ ◦ j , which enter the matrix H . Due to thesimilarity to the lossless case, this work refrains from adetailed analysis of voltage stability and only quotes tworesults from . Corollary 1.
If for all nodes j = 1 , . . . , N ( X j − X ′ j ) − > N X ℓ =1 B j,ℓ , then the matrix H − X − is negative definite. Corollary 2.
If for any subset of nodes
S ⊂{ , , . . . , N } , X j ∈S ( X j − X ′ j ) − ≤ X j,ℓ ∈S B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , then the matrix H − X − is not negative definite and theequilibrium is linearly unstable. C. Rotor angle stability
Criterion I. (a) in Proposition 2 entails the stability ofthe isolated rotor-angle subsystem. Briefly take the loss-less case into consideration, for which rotor-angle stabil-ity is determined by the matrix Λ . If for all connections ( j, ℓ ) in a power grid one has | δ j − δ ℓ | ≤ π π ≤ π , then Λ is a proper Laplacian matrix of a weighted undi-rected graph, which is well known to be positive definiteon D (1) ⊥ . If the condition is not satisfied for a line, the ma-trix Λ rather describes a signed graph, for which positivedefiniteness is more involved. Sufficient and necessary cri-teria have been obtained in .One can generalise the above condition to power gridswith Ohmic losses. Recall that the real parts of the eigen-values of the matrices Λ + Γ and Λ + Γ d are equal. Corollary 3.
If for all connections ( j, ℓ ) in a power grid,one has B j,ℓ cos (cid:0) δ ◦ ℓ − δ ◦ j (cid:1) + G j,ℓ sin (cid:0) δ ◦ ℓ − δ ◦ j (cid:1) > , (16) then the matrix Λ + Γ is positive definite on D (1) ⊥ and theisolated angle subsystem is linearly stable.Proof. The statement can be proved by applying Gerš-gorin’s circle theorem to Λ + Γ . Each eigenvalue λ j isbound to exist in a disk of radius R j = P ℓ = j | Λ j,ℓ + Γ j,ℓ | around the centre Λ j,j + Γ j,j such that | λ j − (Λ j,j + Γ j,j ) | ≤ X ℓ = j | Λ j,ℓ + Γ j,ℓ | . If condition (16) is satisfied, one can simplify this relationto | λ j − (Λ j,j + Γ j,j ) | ≤ X ℓ = j Λ j,ℓ + Γ j,ℓ = (Λ j,j + Γ j,j ) , which directly yields λ j ≥ . The only missing step is to show that the eigenvalue λ =0 . To this end, assume that x is an eigenvector of L = Λ + Γ of the eigenvalue λ = 0 and find its smallestelement, x ℓ = argmin j x j such that ( x ℓ − x j ) ≤ for allnodes j . Thus one has Lx ) ℓ = L ℓ,ℓ x ℓ + X j = ℓ L ℓ,j x j = X j = ℓ L ℓ,j ( x ℓ − x j ) .
8y the assumption (16), the entries L ℓ,j are positive forall nodes j connected to ℓ . Hence one must have x j = x ℓ for all nodes j connected to ℓ . One can now proceed to thenext-to-nearest neighbours of ℓ and then further throughthe power grid to show that x j = x ℓ ∀ j = 1 , . . . , N. Thus, there is a unique eigenvector with eigenvalue λ =0 , which is given by (1 , , . . . , ⊤ . D. Mixed instabilities
One now turns to the interplay of voltage and angle sta-bility, i.e., further investigating criteria I.(b) and II.(b) inProposition 2. Unless stated otherwise, consider an equi-librium such that the criteria I. (a) and II. (a) in Propo-sition 2 are satisfied. Hence, the isolated subsystems arestable, but the full system can still become unstable.To begin, consider the case where the voltage dynamicsare very stiff, i.e., the case where ( X j − X ′ j ) are small. Re-call that in the limit ( X j − X ′ j ) → the voltage dynamicsare trivially stable such that stability is determined solelyby the angular subsystem. One can extend this analysisto the case of small but non-zero ( X j − X ′ j ) and relatestability to the connectivity of the power grid. Recall thatstability of the isolated rotor-angle subsystem is ensuredif (cf. criterion I. (a) in Proposition 2) λ > , where λ is the Fiedler value or algebraic connectivity ofthe Laplacian Λ + Γ . Corollary 4.
A necessary condition for the stability ofan equilibrium point is given by λ > v ⊤ F h A ⊤ XA + 2 A ⊤ XN + N XN i v F (17) + O (( X j − X ′ j ) ) , where v F denotes the Fiedler vector of the Laplacian Λ + Γ for ( X j − X ′ j ) ≡ , as previously indicated in Lemma 1.Proof. The normalised Fiedler vector, at ( X j − X ′ j ) ≡ ,is denoted v F . The actual normalised Fiedler vector, fora particular non-zero value of the ( X j − X ′ j ) , is denoted v ′ F , such that v ′ F = v F + O (( X j − X ′ j ) ) . Take the expansion ( X − − H ) − = ∞ X ℓ =0 X ( XH ) ℓ , such that at lowest order one obtains ( X − − H ) − = X + O (( X j − X ′ j ) ) . Now, Proposition 2, criterion II. (b) can be reformulatedas: ∀ yy ⊤ ( Λ + Γ ) y > y ⊤ ( A + N ) ⊤ ( X − − H ) − ( A + N ) y . For a particular choice of y one obtains a necessary con-dition for stability. Taking y = v ′ F , the above results in λ > v ′ ⊤ F ( A + N ) ⊤ ( X − − H ) − ( A + N ) v ′ F , were applying the aforementioned expansion on the right-hand side, at leading order in ( X j − X ′ j ) , yields λ > v ⊤ F ( A + N ) ⊤ X ( A + N ) v F + O (( X j − X ′ j ) ) . Given now the symmetries of A , N , and X , one canexpand the result λ > v ⊤ F h A ⊤ XA + 2 A ⊤ XN + N XN i v F + O (( X j − X ′ j ) ) . This concludes the proof. This corollary entails a previousresult, vide .Note that each term of matrices on the right hand sideof (17) is symmetric and hence contributes positively,adding to the lower bound on the algebraic connectivity λ of the system. This implies that resistive networks always require a higher degree of connectivity to ensurestability. Claim 1.
A resistive power grid needs to ensure λ > X j ( X j − X ′ j ) v F j N X k = j E ◦ k G j,k , (18) in the limiting case of no power exchange.Proof. If there is a negligible power exchange in the powergrid, all rotor angles δ ◦ j , ∀ j are identical, such that cos( δ ◦ ℓ − δ ◦ j ) = 1 , sin( δ ◦ ℓ − δ ◦ j ) = 0 , for all connections ( j, ℓ ) . This results in A j,ℓ = 0 and Γ dj,j = 0 in (5), and Corollary 4 reads λ > v ⊤ F N XN v F , where all matrices are diagonal matrices. Writing theterms explicitly yields (18), entailing a lower bound to theconnectivity of a power grid with resistive elements.Before proceeding with the two final corollaries, notethat despite the cumbersome mathematical matricial no-tation employed here, one can still extract very usefulinformation—which can easily be computed numericallyif desired—by bestowing several matrix norms .9 emma 4. Let Z ∈ C N × N and W ∈ C N × N be twomatrices, and let k ·k n denote an n -induced matrix norm,one has k ZW k n ≤ k Z k n k W k n , i.e., all induced matrix norms are sub-multiplicative. Furthermore, recall that k · k denotes the ℓ -norm forvectors, also known as spectral norm or Euclidean norm . Corollary 5.
For a positive algebraic connectivity λ > , and all nodes j = 1 , . . . , N , ( X j − X ′ j ) − − N X ℓ =1 B j,ℓ > k ( A + N ) k k ( A + N ) ⊤ k λ , (19) where k · k is the induced ℓ -norm, then an equilibriumpoint is linearly stable.Proof. A positive algebraic connectivity λ > impliesthat both Λ + Γ and Λ + Γ d are positive definite on D (1) ⊥ ,and criterion I. (a) in Proposition 2 is satisfied.Consider now criterion I. (b) in Proposition 2. UsingGeršgorin’s circle theorem, as in the proof of Corollary 3,one finds that condition (19) imply ( X − − H ) − λ − k ( A + N ) k k ( A + N ) ⊤ k , is positive definite. Noting that λ − = k ( Λ + Γ ) + k = k ( Λ + Γ d ) + k , this implies that ∀ yy ⊤ ( X − − H ) y > k A + N k k ( Λ + Γ ) + k k ( A + N ) ⊤ k k y k ≥ k ( A + N )( Λ + Γ ) + ( A + N ) ⊤ k k y k ≥ y ⊤ ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ y . Hence, matrix H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ is negative definite and criterion I. (b) in Proposition 2 issatisfied. The equilibrium is linearly stable. Corollary 6.
If by criterion II. (a) in Proposition 2 thematrix H − X − is negative definite, and if the algebraicconnectivity λ satisfies λ > k ( A + N ) ⊤ ( H − X − ) − ( A + N ) k , (20) where k · k is the induced ℓ -norm, then the equilibriumpoint is linearly stable.Proof. Firstly, take by assumption that H − X − is neg-ative definite as given by criterion II.(a) in Proposition 2.The assumption (20) implies that ∀ y ∈ D (1) ⊥ y ⊤ ( Λ + Γ d ) y ≥ λ k y k > k ( A + N ) ⊤ ( H − X − ) − ( A + N ) k k y k ≥ y ⊤ ( A + N ) ⊤ ( H − X − ) − ( A + N ) y , again noticing that the algebrain connectivity λ isequvalently obtain for ( Λ + Γ d ) or ( Λ + Γ ) . Thus thematrix ( Λ + Γ d ) + ( A + N ) ⊤ ( H − X − ) − ( A + N ) isnegative definite in D (1) ⊥ . Criterion II.(b) in Proposition 2is therefore satisfied and the equilibrium is linearly sta-ble. V. CONCLUSIONS
The third-order model describes the dynamics of syn-chronous machines and takes into account both the rotor-angle and the voltage dynamics. Analytical results for thedynamics on the stability of coupled machines in powergrids with complex topologies are rare, in particular ifOhmic losses are taken into account. In this article, acomprehensive linear stability analysis was carried outand several explicit stability criteria were derivedThe first main result depicts the influence of resistiveterms of the system after the linear stability analysis.Remarkably, these terms enter into the reduced systemJacobian only via the two diagonal matrices Γ d and N ,as shown in Proposition 1. As a second main result, adecomposition of the Jacobian into the rotor-angle andthe voltage subsystems in Proposition 2 is derived. Thisdecomposition reveals clearly how the interplay of bothsubsystems can lead to novel forms of instability and thusrequires additional security margins.Based on this decomposition, several explicit stabilityconditions are uncovered, both for the isolated subsystemas well as for the full systems, including rotor-angle andvoltage dynamics. In particular, one can show that volt-age stability is not affected directly by resistive terms,thus implying that studies on voltage stability can bewithstood in the purely lossless case. Furthermore, Corol-lary 4 and subsequent Claim 1 entail a strict minimumconnectivity of the power-grid network solely by the pres-ence of resistive terms, i.e., a lower bound to possible dy-namics on the system given the presence of losses in thesystem.The analytical insights unveiled here—and in par-ticular the careful mathematical evaluation of lossysystems—can prove relevant to further understand powergrids of all spatial scales, and of general graph construc-tions. By mathematically tackling the presence of lossesin the system, the applicability of the results is now ex-tended from the transmission grid to the distributiongrid, where losses play a fundamentally bigger role. Theseresults can be applied to systems withstanding losses, in-cluding isolated power grids and microgrids. Moreover, itopens the door to further research on higher-order modelsfrom a mathematical point-of-view, and can henceforthbe applied more generally to other power-grid models. VI. SUPPLEMENTAL MATERIAL
Added to this article are three appendices in the Sup-plemental Material. In the Supplemental Material A afull derivation of the linearisation of (2) is made explic-itly. In the Supplemental Material B a full calculation ofthe hermitian Jacobian Ξ † from (13), in Proposition 1,is laid. In the Supplemental Material C the Schur com-plements leading to Proposition 2 is made explicit.10 CKNOWLEDGMENTS
L.R.G. thanks all the help and text revisions byMarieke Helmich. We gratefully acknowledge supportfrom the German Federal Ministry of Education and Re-search (grant no. 03EK3055B) and the Helmholtz Asso-ciation (via the joint initiative “Energy System 2050 –A Contribution of the Research Field Energy” and thegrant no. VH-NG-1025).
REFERENCES R. Sims et al., “Integration of renewable en-ergy into present and future energy systems,” in
IPCC Special Report on Renewable Energy Sources and Climate Change Mitigation ,edited by O. Edenhofer et al. (Cambridge Univer-sity Press, Cambridge, United Kingdom, 2011) pp.609–706. R. H. Lasseter and P. Piagi, “Microgrid: A conceptualsolution,” in
IEEE Power Electronics Specialists Con-ference , Vol. 6 (Citeseer, 2004) pp. 4285–4291. J. Rocabert, A. Luna, F. Blaabjerg, and P. Ro-driguez, “Control of power converters in ac microgrids,”IEEE transactions on power electronics , 4734–4749(2012). P. W. Sauer and M. A. Pai,
Power system dynamicsand stability (Prentice Hall, Upper, Saddle, River, N,J, 1998). J. Machowski, J. Bialek, and J. Bumby,
Power systemdynamics, stability and control (John Wiley & Sons,New York, 2008). P. C. Krause, O. Wasynczuk, S. D. Sudhoff, andS. Pekarek,
Analysis of electric machinery and drivesystems , Vol. 2 (Wiley Online Library, 2002). J. Schiffer, R. Ortega, A. Astolfi, J. Raisch,and T. Sezi, “Conditions for stability ofdroop-controlled inverter-based microgrids,”Automatica , 2457 – 2469 (2014). J. Schiffer, T. Seel, J. Raisch, and T. Sezi,“Voltage stability and reactive power shar-ing in inverter-based microgrids withconsensus-based distributed voltage control,”IEEE Transactions on Control Systems Technology , 96–109 (2016). F. Dörfler and F. Bullo, “Synchroniza-tion and transient stability in power net-works and non-uniform Kuramoto oscillators,”SIAM J. Control Optim. , 1616–1642 (2010). F. Dörfler, M. Chertkov, and F. Bullo, “Synchroniza-tion in complex oscillator networks and smart grids,”Proc. Natl. Acad. Sci. USA , 2005–2010 (2013). K. Schmietendorf, J. Peinke, R. Friedrich, andO. Kamps, “Self-organized synchronization and volt-age stability in networks of synchronous machines,”Eur. Phys. J. Special Topics , 2577 (2013). J. Ma, Y. Sun, X. Yuan, J. Kurths, andM. Zhan, “Dynamics and collapse in a power system model with voltage variation: The damping effect,”PloS ONE , e0165943 (2016). K. Schmietendorf, J. Peinke, and O. Kamps, “On thestability and quality of power grids subjected to in-termittent feed-in,” arXiv preprint arXiv:1611.08235(2016). G. Kron,
Tensor analysis of networks (J. Wiley & Sons,New York, 1939). F. Dörfler and F. Bullo, “Kron reduction ofgraphs with applications to electrical networks,”IEEE Trans. Circuits Syst. I, Reg. Papers , 150–163 (2013). A. J. Korsak, “On the question of uniqueness of sta-ble load-flow solutions,” IEEE Transactions on PowerApparatus and Systems , 1093–1100 (1972). D. Manik, M. Timme, and D. Witthaut, “Cy-cle flows and multistability in oscillatory networks,”Chaos: An Interdisciplinary Journal of Nonlinear Science , 083123 (2017). R. Delabays, T. Coletta, and P. Jacquod, “Multistabil-ity of phase-locking and topological winding numbersin locally coupled kuramoto models on single-loop net-works,” Journal of Mathematical Physics , 032701(2016). S. Jafarpour, E. Y. Huang, K. D. Smith, and F. Bullo,“Multistable synchronous power flows: From geom-etry to analysis and computation,” arXiv preprintarXiv:1901.11189 (2019). S. H. Strogatz,
Nonlinear Dynamics And Chaos: WithApplications To Physics, Biology, Chemistry, And En-gineering (Westview Press, Boulder, 2001). K. Sharafutdinov, L. Rydin Gorjão, M. Matthiae,T. Faulwasser, and D. Witthaut, “Rotor-angle versus voltage instability in the third-order model for synchronous generators,”Chaos: An Interdisciplinary Journal of Nonlinear Science , 033117 (2018). F. R. K. Chung and F. C. Graham,
Spectral GraphTheory , 92 (American Mathematical Soc., 1997). D. Babić, D. Klein, I. Lukovits, S. Nikolić, and N. Tri-najstić, “Resistance-distance matrix: A computationalalgorithm and its application,” International Journal ofQuantum Chemistry , 166–176 (2002). M. Fiedler, “Algebraic connectivity of graphs,”Czechoslovak Math. J. , 298 (1973). M. Fiedler, “A property of eigenvectors of nonnegativesymmetric matrices and its application to graph the-ory,” Czechoslovak Math. J. , 619–633 (1975). H. Chao and Y. Xudong, “Lower bound for the secondsmallest eigenvalue of directed rooted graph laplacian,”in
Proceedings of the 31st Chinese Control Conference (IEEE, 2012) pp. 5994–5997. M. P. Drazin and E. V. Haynsworth, “Criteria forthe reality of matrix eigenvalues,” MathematischeZeitschrift , 449–452 (1962). C. R. Johnson, “Positive Definite Matrices,”The American Mathematical Monthly , 259–264 (1970). F. Zhang,
The Schur complement and its applications ,Vol. 4 (Springer Science & Business Media, 2006). W. Chen, D. Wang, J. Liu, T. Basar, K. H.Johansson, and L. Qiu, “On semidefiniteness of11igned Laplacians with application to microgrids,”IFAC-PapersOnLine , 097–102 (2016). W. Chen, J. Liu, Y. Chen, S. Z. Khong,D. Wang, T. Basar, L. Qiu, and K. H. Johans-son, “Characterizing the positive semidefinitenessof signed Laplacians via effective resistances,”IEEE 55th Conf. Decision and Control , 985 – 990 (2016). Y. Song, D. J. Hill, and T. Liu, “Small-disturbance angle stability analysis of mi-crogrids: A graph theory viewpoint,”IEEE Conf. Control Applications (CCA) , 201 – 206 (2015). D. Zelazo and M. Bürger, “On the def-initeness of the weighted Laplacian andits connection to effective resistance,”IEEE 53rd Conf. Decision and Control , 2895–2900 (2014). S. Gerschgorin, “Über die Abgrenzung der Eigenwerteeiner Matrix,” Izv. Akad. Nauk. USSR Otd. Fiz.-Mat.Nauk , 749–754 (1931). 12 upplemental Materials Appendix A: Derivation of the matrix formulation
The full derivation of the matrices and matricial notation used to examine (2) is made explicit, step-by-step. ˙ δ j = ω j ,M j ˙ ω j = P mj − D j ω j + N X ℓ =1 E j E ℓ [ B j,ℓ sin( δ j − δ ℓ ) + G j,ℓ cos( δ j − δ ℓ )] ,T j ˙ E j = E fj − E j + ( X j − X ′ j ) N X ℓ =1 E ℓ [ B j,ℓ cos( δ j − δ ℓ ) − G j,ℓ sin( δ j − δ ℓ )] , (A1)with M j , P mj , B j,ℓ , G j,ℓ , T j , E fj , X j , and X ′ j parameters of the system. In order to linearise the system, or in other words,to study the effects of a perturbation, the system is assumed to have an equilibrium ( δ ◦ j , ω ◦ j , E ◦ j ) , ∀ j = 1 , , . . . , N .
1. Linearisation δ j ( t ) = δ ◦ j + ξ j ( t ) , ω j ( t ) = ω ◦ j + ν j ( t ) , E j ( t ) = E ◦ j + ǫ j ( t ) . (A2)In vector form, the perturbations are denoted by ( ξ , µ , ǫ ) . By inserting (A2) into (A1) ˙ ξ j = ω ◦ j + ν j ,M j ˙ ν j = P mj − D j ( ω ◦ j + ν j ) + ( E ◦ j + ǫ j ) N X ℓ =1 ( E ◦ ℓ + ǫ ℓ ) B j,ℓ sin( δ ◦ ℓ − δ ◦ j + ξ ℓ − ξ j )+( E ◦ j + ǫ j ) N X ℓ =1 ( E ◦ ℓ + ǫ ℓ ) G j,ℓ cos( δ ◦ ℓ − δ ◦ j + ξ ℓ − ξ j ) ,T j ˙ ǫ j = E fj − ( E ◦ j + ǫ j ) + (cid:0) X j − X ′ j (cid:1) N X ℓ =1 ( E ◦ ℓ + ǫ ℓ ) B j,ℓ cos( δ ◦ ℓ − δ ◦ j + ξ ℓ − ξ j ) − (cid:0) X j − X ′ j (cid:1) N X ℓ =1 ( E ◦ ℓ + ǫ ℓ ) G j,ℓ sin( δ ◦ ℓ − δ ◦ j + ξ ℓ − ξ j ) . (A3)Keeping only first-order terms in ( ξ , µ , ǫ ) , and using the following expansions of the sine and cosine functions sin( x + x ◦ ) ≈ sin( x ◦ ) + x cos( x ◦ ) + O ( x ) , cos( x + x ◦ ) ≈ cos( x ◦ ) − x sin( x ◦ ) + O ( x ) , ˙ ξ j = ν j ,M j ˙ ν j = − D j ν j + (cid:13) ǫ j N X ℓ =1 E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) + E ◦ j N X ℓ =1 ǫ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) + (cid:13) E ◦ j N X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) + (cid:13) ǫ j N X ℓ =1 E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) + E ◦ j N X ℓ =1 ǫ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) − (cid:13) E ◦ j n X ℓ =1 E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) ,T j ˙ ǫ j = − ǫ j + (cid:0) X j − X ′ j (cid:1) (cid:13) n X ℓ =1 ǫ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) − (cid:0) X j − X ′ j (cid:1) (cid:13) n X ℓ =1 E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) − (cid:0) X j − X ′ j (cid:1) (cid:13) n X ℓ =1 ǫ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) − (cid:0) X j − X ′ j (cid:1) (cid:13) n X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) . (A4)
2. Casting into matrix form a. The A matrix First derivation (cid:13) The A matrix can be obtained by correctly understanding the summations involved and applying some trigonometrictricks. Take the terms (cid:13) in (A4) ǫ j N X ℓ =1 E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) + E ◦ j N X ℓ =1 ǫ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) , terms (cid:13) . (A5)The summations can be written in a separate form as ǫ j N X ℓ =1 ℓ = j E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) + E ◦ j N X ℓ =1 ℓ = j ǫ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) , where noticing the first and second term always cancel when ℓ = j , thus the diagonal entries are given entirely by thefirst term here, and the off-diagonal terms are given each element of the second term. Casting this into the matrixform A , associated with a vector ǫ = ( ǫ , ǫ , . . . ) ⊤ . In this case notice that the second term’s indices are exchangedin comparison with the first term, hence one defines by symmetry A = P Nk =1 k =1 E ◦ k B ,k sin( δ ◦ k − δ ◦ ) − E ◦ B , sin( δ ◦ − δ ◦ ) . . . − E ◦ N B ,N sin( δ ◦ − δ ◦ N ) − E ◦ B , sin( δ ◦ − δ ◦ ) P Nk =1 k =2 E ◦ k B ,k sin( δ ◦ k − δ ◦ ) . . . − E ◦ N B ,N sin( δ ◦ − δ ◦ N ) ... ... . . . ... − E ◦ B ,N sin( δ ◦ N − δ ◦ ) . . . . . . P Nk =1 k = N E ◦ k B N,k sin( δ ◦ k − δ ◦ N ) , which rewrites generally as A j,ℓ = (cid:26) − E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) for j = ℓ P Nk E ◦ k B j,k sin( δ ◦ k − δ ◦ j ) for j = ℓ In this case, one has not the application of A onto ξ , but its transpose form, i.e., (A5) (term (cid:13) ) takes the form N X ℓ A j,ℓ ǫ ℓ , or A ⊤ ǫ . econd derivation (cid:13) Identically, A is obtained via the term (cid:13) in (A4), i.e., n X ℓ =1 E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) , terms (cid:13) . Initially one segregates the state variables ξ l and ξ j into two summations and noticing again the terms ℓ = j cancel,one obtains N X ℓ =1 ℓ = j ξ ℓ E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) − ξ j N X ℓ =1 ℓ = j E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j ) . The equivalent matrix A arises, as in (A6). By inverting the indices, applying to the vector ξ = ( ξ , ξ , . . . ) ⊤ andtaking notice of the preceding negative sign in(A4) (preceding the term (cid:13) ), one obtains N X ℓ A ℓ,j ξ ℓ , or Aξ . b. The Λ Laplacian matrix
The Λ Laplacian matrix arises from the term (cid:13) in (A4), E ◦ j N X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) , term (cid:13) . By separating the two terms ( ξ ℓ − ξ j ) into E ◦ j N X ℓ =1 ξ ℓ E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) − ξ j E ◦ j N X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , one notices, for a particular ξ j , the term in the first summation with ℓ = j cancels the equivalent term in the secondsummation, thus having the same mathematical structures as the terms composing the A and F matrices, i.e., theycan be cast into Λ , associated with the vector ξ Λ j,ℓ = (cid:26) − E ◦ j E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P Nk = j E ◦ j E ◦ k B j,k cos( δ ◦ k − δ ◦ j ) for j = ℓ Notice here the that the definition implies a negative sign before, thus one gets, in a matrix form − N X ℓ Λ ℓ,j ξ ℓ , or − Λ ξ . c. The C matrix Proceeding equivalently for the terms associated with the voltage expansion and the conductance G j,ℓ on the secondequation, i.e., terms (cid:13) in (A4), ǫ j n X ℓ =1 E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) + E ◦ j n X ℓ =1 ǫ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , terms (cid:13) , which likewise can take a matrix form C j,ℓ applying on ǫ C j,ℓ = (cid:26) E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P Nk E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓ which is a symmetric matrix, and positive definite for | δ ◦ k − δ ◦ j | < π/ . In the used notation, it takes the form N X ℓ C ℓ,j ǫ ℓ , or Cǫ . . The Γ matrix The term (cid:13) in (A4), associated with the conductance G j,ℓ terms, is cast into a matrix form E ◦ j N X ℓ =1 E ◦ ℓ B j,ℓ sin( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) , terms (cid:13) . Separating the two terms ( ξ ℓ − ξ j ) into E ◦ j N X ℓ =1 ξ ℓ E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) − ξ j E ◦ j N X ℓ =1 E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) , and obtain Γ associated with ξ as (noticing a sign change) Γ j,ℓ = (cid:26) − E ◦ j E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) for j = ℓ P Nk = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) for j = ℓ The final form is − N X ℓ Γ ℓ,j ξ ℓ , or − Γ ξ . e. The F matrix The term (cid:13) in (A4), associated with the voltage expansion and the conductance G j,ℓ , is embodied in the F matrix,derived from N X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j )( ξ ℓ − ξ j ) , term (cid:13) . By segregating the state variables into N X ℓ =1 ξ ℓ E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) − ξ j N X ℓ =1 E ◦ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , The first and second term always cancel when ℓ = j . The matrix form F j,ℓ associated with a vector ξ can be writtenas F j,ℓ = (cid:26) − E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) for j = ℓ P Nk = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓ where − N X ℓ F ℓ,j ξ ℓ , or − F ξ . f. The H matrix The term (cid:13) in (A4) takes the form N X ℓ =1 ǫ ℓ B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , term (cid:13) . and is written as H j,ℓ = B j,ℓ cos( δ ◦ ℓ − δ ◦ j ) , which associated with the vector ǫ yields N X ℓ H ℓ,j ǫ ℓ , or Hǫ . . The K matrix Lastly, the term (cid:13) in (A4) is left to tackle, N X ℓ =1 ǫ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) , term (cid:13) , which is cast into (notice the negative term) − N X ℓ K ℓ,j ǫ ℓ , or − Kǫ . Appendix B: Constructing the hermitian Jacobian Ξ † of the general Jacobian Ξ In order to ensure positiveness, equivalently negativeness, of a non-hermitian matrix Z ∈ C N × N , a decompositioninto an hermitian part Z † and a skew-hermitian part Z s is possible via Z = 12 ( Z + Z † ) + 12 ( Z − Z † )= Z † + Z s , (B1)given as well in (9). Positiveness, equivalently negativeness, of the matrix Z is ensured by the positiveness, equivalentlynegativeness, of the hermitian part Z † , as in Lemma 2.The non-hermitian Jacobian matrix J defined in (6), or in a separated fashion in (7), is cast into hermitian andskew-hermitian parts via (B1). The non-hermitian reduced Jacobian Ξ in (10) in Lemma 3 is permitted an equivalentseparation into hermitian and skew-hermitian parts by (B1).The henceforth referred to as Jacobian matrix Ξ , as in (10), takes the form Ξ = (cid:18) − Λ − Γ A ⊤ + CA + F H + K − X − (cid:19) = (cid:18) − Λ A ⊤ A H − X − (cid:19) + (cid:18) − Γ CF K (cid:19) . From here one can reconstruct the hermitian part of the Jacobian matrix Ξ as Ξ † = 12 "(cid:18) − Λ A ⊤ A H − X − (cid:19) + (cid:18) − Λ A ⊤ A H − X − (cid:19) † + 12 "(cid:18) − Γ CF K (cid:19) + (cid:18) − Γ CF K (cid:19) † , = (cid:18) − Λ A ⊤ A H − X − (cid:19) + 12 (cid:18) − Γ − Γ † C + F † F + C † K + K † (cid:19) . Each matrix in the second term is calculated as − Γ ℓ,j − Γ † ℓ,j = (cid:26) − E ◦ j E ◦ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) − E ◦ j E ◦ ℓ G j,ℓ sin( δ ◦ j − δ ◦ ℓ ) for j = ℓ P Nk = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) + P Nk = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) for j = ℓ = (cid:26) for j = ℓ P Nk = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) for j = ℓC ℓ,j + F † ℓ,j = (cid:26) E ◦ ℓ G j,ℓ cos( δ ◦ ℓ − δ ◦ j ) − E ◦ ℓ G j,ℓ cos( δ ◦ j − δ ◦ ℓ ) for j = ℓ P Nk E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) + P Nk = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓ = (cid:26) for j = ℓ P Nk = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) for j = ℓF ℓ,j + C † ℓ,j = (transpose term above) = ( C ℓ,j + F † ℓ,j ) † = C ℓ,j + F † ℓ,j ,K ℓ,j + K † ℓ,j = N X ℓ =1 ǫ ℓ G j,ℓ sin( δ ◦ ℓ − δ ◦ j ) + N X ℓ =1 ǫ ℓ G j,ℓ sin( δ ◦ j − δ ◦ ℓ ) = 0 . Γ d and N , as Γ d j,j = N X k = j E ◦ j E ◦ k G j,k sin( δ ◦ k − δ ◦ j ) , N j,j = N X k = j E ◦ k G j,k cos( δ ◦ k − δ ◦ j ) , (B2)noticing the half-factor in (B2). As such, the final hermitian formulation of the reduced Jacobian matrix Ξ † , i.e.,inserting (B2) into (B2) yields Ξ † = (cid:18) − Λ A ⊤ A H − X − (cid:19) + (cid:18) Γ d NN (cid:19) . Appendix C: Schur, or Albert, complements of the Jacobian matrix Ξ † To further elucidate the decomposition into Schur or Albert complements of the Jacobian matrix Ξ † evidenced inProposition 2, the explicit matrices are shown here in full glory.
1. Proposition 2, criterion I
The Schur complement takes the form Ξ † = U ⊤ S U as in (14), where here the subscript · is the first Schurcomplement, takes the form Ξ † = (cid:18) − Λ − Γ d ( A + N ) ⊤ A + N H − X − (cid:19) = U ⊤ (cid:18) − Λ − Γ d H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ (cid:19) S U , with U = (cid:18) − ( Λ + Γ d ) + ( A + N ) ⊤ (cid:19) . (C1)It is straightforward to show by insertion that Ξ † ! = U ⊤ S U (cid:20) − ( A + N )( Λ + Γ d ) + (cid:21) (cid:20) − Λ − Γ d H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ (cid:21) (cid:20) − ( Λ + Γ d ) + ( A + N ) ⊤ (cid:21) = (cid:20) − ( A + N )( Λ + Γ d ) + (cid:21) (cid:20) − Λ − Γ d ( A + N ) ⊤ H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ (cid:21) = (cid:20) − Λ − Γ d ( A + N ) ⊤ ( A + N ) H − X − + ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ − ( A + N )( Λ + Γ d ) + ( A + N ) ⊤ (cid:21) = (cid:20) − Λ − Γ d ( A + N ) ⊤ ( A + N ) H − X − (cid:21) .
2. Proposition 2, criterion II
To show now how to obtain criterion II in Proposition 2 take Ξ † = U S U ⊤ where now the subscript · indicatesthe second Schur complement. The formulation takes the form Ξ † = (cid:18) − Λ − Γ d ( A + N ) ⊤ A + N H − X − (cid:19) = U (cid:18) − Λ − Γ d − ( A + N )( H − X − ) − ( A + N ) ⊤ H − X − (cid:19) S U ⊤ , with U = (cid:18) ( A + N ) ⊤ ( H − X − ) − (cid:19) . (C2)18y insertion one shows that Ξ † ! = U S U ⊤ (cid:20) ( A + N ) ⊤ ( H − X − ) − (cid:21) (cid:20) − Λ − Γ d − ( A + N )( H − X − ) − ( A + N ) ⊤ H − X − (cid:21) (cid:20) ( H − X − ) − ( A + N ) (cid:21) = (cid:20) ( A + N ) ⊤ ( H − X − ) − (cid:21) (cid:20) − Λ − Γ d − ( A + N )( H − X − ) − ( A + N ) ⊤ ( A + N ) H − X − (cid:21) = (cid:20) − Λ − Γ d − ( A + N )( H − X − ) − ( A + N ) ⊤ + ( A + N )( H − X − ) − ( A + N ) ⊤ ( A + N ) ⊤ ( A + N ) H − X − (cid:21) = (cid:20) − Λ − Γ d ( A + N ) ⊤ ( A + N ) H − X − (cid:21) ..