Electronics Guide

Power Flow Analysis

Power flow analysis, also called load flow analysis, calculates the steady-state operating point of an electric power network. Given the network's impedances, the real power and voltage of each generator, and the power each load draws, it finds the voltage magnitude and angle at every bus. From those follow the real and reactive power in every line and transformer, the losses, and the reactive output each generator must supply. The results answer the first questions asked of any operating condition: are voltages within limits, is any branch overloaded, and what happens if a line trips?

Loads and generators are specified by power rather than by impedance or current, and power is a product of voltage and current, so the equations are nonlinear and must be solved iteratively. A staff paper from the Federal Energy Regulatory Commission (FERC), History of Optimal Power Flow and Formulations by Mary B. Cain, Richard P. O'Neill, and Anya Castillo (December 2012), records that analog network analyzers solved the problem as early as 1929, that Ward and Hale published the first automated digital solution in 1956, and that Newton's method became the common solution method during the 1960s, after sparsity techniques made it practical.

This article develops the standard formulation and the Gauss–Seidel, Newton–Raphson, fast decoupled, and DC solution methods, and it works a three-bus example by each.

What a Power Flow Study Answers

A power flow solution is a snapshot of one operating condition, one pattern of load and generation on one network topology, in sinusoidal steady state. Engineers solve many snapshots, such as peak load, light load, and planned outages, and each application emphasizes different results.

  • Transmission planning. Planners test whether the network, projected years ahead, stays within ratings and voltage limits both intact and after outages.
  • Interconnection studies. A proposed generator, large load, or transmission link is added to a base case to find the problems it causes and the upgrades that cure them.
  • Operations. Operators study next-day conditions and planned outages, and they run real-time contingency analysis on the network state estimated from telemetry.
  • Industrial and commercial systems. Designers check voltages, equipment loading, tap settings, and power factor correction. IEEE 3002.2-2018, IEEE Recommended Practice for Conducting Load-Flow Studies and Analysis of Industrial and Commercial Power Systems, addresses design considerations for new systems, analytical studies of existing ones, and operational and model validation considerations.
  • Distribution planning. Feeder studies check voltage along radial circuits and, increasingly, the voltage rise and reverse flow that distributed generation causes.

A solved power flow also initializes other studies: stability studies start from its operating point, and short-circuit programs can take its voltages as the prefault condition.

Assumptions of the Standard Formulation

  • The three phases are balanced, so a single-phase positive-sequence equivalent represents the network, as described in Three-Phase Circuits and Power.
  • Voltages and currents are sinusoids at nominal frequency.
  • Loads draw constant real and reactive power unless a voltage-dependent load model is added.
  • Each line is a nominal-π equivalent: a series impedance r + jx with half of the line's total charging susceptance B at each end. Each transformer is its leakage impedance and turns ratio.
  • All quantities are in per unit (pu) on one power base, such as 100 MVA, with each voltage level's nominal voltage as its voltage base, so transformers at nominal ratio drop out of the model.

Bus Types and the Slack Bus

Four quantities describe each bus: the net real power injection P, the net reactive power injection Q, the voltage magnitude |V|, and the voltage angle θ. A net injection is generation minus load, so a bus serving a lagging load has negative P and Q. Two quantities are specified at each bus and the solution finds the other two, which defines three bus types.

Bus Types in the Standard Power Flow
Bus type Typical equipment Specified Found by the solution
PQ (load bus)Loads, and generators run at fixed real and reactive outputP, Q|V|, θ
PV (generator or voltage-controlled bus)Generators, inverter-based plants, synchronous condensers, and static compensators that regulate voltageP, |V|Q, θ
Slack (swing or reference bus)One large generator, or an equivalent of a neighboring system|V|, θP, Q

Why One Bus Must Be Slack

The losses are unknown until the voltages are known, so not every generator's real power can be fixed in advance; one bus must supply whatever the losses turn out to be. And because flows depend only on angle differences, one angle must serve as the reference, conventionally zero. The FERC paper gives the first reason and adds that without a slack bus the problem would have more equations than unknowns. In the example below, the slack bus supplies the 4.86 MW of losses in addition to its share of the load.

A single slack bus assigns all of a case's imbalance to one machine, which is unrealistic when the imbalance is large, as after a generator trips. Some programs therefore distribute the slack among several generators by participation factors.

Counting Equations and Unknowns

With n buses, of which nPV are PV buses and nPQ are PQ buses, the unknowns are the angles at the n − 1 non-slack buses and the magnitudes at the PQ buses. The equations are a real power balance at each non-slack bus and a reactive power balance at each PQ bus. Both counts equal n − 1 + nPQ, which is 2nPQ + nPV. The slack bus P and Q and each PV bus's Q then follow by substitution.

