Grid-forming inverters and power-system stability · Chapter 1

The classical picture: machines that swing

A grid holds its frequency because steel rotors store kinetic energy and release it without being asked. This chapter derives that picture in the smallest setting that contains it: one machine against a stiff grid. It ends with two numbers, a natural frequency and a damping ratio, and with the observation that every term that produced them needs a rotor.

1.1 The problem: why frequency is a shared variable

Switch on a kettle in Manchester and a turbine shaft in Yorkshire slows down. The statement is literal, not a figure of speech, and the mechanism is the subject of this chapter.

An alternating-current grid carries power at one frequency. In the classical grid this chapter describes, that power comes from synchronous generators. A synchronous generator is a rotating machine whose rotor carries a magnetic field, set up by a direct current in a field winding. Its stationary winding, the stator, produces an alternating voltage as that field sweeps past it. The voltage frequency is fixed by the rotor speed and by the number of magnetic poles on the rotor; Definition 1.2 gives the relation. In steady state every such machine on the grid turns at the one speed that produces the grid frequency. That is what synchronous means. The lock is not a control law. It is a physical restoring torque. If one rotor advances in angle against the rest of the system, the electrical power it exports rises. That loads the shaft and pulls the rotor back; §1.4 derives this as equation (1.8). If it falls behind, its export drops, which unloads the shaft and lets it catch up. The machines are coupled like masses on springs, and §1.5 makes the spring constant exact.

The network itself holds no usable energy store. The electric and magnetic fields of its lines and transformers exchange energy every half cycle. They hold no net reserve on the timescale of this chapter, tenths of a second to seconds. At every instant, generation equals demand plus losses. Suppose demand rises and no governor has yet responded. (The governor is the controller that adjusts the mechanical power of a turbine when its speed changes. It acts within seconds; §1.3 lists it among the omissions of the swing equation.) The extra energy must then come from somewhere in the next few hundred milliseconds. The one store of the required size is the kinetic energy of the spinning masses. A 200 MVA machine with the inertia constant of Example 1.1 holds 700 MJ, enough to supply its full rating for 3.5 s (Definition 1.3). Drawing on that store slows them down. Because the machines are locked together, they all slow together, and the system frequency falls. Frequency is therefore not a local measurement. It is a system-wide accounting variable that reports the running balance of power.

That observation drives the whole chapter. It means one scalar, the frequency, tells an operator whether the system is in balance. It means the rate at which frequency moves measures how much kinetic energy is spinning. And it means that protecting the system against a sudden loss of generation is a question about two things. The first is stored energy. The second is how far a rotor angle can swing before the restoring torque runs out.

Two distinct questions follow, and this chapter answers both.

  1. Small disturbances. After a load step of one per cent, does the angle return to its equilibrium, and how does it get there? This is a question about the eigenvalues of a linearised model. §1.5 answers it.
  2. Large disturbances. After a short circuit that removes the electrical load from a machine entirely for a fraction of a second, does the angle come back at all? Linearisation is useless here, because the angle moves through tens of degrees. §1.6 answers it by an energy argument.

The reason to build this picture carefully, in a book about inverters, is that all of it follows from one mechanical fact. A rotor is an integrator whose input is the difference between mechanical power in and electrical power out. Remove the rotor and every result in this chapter must be re-derived or abandoned. Chapter 2 shows which ones are abandoned. Chapter 3 shows which ones a control law can rebuild, and at what cost. So the two numbers computed in Example 1.1, and the critical clearing time computed in Example 1.2, are the yardstick for the rest of the book.

A table of every symbol used in this chapter, with the section that defines it, stands before the Exercises.

1.2 Per unit and base quantities

A transmission system contains transformers, and a transformer changes the numerical value of every voltage, current and impedance that passes through it. Writing the network equations in volts and amperes therefore adds turns ratios to every line. Those ratios carry no physical content. The per-unit system removes them.

Definition 1.1 (per unit and base quantities) R01

Definition — stated here

Choose two independent base quantities for a part of the network: a three-phase base apparent power Sbase in VA, and a line-to-line base voltage Vbase in V. Derive the other two:

Ibase = Sbase / (√3 Vbase)  and  Zbase = Vbase2 / Sbase.
(1.1)

The per-unit value of any quantity q is q/qbase, a dimensionless number written with the tag pu. Time is not normalised in this book; it stays in seconds. Angular frequency uses the base ωbase = ω0 = 2πf0, where f0 is the rated system frequency. Every example in this chapter uses f0 = 60 Hz, so ω0 = 2π × 60 = 376.9911 rad/s.

Why per unit removes the turns ratio

Take an ideal two-winding transformer with N1 turns on the primary and N2 turns on the secondary. An impedance Z2 in ohms connected on the secondary side, seen from the primary side, has the value

Z1 = (N1/N2)2 Z2.

Now pick one Sbase for the whole circuit. Choose the base voltages in the ratio of the turns: Vbase on the primary and (N2/N1)Vbase on the secondary. By (1.1) the base impedance on the secondary side is then (N2/N1)2 times the base impedance on the primary side. Divide each of the two impedances above by its own base:

Z1 / Zbase = (N1/N2)2 Z2 / Zbase = Z2 / [ (N2/N1)2 Zbase ].

The two per-unit impedances are equal. The turns ratio has cancelled. In per unit the ideal transformer becomes a short circuit, and a network of transformers, lines and machines becomes one connected circuit with no ratios in it.

Two further conveniences follow, and both are used without comment for the rest of the book. First, the √3 of a three-phase system disappears: a balanced three-phase power is Sbase times its per-unit value, with no factor. Second, machine parameters in per unit on the machine's own rating fall in narrow numerical ranges, whatever the machine size. A per-unit number therefore carries engineering meaning that an ohm does not. A transient reactance near 0.3 pu and an inertia constant between 2 s and 8 s are typical of large turbo-generators [S01, Kundur 1994, Power System Stability and Control, Ch. 3 — verify].