The Bus Admittance Matrix

Kirchhoff's current law at every bus, with ground as the reference node, gives the nodal equations

Ii = ∑k Yik Vk, for i = 1, 2, …, n

where Ii is the complex current injected at bus i by its generators and loads, Vk is the complex voltage of bus k, and Yik is an element of the bus admittance matrix, Ybus.

Building Ybus by Inspection

  1. For each branch between buses i and k, compute the series admittance yik = 1 / (rik + jxik).
  2. Set Yik = Yki = −yik, summed over parallel branches. Buses with no branch between them have Yik = 0.
  3. Set Yii to the sum of the series admittances of all branches at bus i, plus jB/2 for each line at the bus, plus any shunt admittance: a positive susceptance for a capacitor bank or a negative one for a reactor.

The matrix is symmetric unless the network contains phase-shifting transformers, and it is sparse: each bus connects to only a few others, so in a large network almost every element is zero, and solvers store and factor only the nonzero ones.

Off-Nominal Taps and Phase Shifters

A transformer whose turns ratio differs from the ratio of the voltage bases, as when its tap is off nominal, needs its own model. Represent it as an ideal transformer of ratio t : 1 at bus i in series with the leakage admittance y toward bus k, so that the voltage between the two is Vi / t. The admittance carries y(Vi / t − Vk) toward bus k, and because the ideal transformer conserves complex power, the current entering at bus i is that current divided by t*, the complex conjugate of t. Collecting terms gives the entries the transformer adds to Ybus.

Entries Added to Ybus by a Transformer with Ratio t : 1 at Bus i
Row Column i Column k
Row iy / |t|2−y / t*
Row k−y / ty

With an in-phase tap, t is real and the matrix stays symmetric. A phase-shifting transformer has a complex ratio, t = a e, so the off-diagonal entries differ and the matrix becomes unsymmetric. Programs differ in which end carries the tap and in the sign convention for the shift, so data must follow the program's convention. MATPOWER, for example, places the tap at the from bus of each branch and, for a branch without charging, builds these same entries.

The Power Flow Equations

The complex power injected at bus i is Si = Pi + jQi = Vi Ii*. Substituting the nodal equation gives

Pi + jQi = Vik Yik* Vk*

With Vi = |Vi| ei, Yik = Gik + jBik, and θik = θi − θk, the real and imaginary parts give the polar power flow equations:

Pi = |Vi| ∑k |Vk| (Gik cos θik + Bik sin θik)

Qi = |Vi| ∑k |Vk| (Gik sin θik − Bik cos θik)

Each sum runs over bus i and its direct neighbors. G and B are the real and imaginary parts of matrix elements, not branch values: for an inductive line, the off-diagonal Bik is positive, and the diagonal Bii is normally negative.

The Two-Bus Case

For a lossless line of reactance X between buses 1 and 2, the powers leaving bus 1 reduce to

P12 = |V1| |V2| sin θ12 / X

Q12 = (|V1|2 − |V1| |V2| cos θ12) / X

Real power flows toward the bus with the lagging angle, up to a maximum of |V1| |V2| / X at a 90° difference. For small angles, cos θ12 is close to 1, so reactive power depends mainly on the difference in magnitudes and flows toward the lower voltage. This pairing of real power with angle and reactive power with magnitude holds approximately wherever reactance is much greater than resistance and angle differences are small, and the fast decoupled and DC methods exploit it.

Consequences of Nonlinearity

  • No general closed-form solution exists, so every method iterates.
  • A solution need not exist. If a case asks the network to deliver more power than it can, no set of voltages satisfies the equations.
  • Several solutions can exist. Besides the normal high-voltage solution there can be low-voltage ones; for a single line serving a load, the two lie on the upper and lower halves of its power–voltage curve.
  • The equations do not enforce equipment limits. In the FERC paper's words, power flow methods find “a mathematically but not necessarily physically feasible or optimal solution.”

Gauss–Seidel and Newton–Raphson Solutions

The Gauss–Seidel Iteration

Conjugating the complex power equation gives Pi − jQi = Vi* Ii. Equating the injected current, (Pi − jQi) / Vi*, to the nodal equation and solving for the diagonal term gives the Gauss–Seidel update

Vi(m+1) = [(Pi − jQi) / (Vi(m))* − ∑k<i Yik Vk(m+1) − ∑k>i Yik Vk(m)] / Yii

where the superscript counts iterations. Each new voltage is used at once in the updates that follow, which distinguishes Gauss–Seidel from the Jacobi method. The usual flat start sets every PQ bus to 1.0 pu and every PV bus to its scheduled magnitude, all at zero angle.

A PV bus has an unknown Q and a fixed magnitude, so its update takes three steps:

  1. Compute Qi = −Im{Vi*k Yik Vk} from the latest voltages, and check it against the reactive limits.
  2. Use that Qi in the update to obtain a new complex voltage.
  3. Keep the new angle, and reset the magnitude to the scheduled value.

Iteration stops when the largest voltage change or, more commonly, the largest power mismatch falls below a tolerance. The FERC paper notes that tolerances are generally set on the mismatches, and the example shows why: the largest voltage change fell below 10−4 pu after seven iterations while the largest mismatch was still 2.4 × 10−4 pu. On a 100 MVA base, 10−4 pu is 0.01 MW or Mvar.

Gauss–Seidel is simple and needs little memory, but it converges linearly, and a correction at one bus reaches distant buses only over many sweeps. Over-relaxation scales each voltage change by an acceleration factor α greater than 1:

Vi,acc(m+1) = Vi(m) + α (Vi(m+1) − Vi(m))

The best factor depends on the network and is found by trial. In the example, with acceleration applied at both buses, α = 1.15 cut the iterations needed to reach a mismatch of 10−6 pu from 12 to 7, while α = 1.6 raised them to 26.

The Newton–Raphson Method

Newton–Raphson linearizes the equations around the present estimate and corrects all the unknowns at once. With the angles θ at the non-slack buses and the magnitudes |V| at the PQ buses as unknowns, the mismatches are

ΔPi = Pi,spec − Pi(θ, |V|), at every PV and PQ bus

ΔQi = Qi,spec − Qi(θ, |V|), at every PQ bus

and a first-order Taylor expansion relates them to the corrections through the Jacobian matrix of partial derivatives, which divides into four blocks:

ΔP = H Δθ + N Δ|V|

ΔQ = M Δθ + L Δ|V|

Differentiating the polar equations gives the elements in the table, with Pi and Qi evaluated at the present estimate.

Elements of the Polar Newton–Raphson Jacobian
Block Off-diagonal element, k ≠ i Diagonal element
H = ∂Pi/∂θk|Vi| |Vk| (Gik sin θik − Bik cos θik)−Qi − Bii |Vi|2
N = ∂Pi/∂|Vk||Vi| (Gik cos θik + Bik sin θik)Pi / |Vi| + Gii |Vi|
M = ∂Qi/∂θk−|Vi| |Vk| (Gik cos θik + Bik sin θik)Pi − Gii |Vi|2
L = ∂Qi/∂|Vk||Vi| (Gik sin θik − Bik cos θik)Qi / |Vi| − Bii |Vi|
  1. Compute Pi and Qi from the present estimate and form the mismatches; stop if the largest is below the tolerance.
  2. Evaluate the Jacobian at the present estimate.
  3. Solve for Δθ and Δ|V| by sparse LU factorization; the Jacobian is never inverted.
  4. Add the corrections to the estimate, and return to step 1.

The Jacobian has 2nPQ + nPV rows, and its four blocks share the sparsity pattern of Ybus, but it is unsymmetric in value and changes at every iteration. Ordering the equations to limit fill-in during factorization keeps storage and time roughly proportional to network size. Tinney and Hart described the method in “Power Flow Solution by Newton's Method” (IEEE Transactions on Power Apparatus and Systems, November 1967). They reported that five iterations, each equivalent to about seven Gauss–Seidel iterations, gave an exact solution, and they credited optimally ordered Gaussian elimination and special programming techniques with making the method practical.

Near a solution at which the Jacobian is nonsingular, Newton's method converges quadratically: each error is roughly proportional to the square of the previous one, so the number of correct digits about doubles per iteration. Far from a solution it can overshoot, and near the network's loading limit, where the Jacobian approaches singularity, its corrections become erratic. From a flat start, the example converges in four iterations.

Fast Decoupled and DC Power Flow

Decoupling

On transmission networks the coupling blocks N and M are small compared with H and L. At the example's flat start, the diagonal entries of H are 17.17 and 11.71 and that of L is 10.62, while no entry of N or M exceeds 2.35 in magnitude. Neglecting N and M splits the Newton step into two smaller problems, ΔP = H Δθ and ΔQ = L Δ|V|.

The Fast Decoupled Load Flow

Stott and Alsaç published the fast decoupled load flow in IEEE Transactions on Power Apparatus and Systems in May 1974. Its usual derivation adds three approximations that hold on typical transmission networks: cos θik ≈ 1; |Gik sin θik| is much smaller than |Bik|; and |Qi| is much smaller than |Bii| |Vi|2. Then Hik ≈ −|Vi| |Vk| Bik and Lik ≈ −|Vi| Bik, diagonal elements included. Dividing each mismatch by the voltage magnitude at its own bus, and setting the remaining |Vk| to 1 in the angle equation, gives