One warning. A per-unit value is meaningless without its base. A reactance of 0.2 pu on a 100 MVA base and a reactance of 0.2 pu on a 500 MVA base differ by a factor of five as impedances. Every per-unit number in this book states its base, or takes the base of the unit it belongs to.

1.3 The swing equation

This section derives the equation that governs every result in the chapter. The derivation starts at Newton's second law for rotation and ends with a two-state model in per unit.

The rotor angle and the speed deviation

Definition 1.2 (rotor angle and per-unit speed deviation) R02

Definition — stated here

Let ωm be the mechanical angular speed of the rotor in rad/s. A machine with p magnetic poles produces an electrical angular frequency

ω = (p/2) ωm.

Define the per-unit speed deviation as the single symbol

Δω = (ω − ω0) / ω0.
(1.2)

Define the rotor angle δ as the electrical angle of the rotor measured against a reference frame that rotates at exactly ω0. The rotor angle therefore moves only when the machine is off rated speed, and

dδ/dt = ω − ω0 = ω0 Δω.
(1.3)

A machine running exactly at rated speed holds δ constant at whatever value it had. The angle is a position, not a rate. That is what makes it the state variable of the system.

The inertia constant

Definition 1.3 (inertia constant and stored energy) R03

Definition — stated here

A rotor of moment of inertia J turning at rated mechanical speed stores kinetic energy. By Definition 1.2 the rated mechanical speed of a machine with p poles is 2ω0/p, so

Ekin = ½ J (2ω0/p)2.
(1.4)

Ekin is a constant of the machine, evaluated at rated speed. It is not the instantaneous kinetic energy, which varies with ωm.

The inertia constant is that energy normalised by the machine rating:

H = Ekin / Sbase,  so  Ekin = H Sbase.
(1.5)

Its unit is joule per volt-ampere, which is the second. That is not an accident of algebra; it has a direct reading. H is the time for which the machine could supply its own rated power from its stored kinetic energy alone, before the rotor stopped. Large turbo-generators have H between 2 s and 8 s [S01, Kundur 1994, Power System Stability and Control, Ch. 3 — verify]. Exercise 2 computes 0.444 s for a generator rotor alone, which is 22 % of the 2 s lower end. The turbine stages carry the remainder.

From Newton to the swing equation

Newton's second law for a rotating shaft states that the net torque equals the moment of inertia times the angular acceleration:

J dωm/dt = Tm − Te,
(1.6)

where Tm is the mechanical torque applied by the turbine and Te is the electrical torque opposing it. Three steps convert (1.6) into the working form.

Step 1 — multiply by speed to get power. Multiply both sides by ωm. The right-hand side becomes ωm(Tm − Te), which is mechanical power minus electrical power. Here the first approximation enters. Power equals torque times speed exactly. Take the torque base as Sbase divided by the rated mechanical speed 2ω0/p. Then per-unit torque equals per-unit power exactly at rated speed, and differs from it by the factor 1 + Δω off rated speed. In every case in this book the speed deviation is a few parts in a thousand, so we set Tm = Pm and Te = Pe in per unit. Assumption: |Δω| ≪ 1, so that ωm may be replaced by its rated value wherever it multiplies a torque.

Step 2 — normalise by the rating. Divide by Sbase. The coefficient on the left is Jωm/Sbase. Under the assumption of Step 1, replace the ωm in it by the rated value 2ω0/p. By (1.4), J = 2Ekin/(2ω0/p)2, so the coefficient becomes 2Ekin/[(2ω0/p) Sbase] = 2H/(2ω0/p) by (1.5). The left side is then 2H times the time derivative of ωm/(2ω0/p) = ω/ω0 = 1 + Δω, whose derivative is dΔω/dt by (1.2). The pole factor p/2 has cancelled between mechanical and electrical speed.

Step 3 — add the damping term. Real machines dissipate power in damper windings. These are short-circuited conductors on the rotor; current flows in them only when the rotor speed differs from the speed of the stator field. Real load also falls when frequency falls. The classical model represents both by one torque proportional to the speed deviation, with coefficient KD.

Model 1.4 (the swing equation, per unit) R04

Derived here from (1.6)

2H dΔω/dt = Pm − Pe − KD Δω,    dδ/dt = ω0 Δω.
(1.7)

Assumptions. (i) |Δω| ≪ 1, so per-unit torque equals per-unit power. (ii) H is constant; the rotor does not change shape. (iii) Damping is linear in the speed deviation. (iv) Pm is an external input; governor dynamics are not modelled.

Omissions, named. Turbine and governor response, which acts on a timescale of seconds and is therefore absent from the sub-second events of §1.6. Shaft torsional modes, which are reported above 5 Hz [S01, Kundur 1994, Power System Stability and Control, Ch. 15 — verify]. Any speed dependence of H.

Reading. The rotor is an integrator. Its input is the power imbalance Pm − Pe. Its output is a speed deviation, which (1.3) integrates a second time into an angle. Everything in this chapter follows from that double integration and from what Pe depends on.

Notation decision: KD, never D

The damping coefficient of (1.7) is written KD throughout this book. The bare symbol D is reserved for the state-space feedthrough matrix of Chapter 4, and Dload for load damping in MW/Hz. Virtual and equivalent damping in Chapter 3 are written KD,v and KD,eq, never Dv or Deq.

The convention that is excluded

A second convention for (1.7) is common. It writes the inertia coefficient as M = 2H/ω0 and lets the damping act on a speed deviation in rad/s rather than in per unit:

M dω/dt = Pm − Pe − K′D(ω − ω0).