ΔP / |V| = B′ Δθ

ΔQ / |V| = B″ Δ|V|

Both matrices derive from the negative of the imaginary part of Ybus and do not depend on the solution, so each is factored once and reused while the iterations alternate angle and magnitude corrections. Implementations refine them. MATPOWER, for example, builds B′ without shunts or line charging and with taps at nominal, and B″ without phase shifts; its XB version also drops series resistance from B′, and its BX version drops it from B″ instead.

Because the mismatches are always computed from the full equations, the approximations change only the path to the solution: when the method converges, it reaches the Newton–Raphson solution. Stott and Alsaç reported that on practical networks of up to 1,080 buses a solution to within 0.01 MW or Mvar maximum mismatch normally took 4 to 7 iterations, each equal in speed to one and a half Gauss–Seidel iterations or one-fifth of a Newton iteration, and they described the method as attractive for routine and contingency calculations.

DC Power Flow

The DC power flow drops reactive power and voltage magnitude altogether. It assumes that all voltage magnitudes are 1.0 pu, that branch resistance is negligible, that angle differences are small enough for sin θik ≈ θik in radians, and that shunts can be ignored. The flow from bus i to bus k becomes

Pik = (θi − θk) / xik

and the injections satisfy the linear equation P = B′θ, with B′ built from branch reactances as in the XB method and the slack bus row and column removed. One sparse linear solve gives every angle and flow, with no iteration, provided the network is connected. The name is an analogy: the equations match those of a resistive DC circuit, with power in place of current, angle in place of voltage, and reactance in place of resistance.

Linearity allows superposition, which yields sensitivity factors. A power transfer distribution factor gives the change in a branch flow per megawatt transferred between two buses, and a line outage distribution factor gives the change in a branch flow, as a fraction of another branch's pre-outage flow, when that branch trips. The model is lossless, so the slack bus supplies only the difference between load and scheduled generation. Stott, Jardim, and Alsaç, re-examining the model in “DC Power Flow Revisited” (IEEE Transactions on Power Systems, August 2009), described such models as in widespread and increasing use, particularly in congestion-constrained market applications, and cautioned that their MW flows are by no means assured to be reasonably correct.

Where the Approximations Fail

These approximations weaken where resistance is comparable to reactance, as on distribution feeders, and where angle differences are large: at 30°, sin θ is 4.5 percent smaller than θ in radians. The DC model also cannot see voltages, reactive shortfalls, or losses. With every resistance in the example raised to equal its reactance, Newton–Raphson still converges in four iterations, but the fast decoupled method needs 24 iterations to reach 10−8 pu instead of 6. The DC solution does not change, because it ignores resistance, yet the true losses rise to 33.18 MW, and the DC slack generation understates the real one by that amount. The DC model also misses the fall of the bus 3 voltage to 0.835 pu and gives the wrong direction for the small flow on line 1–2. With every resistance raised to twice its reactance, the base-case load has no solution: the network can deliver only about 87 percent of it.

Reactive Limits and Voltage Control

A PV bus holds its voltage only while its equipment can supply or absorb the reactive power needed. A generator's range depends on its real power output and terminal voltage, as its capability curve shows; a static compensator's depends on its rating; and an inverter-based plant's depends on its converter current rating and interconnection requirements. A power flow that ignores these limits can report a voltage profile the equipment cannot deliver.

Switching PV Buses to PQ

  1. After an iteration, compute the reactive generation at each PV bus: its computed Qi plus any reactive load at the bus.
  2. If the generation exceeds Qmax or falls below Qmin, fix it at that limit and treat the bus as a PQ bus, which releases its voltage.
  3. Continue iterating. A bus held at Qmax settles below its voltage setpoint, and a bus held at Qmin settles above it.
  4. Check for back-switching. A bus held at Qmax whose voltage rises above its setpoint, or held at Qmin whose voltage falls below it, can regulate again, so return it to PV.

Early iterations estimate Q poorly, and neighboring PV buses can trade the reactive burden and switch back and forth, so deferring limit checks until the mismatches are small, and capping the number of switches, helps convergence. In the example, generator 2 supplies 15.23 Mvar to hold 1.04 pu. With a 10 Mvar limit, bus 2 becomes a PQ bus, the case re-solves at |V2| = 1.0360 pu and |V3| = 0.9743 pu, and the slack bus's reactive output rises from 35.02 to 40.57 Mvar. Because 1.0360 pu lies below the setpoint, the bus correctly stays at its limit.