This book does not use that form. The two forms are not interchangeable. The two speed deviations differ by the factor ω0 = 376.9911 s−1, so mixing them changes the computed damping ratio by exactly that factor. Exercise 3 carries the mixture through Example 1.1 and obtains a damping ratio (the symbol ζ, defined in §1.5) of 6.01 instead of 0.01594. A machine whose swing decays by the factor e in 7.0 s is then reported as overdamped, with no oscillation at all. From here on, every equation uses (1.7) with Δω in per unit.

1.4 The classical machine against an infinite bus

The swing equation is not yet closed, because Pe is not yet a function of the state. This section supplies that function under four assumptions, listed in Model 1.5. They are the fewest that make Pe a function of δ alone.

Two modelling choices do the work. The first is the infinite bus: a node whose voltage magnitude Vinf and frequency are fixed, whatever current is drawn from it. It is a Thevenin source of zero impedance. It idealises a machine whose own swing does not measurably move the system voltage or frequency. Chapter 4 replaces it with a network of finite size. The second is the classical machine model.

Model 1.5 (classical machine: constant emf behind transient reactance) R05

Derived here; its omissions are listed below and repaid in §1.6.1

Represent the machine as a voltage source of constant magnitude E′ behind a reactance X′d. Here E′ is the internal voltage that the rotor field induces in the stator winding. The prime marks it as the value fixed by the field flux in the first second after a disturbance. The reactance X′d is the transient reactance. It is the reactance the machine presents between that internal voltage and its terminals on the same timescale, before the field flux has changed. The qualifier direct-axis names the rotor axis that lies along the field winding. Under assumption (iii) below the two rotor axes have equal reactance, so the qualifier carries no content in this chapter. The symbol keeps it because it is standard. Typical values lie near 0.3 pu on the machine base [S01, Kundur 1994, Power System Stability and Control, Ch. 3 — verify]. The angle of that source is the rotor angle δ of Definition 1.2. Let X be the total reactance from the source to the infinite bus, that is X′d plus the transformer reactance Xt plus the reactance of the transmission path. That path is built from transmission lines, each of reactance XL. Then the active power delivered to the bus is

Pe = (E′ Vinf / X) sin δ = Pmax sin δ.
(1.8)

Assumptions. (i) The path between the source and the bus is purely reactive; resistance is zero. (ii) E′ is constant over the time of interest. (iii) The rotor is cylindrical (a round rotor), so the reactance the machine presents does not depend on the rotor position. (iv) Quantities are balanced and at fundamental frequency.

Omissions, named. Flux decay: the field flux linkage actually falls on a time constant reported in the range of seconds [S01, Kundur 1994, Power System Stability and Control, Ch. 3 — verify], so E′ holds constant only for about the first second. Excitation control: an automatic voltage regulator acts on E′ and changes both damping and synchronizing torque. Saliency: a rotor with projecting poles, as in hydro machines, presents a reactance that varies with rotor position. It adds a term in sin 2δ to (1.8). Stator resistance: including it adds a small cosine term and a small loss. Each omission is repaid by the classification of §1.6.1 and by the source cited there [S01, Kundur 1994, Power System Stability and Control, Ch. 3 and Ch. 12 — verify].

The derivation of (1.8) is four lines of phasor circuit theory. Take the infinite-bus voltage as the phase reference, Vinf∠0, and the internal voltage as E′∠δ = E′(cos δ + j sin δ). The two sources are joined by the reactance X, whose impedance is jX. The current I from the machine toward the bus is

I = (E′∠δ − Vinf∠0) / (jX).

The complex power S delivered to the bus is S = Vinf I*, where the asterisk denotes the complex conjugate. Substitute I and use 1/(−j) = j:

S = jVinf (E′ cos δ − jE′ sin δ − Vinf) / X.

Multiply out. The real part is Pe = (E′Vinf/X) sin δ, which is (1.8). The imaginary part is the reactive power (E′Vinf cos δ − Vinf2)/X, which this chapter does not use. In per unit S = VinfI* carries no factor 3, by the first convenience of §1.2.

Equation (1.8) is the whole electrical side of the classical picture, and it has three properties worth naming before they are used.

  1. It is bounded. No angle delivers more than Pmax = E′Vinf/X. A machine asked for more than that has no equilibrium at all.
  2. It is non-monotonic. Beyond δ = 90° the power falls as the angle rises. This is the source of every stability limit in the chapter.
  3. It depends on the network. Tripping a line raises X and lowers Pmax. Fig. 1.1 shows exactly that.

Take a mechanical power Pm < Pmax. The equilibrium angle δ0 is the angle at which electrical output equals mechanical input: Pmax sin δ0 = Pm. There are two solutions in (0, π), one below 90° and one above; §1.5 shows why only the first is usable. Fig. 1.1 marks δ0 for the pre-fault and post-fault networks of Example 1.2.

0 0.5 1.0 1.5 P (pu) 0 45 90 135 180 rotor angle δ (degrees) P_m = 0.8 pu pre-fault, P_max = 1.6923 post-fault, P_max = 1.2941 fault-on, P_e = 0 28.21° 38.18° 141.82°
Fig. 1.1 — The power-angle curve moves when the network changes. Equation (1.8) plotted for the three networks of Example 1.2, all with E′ = 1.1 pu and Vinf = 1.0 pu. Orange: before the fault, X = 0.65 pu, Pmax,pre = 1.6923 pu. Grey dashes along the axis: during the three-phase fault, Pe = 0. Teal: after line 2 is tripped, X = 0.85 pu, Pmax,post = 1.2941 pu. The white dashed line is Pm = 0.8 pu. Marked: the pre-fault equilibrium at 28.21°, the post-fault equilibrium at 38.18°, and the far intersection of the post-fault curve with Pm at 141.82°. §1.6 names that angle δmax; beyond it the post-fault curve lies below Pm again. The fault-on curve never meets Pm, so while the fault is on the rotor can only accelerate. Illustrates Model 1.5 (R05) and sets up Theorem 1.8 (R09).