Other Controls

Production programs also move transformer taps to hold voltages, switch shunt capacitors and reactors, set phase-shifting transformers to hold scheduled flows, and adjust generation to hold interchange between control areas. They handle these controls in outer loops around the network solution or as extra variables in the Newton iteration. Discrete taps and capacitor steps can hunt between positions.

Worked Example: A Three-Bus System

Three buses are joined in a triangle by three lines. Bus 1 is the slack bus, held at 1.05 pu and zero angle. A generator at bus 2 produces 100 MW and holds 1.04 pu. Bus 3 serves a load of 160 MW and 40 Mvar. Values are in per unit on a 100 MVA base, and each line is a nominal-π model. The lines share one design, so their impedance and charging scale with length, and all three have an X/R ratio of 5. Results are rounded separately, so a sum of rounded entries can differ from a rounded total by 0.01.

Line Data for the Three-Bus Example, per Unit on 100 MVA
Line Resistance r Reactance x Total charging B Series admittance 1/(r + jx)
1–20.020.100.0301.9231 − j9.6154
1–30.040.200.0600.9615 − j4.8077
2–30.030.150.0451.2821 − j6.4103

The Admittance Matrix

For bus 1, Y11 = y12 + y13 + j(0.015 + 0.030) = 2.8846 − j14.3781: the two series admittances plus half of each line's charging susceptance. Each off-diagonal element is the negative of one line's series admittance.

Bus Admittance Matrix of the Example, per Unit
Row Column 1 Column 2 Column 3
Bus 12.8846 − j14.3781−1.9231 + j9.6154−0.9615 + j4.8077
Bus 2−1.9231 + j9.61543.2051 − j15.9881−1.2821 + j6.4103
Bus 3−0.9615 + j4.8077−1.2821 + j6.41032.2436 − j11.1654

One Gauss–Seidel Iteration

  1. Start flat: V2 = 1.04 pu and V3 = 1.00 pu, both at zero angle, with V1 = 1.05 pu fixed.
  2. Bus 2 is a PV bus, so first compute its reactive injection from these voltages: Q2 = −Im{V2*(Y21V1 + Y22V2 + Y23V3)} = 0.1261 pu.
  3. Update bus 2. With Y21V1 + Y23V3 = −3.3013 + j16.5064, V2 = [(1.00 − j0.1261) / 1.04 − (−3.3013 + j16.5064)] / Y22 = 1.0512 + j0.0559, a magnitude of 1.0527 pu at 3.04°.
  4. Keep the angle and reset the magnitude: V2 = 1.04 ej3.04° = 1.0385 + j0.0552 pu.
  5. Update bus 3 with the new V2. With Y31V1 + Y32V2 = −2.6950 + j11.6346, V3 = [(−1.60 + j0.40) / 1.00 − (−2.6950 + j11.6346)] / Y33 = 0.9861 − j0.1001, or 0.9911 pu at −5.79°.

After this sweep the largest mismatch is 0.64 pu. Further sweeps reach 10−4 pu after 8 iterations and 10−8 pu after 15.

Newton–Raphson Iterations

The unknowns are θ2, θ3, and |V3|, so the Jacobian is 3 × 3. At the flat start, the mismatches ΔP2, ΔP3, and ΔQ3 are 0.9667, −1.5006, and 0.1493 pu, and the Jacobian takes the values below.

Jacobian of the Example at the Flat Start
Mismatch ∂/∂θ2 ∂/∂θ3 ∂/∂|V3|
ΔP217.1667−6.6667−1.3333
ΔP3−6.666711.71472.1442
ΔQ31.3333−2.342910.6162
Newton–Raphson Iterations from a Flat Start
Iteration θ2 (deg) θ3 (deg) |V3| (pu) Largest mismatch (pu)
00.00000.00001.000001.5
10.4761−6.92530.986348.8 × 10−2
20.3771−7.16930.976949.3 × 10−4
30.3762−7.17240.976831.2 × 10−7
40.3762−7.17240.976834.0 × 10−15

At iterations 2 and 3, each largest mismatch is about 0.12 to 0.13 times the square of the one before, the signature of quadratic convergence. The last value sits at the limit of double-precision arithmetic.

Converged Solution

Converged Bus Voltages and Net Injections
Bus |V| (pu) θ (deg) Net P (MW) Net Q (Mvar)
1 (slack)1.05000.0064.8635.02
2 (PV)1.04000.38100.0015.23
3 (PQ)0.9768−7.17−160.00−40.00
Line Flows and Losses
Line P leaving first bus (MW) Q leaving first bus (Mvar) P leaving second bus (MW) Q leaving second bus (Mvar) Real loss (MW) Reactive loss (Mvar)
1–2−4.879.844.90−12.980.03−3.13
1–369.7325.17−67.67−21.052.064.12
2–395.1028.21−92.33−18.952.779.26