1.5 Small-signal stability of one machine

Close the loop. Substitute (1.8) into the swing equation (1.7), and the system becomes two coupled first-order equations in the two states δ and Δω:

dδ/dt = ω0 Δω,     2H dΔω/dt = Pm − Pmax sin δ − KD Δω.

An equilibrium needs both derivatives to be zero, so Δω = 0 and Pmax sin δ0 = Pm: the equilibrium angle δ0 of §1.4. The rest of this section shows why only the solution below 90° is usable.

The synchronizing torque coefficient

Lemma 1.6 (synchronizing torque coefficient) R06

Proved here

Define the synchronizing torque coefficient as the slope of the power-angle curve at the equilibrium:

Ks = dPe/dδ at δ0 = (E′ Vinf / X) cos δ0.
(1.9)

Then Ks > 0 exactly when |δ0| < 90°.

Proof

Differentiate (1.8) with respect to δ. The derivative of Pmax sin δ is Pmax cos δ. Evaluate at δ0 to obtain (1.9). Now Pmax = E′Vinf/X is positive for a positive reactance and positive voltages. So the sign of Ks is the sign of cos δ0, which is positive exactly on (−90°, 90°). ∎

Reading. Ks is a spring constant, in units of power per radian. Push the rotor forward by one radian and the network pulls it back with an extra Ks per unit of load. Past 90° the spring reverses sign: pushing the rotor forward reduces its export, which accelerates it further. That is a runaway, not an oscillation.

The two numbers

Theorem 1.7 (small-signal single-machine infinite-bus result) R07

Proved here — the central result of Chapter 1

Linearise (1.7) with (1.8) about an equilibrium (δ0, 0). Write Δδ for the perturbation of the angle and Δω for the speed deviation. Then the swing mode has undamped natural frequency and damping ratio

ωn = √( Ks ω0 / (2H) ),    ζ = KD / (4 H ωn),
(1.10)

and the equilibrium is asymptotically stable exactly when Ks > 0 and KD > 0.

Proof

The only non-linear term is Pmax sin δ. Its first-order expansion about δ0 is Pmax sin δ0 + Pmax cos δ0 Δδ, that is Pm + Ks Δδ, using the equilibrium condition and (1.9). Substitute and cancel Pm:

2H dΔω/dt = −Ks Δδ − KD Δω,    dΔδ/dt = ω0 Δω.

Eliminate Δω. Differentiate the second equation once and substitute the first:

d2Δδ/dt2 + (KD/(2H)) dΔδ/dt + (Ks ω0/(2H)) Δδ = 0.

The characteristic polynomial is

λ2 + (KD/(2H)) λ + Ks ω0/(2H) = 0.
(1.11)

Compare with the standard second-order form λ2 + 2ζωnλ + ωn2 = 0. Matching the constant term gives ωn2 = Ksω0/(2H), which is the first half of (1.10). Matching the linear term gives 2ζωn = KD/(2H), so ζ = KD/(4Hωn), which is the second half.

For the stability claim, call the equilibrium asymptotically stable when every solution that starts near it returns to it as t grows without bound. For a linear second-order system that holds exactly when both roots of (1.11) have negative real part. By the Routh–Hurwitz condition for a monic quadratic, both roots have negative real part exactly when both of the remaining coefficients are positive. Here H > 0 and ω0 > 0 always, so the two conditions reduce to KD > 0 and Ks > 0. ∎

Boundary cases, stated because they matter later. If KD = 0 the roots are ±jωn: purely imaginary. The angle then oscillates forever at ωn, neither growing nor decaying. That is marginal stability, not asymptotic stability, and §1.6 uses exactly this undamped case. If Ks < 0 the constant term of (1.11) is negative. The product of the two roots is then negative, so one root is real and positive. The angle runs away.

Reading. The machine is a mass, a spring and a dashpot. To see the three constants, write the linearised pair with Δδ as the position. From dΔδ/dt = ω0Δω, the speed deviation is Δω = (1/ω0) dΔδ/dt. Substitute it into the first equation of the proof:

(2H/ω0) d2Δδ/dt2 + (KD/ω0) dΔδ/dt + Ks Δδ = 0.

The mass is 2H/ω0, the dashpot is KD/ω0, and the spring is Ks. Check against (1.10): (dashpot)2/(mass × spring) = KD2/(2Hω0Ks) = 4ζ2. Do not read the dashpot as KD alone. That pairs a per-unit damping coefficient with a position in radians, which is the convention mixture excluded in §1.3. Both numbers in (1.10) are inherited from hardware: H from the steel, Ks from the network reactance and the operating angle. Neither is chosen by a control engineer. Chapter 3 asks what happens when both become settings in a firmware file.

Example 1.1 — a 60 Hz machine against an infinite bus R08

Computed here; every digit is reproduced by the chapter's own arithmetic

Data. f0 = 60 Hz, so ω0 = 2π × 60 = 376.9911 rad/s. H = 3.5 s. KD = 2 pu. E′ = 1.1 pu. Vinf = 1.0 pu. Pm = 0.8 pu. Network: X′d + Xt = 0.45 pu, then two parallel transmission lines of XL = 0.40 pu each. All values on the machine base.

Step 1 — total reactance. Two equal reactances in parallel give half of one: 0.40/2 = 0.20 pu. So

X = 0.45 + 0.20 = 0.65 pu,    Pmax,pre = E′Vinf/X = 1.1 × 1.0 / 0.65 = 1.6923 pu.

Step 2 — equilibrium angle. Set Pe = Pm in (1.8):

sin δ0 = Pm X / (E′Vinf) = 0.8 × 0.65 / 1.1 = 0.52 / 1.1 = 0.472727,
δ0 = arcsin(0.472727) = 0.492383 rad = 28.211°,  cos δ0 = 0.881209.

Step 3 — synchronizing torque coefficient, from (1.9):

Ks = 1.692308 × 0.881209 = 1.4913 pu/rad.

Step 4 — natural frequency, from (1.10):

Ksω0/(2H) = 1.491276 × 376.9911 / 7 = 562.198 / 7 = 80.3140 s−2,
ωn = √80.3140 = 8.9618 rad/s = 8.9618 / (2π) = 1.4263 Hz.

Step 5 — damping ratio, from (1.10):

ζ = 2 / (4 × 3.5 × 8.9618) = 2 / 125.465 = 0.015941.

Step 6 — eigenvalues. The roots of (1.11) are −ζωn ± jωn√(1 − ζ2):

−0.015941 × 8.9618 = −0.142857 s−1, which equals −KD/(4H) = −2/14 exactly;
8.9618 × √(1 − 0.000254) = 8.9607 rad/s;
λ = −0.1429 ± j8.9607 s−1.

Interpretation. The oscillation period is 2π/8.9607 = 0.7012 s. The amplitude decays with time constant 1/0.142857 = 7.000 s, so the disturbance falls to 36.8 % of its initial amplitude after 7.000 s, which is 9.98 periods. The frequency, 1.43 Hz, sits inside the 0.7 Hz to 2 Hz band reported for local plant modes [S01, Kundur 1994, Power System Stability and Control, Ch. 12 — verify]. The damping ratio is 0.015941, that is 1.6 % of critical damping. After a load step this mode is still at 36.8 % of its initial amplitude ten swings later. Power system stabilisers exist to raise that damping. They add a supplementary signal to the excitation system whose purpose is to increase the damping of this mode [S01, Kundur 1994, Power System Stability and Control, Ch. 12 — verify].

0 −1 −2 −3 real part of λ (s⁻¹) +8.96 −8.96 0 imaginary part (rad/s) K_D = 0 2 10 30 |λ| = ω_n = 8.9618 s⁻¹
Fig. 1.2 — Damping moves the pair to the left; it does not move the natural frequency. Eigenvalues of Example 1.1 computed from (1.11) with H = 3.5 s and Ks = 1.4913 pu/rad held fixed, and KD swept over 0, 2, 10 and 30 pu. The pairs are 0 ± j8.9618, −0.1429 ± j8.9607, −0.7143 ± j8.9333 and −2.1429 ± j8.7018 s−1, with ζ = 0, 0.01594, 0.07970 and 0.23911. Every pair has the same magnitude, |λ| = ωn = 8.9618 s−1, so the locus is a circle of that radius. The white dashed arc is that circle. It looks nearly flat because one unit on the horizontal axis spans 10.2 times the length of one unit on the vertical axis. The teal pair at KD = 0 sits on the imaginary axis: marginally stable, not asymptotically stable, as Theorem 1.7 states. Illustrates Theorem 1.7 (R07) and Example 1.1 (R08).

1.6 Large-signal stability: the equal-area criterion

A short circuit on the transmission system is not a small disturbance. While the fault is on, the machine's electrical export collapses, the angle moves through tens of degrees, and the linearisation of §1.5 is void. The question changes from "how does it come back" to "does it come back at all". Protection engineers need the answer as a time: how long may a circuit breaker take to clear the fault before the machine is lost? That time is the critical clearing time tcr.

The tool is an energy argument. It comes out of one integration of the undamped swing equation.

Theorem 1.8 (equal-area criterion) R09

Proved here; the standard treatment is [S04, Anderson and Fouad 2003, Power System Control and Stability, 2nd ed., Ch. 2 — verify]

Take (1.7) with KD = 0 and constant Pm. Suppose the machine starts at the equilibrium δ0 at synchronous speed, so Δω = 0. Suppose a fault holds Pe = 0 from δ0 up to a clearing angle, and that after clearing Pe = Pmax,post sin δ. Write δmax = π − arcsin(Pm/Pmax,post) for the far intersection of the post-fault curve with Pm. Then the rotor first stops swinging at the angle δ at which

∫ from δ0 to δ of (Pm − Pe) dδ = 0.
(1.12)

Split that integral at the clearing angle. Call Aacc the integral of Pm − Pe over the angles where that difference is positive: the rotor gains kinetic energy there. Call Adec the integral of Pe − Pm over the angles after clearing where that difference is positive, up to δmax. The rotor loses kinetic energy there. When the clearing angle lies above the post-fault equilibrium arcsin(Pm/Pmax,post), Aacc is the rectangle of height Pm between δ0 and the clearing angle. Example 1.2 is such a case. The machine keeps synchronism exactly when Adec ≥ Aacc. The critical clearing angle δcr is the clearing angle at which the two areas are equal. It satisfies (1.13).

Proof

With KD = 0, (1.7) reads 2H dΔω/dt = Pm − Pe. Multiply both sides by dδ/dt = ω0Δω:

2H ω0 Δω dΔω/dt = (Pm − Pe) dδ/dt.

The left side is Hω0 d(Δω2)/dt. Integrate both sides in time from the start of the disturbance to the instant the angle reaches δ. On the right, change the variable of integration from t to δ. The change is valid because dδ/dt = ω0Δω > 0 until the rotor first stops, so δ(t) is strictly increasing on that interval:

Hω0 [Δω2] = ∫ from δ0 to δ of (Pm − Pe) dδ.

The left side is proportional to the kinetic energy the rotor has gained. The machine starts at rest, so Δω = 0 there, and it stops swinging when Δω returns to zero. The left side is therefore zero at both ends, which is (1.12).

Now let the fault be on from δ0 to the clearing angle, and be cleared from there onward. While the fault is on, Pe = 0, so the integrand is Pm > 0 and the integral accumulates a positive area: that is Aacc. After clearing, Pe follows the post-fault curve, which exceeds Pm between the post-fault equilibrium and δmax. There the integrand is negative and the integral subtracts: that is Adec. Condition (1.12) is therefore Aacc = Adec.