Each loss is the sum of the powers leaving the two ends of the line. The reactive loss of line 1–2 is negative because this lightly loaded line generates more reactive power in its charging capacitance than its series reactance absorbs. Across all three lines, the series reactances absorb 24.28 Mvar and the charging supplies 14.03 Mvar, a net reactive loss of 10.25 Mvar. The power balances close:

  • Real power: generation of 64.86 + 100.00 = 164.86 MW equals the 160.00 MW load plus losses of 0.03 + 2.06 + 2.77 = 4.86 MW.
  • Reactive power: generation of 35.02 + 15.23 = 50.25 Mvar equals the 40.00 Mvar load plus the net reactive loss of 10.25 Mvar.

Line 1–2 also shows decoupling at work. Real power flows from bus 2 to bus 1, because the bus 2 angle leads by 0.38°, while reactive power flows from bus 1 to bus 2, because bus 1 has the higher voltage.

Comparing the Methods

Iterations to Reach a Given Largest Mismatch
Method 10−4 pu (0.01 MW) 10−8 pu
Gauss–Seidel, α = 1815
Newton–Raphson34
Fast decoupled, XB version36
AC and DC Solutions of the Example
Quantity AC solution DC solution
θ3 (deg)−7.17−7.64
Line 1–3 flow at bus 1 (MW)69.7366.67
Line 2–3 flow at bus 2 (MW)95.1093.33
Line 1–2 flow at bus 1 (MW)−4.87−6.67
Slack generation (MW)64.8660.00
|V3| (pu)0.9768Not modeled (1.0 assumed)

The DC model puts the two heavily loaded lines within 5 percent of their sending-end flows but omits the losses, so its slack generation is 4.86 MW low. Its largest relative error, 37 percent, is on the lightly loaded line 1–2, where an absolute error of 1.8 MW is a large fraction of the flow.

Convergence and Practical Difficulties

Nonconvergence can come from the operating condition, the data, or the controls.

  • No solution exists. When a case demands more than the network can deliver, the mismatches stall or grow while voltages sag. When the example's bus 3 load is raised at constant power factor, with generator 2 held at 100 MW and no reactive limits, Newton–Raphson converges up to about 2.64 times the base load, the nose of the power–voltage curve. There the bus 3 voltage has fallen to about 0.62 pu and the Jacobian becomes singular: its smallest singular value falls from 6.83 at base load to 0.32 at 2.64 times base load and reaches zero at the nose. Beyond the nose no solution exists. Continuation power flow, which MATPOWER includes, traces the solution as load increases to locate such limits.
  • A poor starting point. After large changes, starting from a previous solution or from DC power flow angles helps, and limiting the size of each Newton correction tames early overshoot.
  • Islands. An outage or switching error can split the network. Each island with generation needs its own slack bus; without one, its angles are undefined and the Jacobian is singular.
  • Data errors. Impedances in ohms or on the wrong base, line charging entered as half where the program expects the total, a tap on the wrong side, a reversed phase shift, or a PV bus with no generator can prevent convergence or produce a plausible wrong answer. The bus with the largest remaining mismatch is a good place to start looking.
  • Very small impedances. Bus ties entered as short branches create very large admittances that degrade numerical conditioning; merging the buses avoids the problem.
  • Low-voltage solutions. A case can converge to a mathematically valid solution with voltages far below normal. Such results should be rerun from a better starting point before they are trusted.
  • Control interactions. Limit switching, tap changers, and switched shunts can oscillate. Freezing the controls and releasing them group by group isolates the cause.
  • High resistance. The fast decoupled method slows or fails where resistance is comparable to reactance; Newton–Raphson or a distribution method suits such networks better.

Reading the Results

A review of a converged solution covers the following.

  • Convergence. Confirm the final largest mismatch and where it occurred, and check that every island solved.
  • Voltages. Compare each bus with the criteria set by the planning or operating entity. Low voltages at load buses indicate a reactive shortfall; high voltages at light load often come from line charging.
  • Branch loading. Compare the apparent power or current at both ends with the rating. The ends differ because of losses and charging, so the larger value governs.
  • Reactive resources. Note which generators and compensators sit at their limits, since they have no voltage-control margin left.
  • Slack output. Check that the slack generator's output lies within its capability; a large imbalance there signals an unbalanced dispatch.
  • Angles. Large angle differences across a path mark heavy transfers; stability studies, starting from this operating point, quantify the remaining margin.
  • Losses. Compare real losses between alternatives, and use reactive losses to see where the network consumes reactive power.