Beyond δmax, the far intersection of the post-fault curve with Pm, the integrand turns positive again and no further deceleration is available. If the rotor is still moving forward at δmax, it accelerates from there without bound. The critical case is therefore the one in which the rotor arrives at δmax with exactly zero speed. Write (1.12) out between δ0 and δmax, with Pe = 0 during the fault and Pe = Pmax,post sin δ after it, and solve for the clearing angle:

cos δcr = [ Pm(δmax − δ0) + Pmax,post cos δmax ] / Pmax,post. ∎
(1.13)

Assumptions. (i) No damping. (ii) Pm constant: the governor does not act inside the swing. (iii) E′ constant, as in Model 1.5. (iv) One machine against an infinite bus, so there is a single angle to integrate. The criterion does not extend to three or more machines without further argument. (v) Pe = 0 while the fault is on. A fault at a bus further from the machine leaves a reduced power-angle curve in place during the fault. Equation (1.12) still holds, because it holds for any Pe(δ). Equation (1.13) then gains the extra term −(1/Pmax,post) × (integral of the fault-on power from δ0 to δcr). Re-derive it from (1.12).

Why neglecting damping is the safe direction. Damping removes energy from the swing. A machine that survives with no damping survives with damping. The criterion therefore returns a conservative clearing time, and Exercise 4 checks that a numerical integration with no damping reproduces it.

Once δcr is known, the clearing time follows from the fault-on motion alone. While the fault is on, Pe = 0, so (1.7) with KD = 0 has a constant right-hand side and the angle follows a parabola in time. Integrating d2δ/dt2 = ω0Pm/(2H) twice from rest gives δ = δ0 + ω0Pmt2/(4H), which inverts to

tcr = √( 4H(δcr − δ0) / (ω0 Pm) ).
(1.14)

Example 1.2 — critical clearing time for the machine of Example 1.1 R10

Computed here

Event. Same machine, same network, same operating point: H = 3.5 s, Pm = 0.8 pu, δ0 = 0.492383 rad = 28.211°, f0 = 60 Hz. A three-phase fault occurs at the sending-end bus of line 2. It holds the voltage of that bus at zero. The path from E′ to that bus is a pure reactance, so no active power leaves the machine. Pe = 0 while the fault is on. The fault is cleared by tripping line 2, which leaves one line in service. Damping is neglected.

Step 1 — the three power-angle curves.

pre-fault: X = 0.45 + 0.20 = 0.65 pu, Pmax,pre = 1.1 / 0.65 = 1.6923 pu;
fault-on: Pe = 0;
post-fault: X = 0.45 + 0.40 = 0.85 pu, Pmax,post = 1.1 / 0.85 = 1.2941 pu.

Step 2 — the far intersection. The post-fault curve crosses Pm twice. The near crossing is arcsin(0.8 / 1.294118) = arcsin(0.618182) = 0.666427 rad = 38.183°. The far crossing is

δmax = π − 0.666427 = 2.475165 rad = 141.817°.

Step 3 — critical clearing angle, from (1.13):

Pm(δmax − δ0) = 0.8 × (2.475165 − 0.492383) = 0.8 × 1.982782 = 1.586226;
cos δmax = cos(2.475165) = −0.786035, so Pmax,post cos δmax = 1.294118 × (−0.786035) = −1.017222;
numerator = 1.586226 − 1.017222 = 0.569004;
cos δcr = 0.569004 / 1.294118 = 0.439685;
δcr = arccos(0.439685) = 1.115549 rad = 63.916°.

Step 4 — area check. With δcr − δ0 = 1.115549 − 0.492383 = 0.623166 rad,

Aacc = Pm(δcr − δ0) = 0.8 × 0.623166 = 0.498533 pu rad;
Adec = Pmax,post(cos δcr − cos δmax) − Pm(δmax − δcr)
= 1.294118 × (0.439685 + 0.786035) − 0.8 × 1.359616 = 1.586226 − 1.087693 = 0.498533 pu rad.

The two areas agree to six digits, as (1.12) requires.

Step 5 — critical clearing time, from (1.14):

4H(δcr − δ0) = 14 × 0.623166 = 8.724324;
ω0Pm = 376.9911 × 0.8 = 301.5929;
tcr = √(8.724324 / 301.5929) = √0.0289275 = 0.170081 s.

At 60 Hz one cycle is 1/60 = 0.016667 s, so 0.170081 × 60 = 10.205 cycles.

Interpretation. Transmission protection with a modern relay and circuit breaker clears a three-phase fault in about 4 to 6 cycles [S04, Anderson and Fouad 2003, Power System Control and Stability, 2nd ed., Ch. 2 — verify]. This machine therefore has 10.2 cycles against 6, a ratio of 1.7. A slower breaker, a heavier pre-fault loading, or a network with a larger X would consume that ratio. Note also that tcr scales as √H in (1.14). Half the inertia gives 0.7071 of the clearing time. That is 0.1203 s, or 7.2 cycles, against the 4 to 6 cycles a breaker needs.

0 0.5 1.0 1.5 P (pu) P_m = 0.8 δ₀ 28.21° δ_cr 63.92° δ_max 141.82° A_acc 0.4985 A_dec 0.4985 pre-fault post-fault
Fig. 1.3 — The critical clearing angle is the angle at which the two areas are equal. The same curves as Fig. 1.1, with the areas of Theorem 1.8 shaded. Orange rectangle: the accelerating area, accumulated from δ0 = 28.21° to δcr = 63.92° while Pe = 0. Its height is Pm = 0.8 pu and its width is 0.623166 rad, so Aacc = 0.498533 pu rad. Teal region: the decelerating area between the post-fault curve and Pm, from δcr to δmax = 141.82°, computed in Step 4 of Example 1.2 as Adec = 0.498533 pu rad. The two numbers agree to six digits, and that equality is what defines δcr. Clearing later than 63.92° makes the orange area larger than any teal area that remains, and the machine is lost. Illustrates Theorem 1.8 (R09) and Example 1.2 (R10).
0 50 100 150 200 δ (degrees) 0 0.3 0.6 0.9 1.2 time t (s) 0.15 s 0.19 s clear at 0.15 s — stable clear at 0.19 s — loses synchronism
Fig. 1.4 — The stability boundary sits between 0.15 s and 0.19 s. Numerical integration of the undamped swing equation (1.7) with (1.8), for the machine and the event of Example 1.2, by second-order Runge–Kutta at a step of 10 μs, starting from δ = 28.211° at rest. Teal, cleared at 0.15 s: the angle reaches a first peak of 106.7° and swings back, then oscillates undamped. Orange, cleared at 0.19 s: the angle passes δmax = 141.82° with speed still positive and runs away; the trace crosses 200° at 0.5662 s, and reaches 1206° at 1.2 s. A bisection on clearing time over the same integration, with failure declared when δ passes 200° inside the simulated window, depends on the window length. A 1.2 s window gives 0.1703 s; a 1.5 s window gives 0.1701 s; a 3 s window gives 0.170085 s. The last value sits 4 μs above the 0.170081 s of (1.14). That residual is the 10 μs step: an exact fault-on solution followed by a fourth-order integrator converges to 0.170081 s. The window matters because a rotor cleared just after tcr lingers near δmax before it runs away. Illustrates Theorem 1.8 (R09) and Example 1.2 (R10).

1.6.1 Where this chapter sits: the stability classification

The chapter has now produced a complete account of one machine. It is worth naming where that account sits in the wider subject, and then naming the single physical assumption it rests on.

Remark 1.1 (classification of power-system stability) R11

Cited, not proved here

The accepted classification divides power-system stability into three classes [S02, Kundur et al. 2004, "Definition and Classification of Power System Stability", IEEE Trans. Power Systems 19(3), 1387–1401, §II — verify]:

The 2021 revision of that classification adds two classes that the 2004 classification did not contain [S03, Hatziargyriou et al. 2021, "Definition and Classification of Power System Stability — Revisited and Extended", IEEE Trans. Power Systems 36(4), 3271–3281, §III — verify]:

That two classes had to be added in 2021 is itself the evidence that the classical picture is incomplete for a converter-dominated system.

1.6.2 The single assumption

Look back at what produced each result. One object appears in every line.

Every result of Chapter 1 traced back to the rotor.
ResultWhat it needsWhere that comes from
Swing equation (1.7), Model 1.4a mass that integrates power imbalancethe rotor's moment of inertia J
H, Definition 1.3stored kinetic energy at rated speedthe rotor's mass and radius
Ks, Lemma 1.6a voltage source of fixed magnitude behind a reactance, whose angle is a physical positionthe rotor's angular position, not a measurement
ωn and ζ, Theorem 1.7both of the abovethe rotor
tcr, Theorem 1.8the kinetic energy the rotor absorbs during the faultthe rotor

An inverter has no rotor. It has a direct-current capacitor whose stored energy is smaller by a factor that Chapter 3 computes exactly. It also has an angle that conventional control takes from a measurement, not from a state. None of the five rows above therefore transfers automatically. Chapter 2 shows what happens to an inverter that takes its angle from a measurement as the machines around it retire. It derives a limit that is the structural twin of Lemma 1.6. That limit is a term of the form (voltage) × cos(angle). It falls to zero, and it takes the loop's stability with it. Chapter 3 asks which of the five rows a control law can rebuild. Chapter 4 asks how much of it a system needs.

What you can now do. (1) Write the swing equation of a machine in per unit from its moment of inertia and rating, Model 1.4. (2) Compute the synchronizing torque coefficient of a machine against an infinite bus, Lemma 1.6. (3) From it, compute the natural frequency and damping ratio of the swing mode, Theorem 1.7. (4) Decide small-signal stability from the sign of two coefficients. (5) Compute a critical clearing angle and time by the equal-area criterion, Theorem 1.8. Chapter 2 asks which of these survive when the rotor is removed.

1.6.3 What this chapter does not cover

Notation of this chapter

Every symbol below is defined in the prose before its first use. This list exists so that a reader can check a symbol without hunting for it. The first block is book-wide notation. The second block is local to this chapter.

Book-wide symbols used in Chapter 1.
SymbolMeaningUnitDefined
Sbasethree-phase base apparent power of a unitVA (MVA)§1.2
Vbaseline-to-line base voltageV (kV)§1.2
Ibasebase currentA§1.2
Zbasebase impedanceohm§1.2
f0rated system frequency; 60 Hz in every example hereHz§1.2
ω0rated electrical angular frequency, 2πf0rad/s§1.2
ωbasebase angular frequency; equal to ω0rad/s§1.2
δrotor angle against the synchronous reference framerad (printed in degrees where noted)§1.3
δ0equilibrium value of δrad§1.4
ωelectrical angular frequency of the unitrad/s§1.3
ωmmechanical angular speed of the rotorrad/s§1.3
Δωper-unit speed deviation (one symbol, not a perturbation)dimensionless§1.3
Jrotor moment of inertiakg m2§1.3
Hinertia constants§1.3
Ekinkinetic energy stored by a unit at rated speedJ (MJ)§1.3
KDdamping coefficient of the per-unit swing equationpu§1.3
Pmmechanical power inputpu§1.3
Peelectrical power outputpu§1.3
Tm, Temechanical and electrical torquepu§1.3
E′constant emf magnitude behind transient reactancepu§1.4
X′ddirect-axis transient reactancepu§1.4
Xttransformer reactancepu§1.4
XLreactance of one transmission linepu§1.4
Xtotal reactance between E′ and the infinite buspu§1.4
Vinfinfinite-bus voltage magnitudepu§1.4
Kssynchronizing torque coefficientpu/rad§1.5
ωnundamped natural frequency of the swing moderad/s§1.5
ζdamping ratio of the swing modedimensionless§1.5
δcrcritical clearing anglerad§1.6
δmaxlargest angle with decelerating area availablerad§1.6
tcrcritical clearing times§1.6
Aacc, Adecaccelerating and decelerating areaspu rad§1.6
Local notation — symbols outside the book-wide table