In the example, bus 3 sits at 0.9768 pu, the most heavily loaded line, 2–3, carries 99.2 MVA at bus 2, generator 2 supplies 15.23 Mvar, and losses are 4.86 MW, about 3 percent of the load. Whether those values are acceptable depends on criteria and ratings, which contingency analysis applies.

Contingency Analysis and N−1 Security

The NERC Glossary of Terms defines a contingency as “the unexpected failure or outage of a system component, such as a generator, transmission line, circuit breaker, switch or other electrical element.” The N−1 criterion asks whether a system with N major elements can lose any one of them and still operate acceptably: no branch above its emergency rating, no voltage outside its post-contingency limits, and no instability, cascading, or uncontrolled islanding.

NERC reliability standard TPL-001-5.1, Transmission System Planning Performance Requirements, defines the planning events. Category P0 is the normal system with no contingency. Category P1, single contingency, is the loss of a generator, a transmission circuit, a transformer, a shunt device, or a single pole of a DC line. Later categories add the opening of a line section without a fault, bus section and breaker faults, stuck breakers, protection system failures, sequences of two outages with system adjustments between them, and the loss of two adjacent circuits on a common structure or of a bipolar DC line. Steady-state results must respect applicable facility ratings and the voltage limits established by the Planning Coordinator and the Transmission Planner, and in every category the system must remain stable, without cascading or uncontrolled islanding.

The Procedure

  1. Solve and review the base case.
  2. List the contingencies that the criteria require.
  3. For each, remove the element or elements, including any that protection would automatically disconnect, and re-solve from the base-case voltages.
  4. Compare flows with emergency ratings and voltages with post-contingency limits. Record the violations, and treat nonconvergence as a finding, since it can mean that no post-contingency operating point exists.
  5. Develop corrective actions, such as redispatch, switching, reactive support, or new facilities, and verify them with further solutions.

Long contingency lists are screened first. Line outage distribution factors from the DC model estimate every post-outage flow with a few multiplications, fast decoupled solutions give approximate AC answers, and full AC solutions then examine the flagged cases. In operations, an energy management system repeats the analysis periodically, and Artificial Intelligence in Power Control describes learning-based screening intended to speed it further.

N−1 in the Example

Removing each line in turn and re-solving gives the results below, with loading and real power taken at the sending end, where they are largest.

Single-Line Outages in the Three-Bus Example
Outage |V3| (pu) Most heavily loaded line Its loading (MVA) Its real power, AC (MW) Its real power, DC estimate (MW) Total losses (MW)
None (base case)0.9772–399.295.193.34.86
Line 1–20.9772–3103.8100.0100.04.86
Line 1–30.8892–3191.3170.3160.011.18
Line 2–30.8061–3212.5176.7160.018.64

Suppose, for illustration, that each line has an emergency rating of 150 MVA and that post-contingency voltages must stay at or above 0.90 pu. Losing line 1–2 is harmless. Losing line 1–3 overloads line 2–3 and pulls bus 3 below 0.90 pu, and losing line 2–3 does both more severely, so the system is not N−1 secure at this load. The DC estimates flag the same overloads but understate each overloaded line's sending-end flow by that line's own losses, and they cannot see the low bus 3 voltages at all. Reactive support at bus 3, a second circuit, or a lower load would be candidate remedies, each tested with further solutions.

Optimal Power Flow

A power flow evaluates a dispatch that has already been chosen; optimal power flow (OPF) chooses it. OPF minimizes an objective, usually generation cost but sometimes losses, by adjusting controls such as generator real power and voltage setpoints, transformer taps, phase shifters, and switched shunts. The power flow equations enter as equality constraints, and limits on generator output, bus voltages, and branch flows enter as inequality constraints. Security-constrained OPF adds constraints so that the operating point also satisfies the limits after each listed contingency.

The FERC staff paper traces the first full formulation to Carpentier in 1962 and describes the AC OPF as solved in some form every year for system planning, every day for day-ahead markets, every hour, and even every five minutes. Two families of formulation dominate:

  • AC OPF uses the exact AC power flow equations and so represents voltage and reactive power limits directly. It is nonconvex, and in 2012 the FERC authors wrote that there was still no fast, robust solution technique for the full AC OPF.
  • DC OPF substitutes the DC power flow, fixing voltage magnitudes and treating angle differences as small, which gives a linear or quadratic program that solves quickly and reliably but cannot represent voltage limits directly. The FERC paper notes that most implementations of security-constrained economic dispatch include only thermal limits and proxies for voltage limits.

Economic dispatch is the simpler relative. It finds the least-cost generation to serve a given load with a reserve margin, and in the FERC paper's words its formulations “simplify or sometimes altogether ignore power flow constraints.” OPF is, in effect, economic dispatch with the network added. Neural networks that approximate OPF solutions or warm-start its solvers are covered in Artificial Intelligence in Power Control.

Distribution Feeders, Microgrids, and Time-Series Studies

Distribution systems differ from transmission networks in three ways that matter here. Unbalanced loads and single-phase laterals require models of each phase and of the mutual coupling between phases, rather than a positive-sequence equivalent. Resistance comparable to or greater than reactance defeats the decoupling that the fast decoupled and DC methods need. And mostly radial topology permits a forward–backward sweep, which computes branch currents from the feeder ends back toward the substation, then voltage drops from the substation outward, and repeats until the voltages settle, with no Jacobian at all. OpenDSS, which EPRI describes as a distribution system simulator for supporting distributed resource integration and grid modernization efforts, and Eaton's CYME, whose calculation engines Eaton says support balanced or unbalanced models in radial, looped, or meshed configurations, are built for these networks.

Distribution studies also add time, because a snapshot at peak load misses conditions such as overvoltage at minimum daytime load. Quasi-static time-series analysis solves a sequence of steady-state power flows at fixed time steps across a representative period, carrying regulator tap positions and capacitor states from each step to the next. Hosting Capacity and Distribution Impacts explains how utilities use such studies to find how much generation and load a feeder can absorb.

Islanded microgrids break another assumption: that one bus holds frequency and absorbs every imbalance. When grid-forming inverters and generators share load through frequency and voltage droop, no single unit is a slack bus and the steady-state frequency depends on the load, so a suitable power flow treats frequency as an unknown and represents each grid-forming unit by its droop characteristic. Microgrid Power Electronics covers islanded operation and power sharing, and Grid Synchronization and Control covers the droop control laws.

Study Software

Study programs include the following, described in their owners' terms.

  • PSS/E (Siemens): transmission planning and analysis software with power flow analysis for up to 200,000 buses, contingency and fault analysis, and PV/QV studies.
  • PowerWorld Simulator (PowerWorld Corporation): an interactive package that simulates high-voltage power system operation over time frames from several minutes to several days, with a power flow that solves systems of up to 250,000 buses.
  • PowerFactory (DIgSILENT): power system analysis software for generation, transmission, distribution, and industrial systems.
  • ETAP (Operation Technology, Inc.): load flow tools that calculate bus voltages, branch power factors, currents, and power flows throughout the electrical system.
  • CYME (Eaton) and OpenDSS (EPRI): the distribution programs described above. OpenDSS is open source under a BSD license.
  • MATPOWER: free, open-source MATLAB-language code initially developed at Cornell University within the Power Systems Engineering Research Center (PSERC). It runs in MATLAB or GNU Octave, solves power flow, continuation power flow, and optimal power flow problems, and is intended as “a simulation tool for researchers and educators that is easy to use and modify.”
  • pandapower: an open-source Python tool developed jointly by the University of Kassel and Fraunhofer IEE. It builds on the pandas data analysis library, can use the PYPOWER solvers among others, and supports power flow, optimal power flow, state estimation, and short-circuit calculation.

Whatever the program, a study is only as good as its data: impedances, ratings, load forecasts and power factors, reactive limits, and control settings. Published test networks, such as the IEEE 14-bus and 118-bus systems, give standard cases for checking a program before trusting it with real data.

Summary

Power flow analysis finds a network's steady-state voltages from specified powers, and from them the flows, losses, and reactive requirements that decide whether an operating condition is acceptable. PQ, PV, and slack buses define what is known; the slack bus supplies the losses and fixes the angle reference; and the sparse bus admittance matrix turns Kirchhoff's current law into the nonlinear power flow equations.

Gauss–Seidel solves those equations by simple sweeps but converges linearly. Newton–Raphson converges quadratically and, with ordered sparse factorization, handles large networks. The fast decoupled method exploits the weak coupling between real power and voltage magnitude, and between reactive power and angle, to reuse constant matrices, while the DC power flow linearizes the problem, giving up losses, voltages, and some accuracy for speed. Both approximations degrade where resistance approaches reactance.

Practical studies add reactive limits and controls, diagnose nonconvergence, and feed contingency analysis under criteria such as N−1, while optimal power flow wraps the same equations in an optimization that chooses the dispatch. Distribution feeders and islanded microgrids need formulations that handle unbalance, high resistance, and droop-controlled frequency.

Related Topics