These symbols are local to Chapter 1. They are listed here because the book-wide notation table does not contain them.


Exercises

Exercise 1.conceptual

A machine operates at δ0 = 95°. State whether it is small-signal stable, and name the term that decides.

Answer target

Not stable. The deciding term is the synchronizing torque coefficient of Lemma 1.6: Ks = (E′Vinf/X) cos 95°, and cos 95° = −0.087156, so Ks < 0. By Theorem 1.7 the constant term of (1.11) is then negative, the product of the two roots is negative, and one root is real and positive.

With the numbers of Example 1.1 (Pmax,pre = 1.692308 pu, H = 3.5 s, KD = 2 pu): Ks = 1.692308 × (−0.0871557) = −0.147494 pu/rad, the polynomial is λ2 + 0.285714λ − 7.943436 = 0, and the roots are λ = +2.679 s−1 and −2.965 s−1. The positive root means the angle runs away. Damping cannot save it: KD appears only in the linear term, and no value of it makes a negative constant term positive.

Exercise 2.derivational

Derive H for a rotor of J = 5000 kg m2 on a four-pole machine at 60 Hz rated 200 MVA.

Answer target

From Definition 1.2, ωm = 2ω/p with p = 4, so ωm = 2π × 60 / 2 = 188.4956 rad/s, that is 1800 rpm. From (1.4), Ekin = ½ × 5000 × 188.49562 = 88.826 MJ. From (1.5), H = 88.826 MJ / 200 MVA = 0.44413 s, that is 0.4441 s.

Comment. 0.4441 s is 22 % of the lower end of the 2 s to 8 s range quoted in Definition 1.3 [S01, Kundur 1994, Power System Stability and Control, Ch. 3 — verify]. Against the H = 3.5 s of Example 1.1, this rotor alone supplies 0.4441/3.5 = 12.7 % of the inertia. The turbine stages, the exciter and the shaft carry the rest. The number to quote for a unit is always the inertia of the complete turbo-set, never of the generator rotor alone.

Exercise 3.derivational

Show that writing the damping term as K′D(ω − ω0), with ω in rad/s, while keeping 2H dΔω/dt on the left-hand side of (1.7), multiplies the computed ζ by ω0.

Answer target

By (1.2), ω − ω0 = ω0Δω. So the mixed equation reads 2H dΔω/dt = Pm − Pe − K′Dω0Δω, which is (1.7) with KD replaced by K′Dω0. Carry that through Theorem 1.7. It leaves ωn unchanged, because ωn does not contain the damping coefficient, and it multiplies ζ = KD/(4Hωn) by ω0.

For Example 1.1: ζ = 376.9911 × 0.0159407 = 6.0095. That would report a machine which rings for seven seconds as strongly overdamped, with no oscillation at all. The two conventions differ by a factor of 377 in this one number, which is why §1.3 fixes one convention and excludes the other.

Exercise 4.computational

Integrate the swing equation for the event of Example 1.2. Find the clearing time at which the machine first fails to return. Use a resolution of 1 ms.

Answer target

The answer lies between 0.169 s and 0.171 s. Use (1.7) with KD = 0, Pe = 0 before the clearing time and Pe = 1.294118 sin δ after it, starting from δ = 0.492383 rad with Δω = 0. Declare failure when δ exceeds a large threshold, for example 200°, and bisect on the clearing time. The reference integration behind Fig. 1.4 uses second-order Runge–Kutta, step 10 μs, failure threshold 200°. It returns 0.1701 s with a 1.5 s window and 0.170085 s with a 3 s window. Report the value your run returns, with the step size, the threshold and the window you used. A short window biases the result upward: a rotor cleared just after tcr lingers near δmax before it runs away. A finite step leaves a residual of a few microseconds.

Why the two values agree. The equal-area result 0.170081 s from (1.14) and the integrated value are the same calculation done two ways, because damping is neglected in both. Adding KD = 2 pu to the integration would push the boundary slightly later, since damping removes energy from the swing. The equal-area value is then a conservative estimate, as Theorem 1.8 states.

Exercise 5.computational

Repeat Example 1.1 for H over 2, 3.5, 6 and 8 s, and plot ωn and ζ. Verify the scaling exponent numerically before you state it.

Answer target

Hold Ks = 1.491276 pu/rad and KD = 2 pu, and use (1.10):

H (s)ωn (rad/s)ωn (Hz)ζ
2.011.85541.88680.02109
3.58.96181.42630.01594
6.06.84471.08940.01217
8.05.92770.94340.01054

Both quantities scale as H−0.5. The numerical check over the two endpoints confirms the exponent: log(5.9277/11.8554)/log(8/2) = −0.50000 for ωn, and log(0.0105438/0.0210875)/log(8/2) = −0.50000 for ζ, using the unrounded values 0.0210875 and 0.0105438 that the table prints as 0.02109 and 0.01054. The rounded table values give −0.5003.

Why. ωn ∝ H−0.5 directly from (1.10). For the damping ratio, ζ = KD/(4Hωn) ∝ H−1 × H+0.5 = H−0.5. Adding inertia therefore slows the oscillation and damps it less, in the same proportion. Inertia is not a substitute for damping.


Every source identifier above refers to §5.2 of the book plan. Every citation carries the mark verify, because no source page was read while this chapter was drafted. Only a human who has read the cited page may remove that mark.