Gauss–Seidel Load Flow

The definition. Every unknown voltage is set to \(1.0\angle0^\circ\) pu:

The slack keeps its specified value; every other bus starts at nominal magnitude and zero angle.

The justification is empirical and strong. A transmission network in normal operation has:

In this case the true solution is 1.018 to 1.047 pu in magnitude and 0 to \(-6.15^\circ\) in angle. The flat start is within 5% and 6° — a good initial guess by any standard.

It also selects the right root. The load-flow equations have multiple solutions, and an iterative method converges to whichever one it starts nearest. The flat start is close to the high-voltage operating solution and far from the low-voltage ones, so it finds the physical answer.

When it is a poor choice, and there are four cases:

The last is not a failure of the flat start but a missed opportunity: a contingency study solving a thousand near-identical cases should start each from the base-case solution, which typically halves the iteration count.

The alternative for hard cases is the DC power flow of Set 17. Solving \(\mathbf{P} = \mathbf{B}\boldsymbol{\theta}\) costs one factorisation, always converges, and gives angles typically within 1° of the AC answer. Starting the AC iteration from those angles with flat magnitudes is standard practice for difficult cases and is sometimes called a "DC start".

What the flat start is not. It is not a claim that the voltages are 1.0; it is a starting point for an iteration whose fixed point is the answer. If the method converges, the starting point affects only how long it took — and, when there are multiple roots, which one was found.

Bus 2, in full. The row of \(\mathbf{Y}_{\text{bus}}\) is

The injection term, with \(S_2 = 0.20+j0.20\) and \(V_2^{(0)} = 1.0\):

The neighbour sum, with \(V_1 = 1.06\) and \(V_3 = V_4 = V_5 = 1.0\):

The update:

A rise of nearly 4% at the first touch, and a small positive angle — bus 2 is net-exporting 0.20 pu, so it leads.

Buses 3, 4 and 5 follow, each using the values updated so far in this same sweep:

Note that bus 3's sum already contains the new \(V_2 = 1.03753\angle0.159^\circ\), not the flat 1.0 — that is Seidel at work.

What the first sweep has achieved. All four angles have moved in the right direction and are ordered correctly by distance from the slack. The magnitudes are close to their final values at bus 2 (1.0375 against 1.0474) and much less so at bus 5 (1.0073 against 1.0179). The information has travelled outward one bus per sweep, and bus 5 is furthest out.

The injection term has changed, because \(V_2\) is no longer 1.0:

Dividing by a larger voltage reduces the current, as it must for fixed power. This term is the only place the nonlinearity enters.

The neighbour sum, now using the first sweep's values:

Against \(-11.13333+j33.40000\) before — a change of about 4% in the real part.

The update:

The full sweep:

The pattern of movement. The magnitude changes have collapsed — from 4% in the first sweep to 0.5% in the second — while the angle changes are almost unchanged, about \(-0.9^\circ\) per bus per sweep. That asymmetry persists throughout:

Angle is a global quantity — bus 5's angle depends on every bus between it and the slack — and information travels one bus per sweep.

Bus 2's angle has crossed zero, from \(+0.159^\circ\) to \(-0.702^\circ\), on its way to \(-2.806^\circ\). The first sweep saw bus 2 exporting into a network still at zero angle and put it ahead; the second sweep, with the load buses now behind, corrects it. An intermediate iterate is not an approximate answer — it is a point on the path, and it may be on the wrong side of the target.

After 43 sweeps at a tolerance of \(10^{-6}\):

The check. Compute the injections back from the solution, \(S_i = V_i(\mathbf{Y}\mathbf{V})_i^{*}\):

Exact to five decimals at every specified bus. That is the only verification a load flow admits, and it must be done — a converged iteration proves the fixed point was reached, not that the fixed point is the right one.

Reading the profile. The voltages fall monotonically with distance from the slack:

A total drop of 4% across the network, comfortably inside a \(\pm5\%\) band. Buses 3 and 4 are within 0.06% of each other — the stiff \(0.01+j0.03\) tie again, and the same near-identity the impedance matrix of Set 18 showed.

And the angles spread over 6.15°, the largest branch difference being \(2.81^\circ\) across line 1–2. Every angle difference is well inside the range where \(\sin\delta \approx \delta\) holds to 1%, which is why the DC flow of Set 17 did as well as it did.

Bus 2 is the highest non-slack bus despite being nearest to the loads, because it is the only other bus with generation. Its \(+0.20\) pu of reactive injection holds it up; without it the voltage there would fall below bus 3's.

The convergence history:

The order. Take the ratio of successive errors after the transient dies:

A constant ratio near 0.78 — the signature of linear (first-order) convergence: \(e_{k+1} = \mu e_k\) with \(\mu \approx 0.78\). Each sweep removes 22% of the error.

The sweep count follows directly. To reduce the error by a factor \(F\):

Each extra decimal place costs about nine more sweeps, forever. There is no acceleration in the tail.

Where \(\mu\) comes from. It is the spectral radius of the iteration matrix, and for a load flow it is governed by how strongly the diagonal of \(\mathbf{Y}_{\text{bus}}\) dominates the rest of its row. Set 16 found that dominance to be marginal — equality without charging, slight deficiency with it. A ratio of 0.78 is exactly what a barely-dominant matrix produces.

And \(\mu\) grows with system size. Information travels one bus per sweep, so a network of diameter \(d\) needs at least \(d\) sweeps for the slack's influence to reach the far end:

The count grows roughly with \(n\), and the cost per sweep also grows with \(n\) — so the total cost is \(O(n^{2})\) against Newton's \(O(n^{1.4})\). That is the whole argument for Set 20.

The modification is one line. Compute the update as usual, then over-shoot it:

With \(\alpha = 1\) this is the plain iteration. With \(\alpha > 1\) the correction is exaggerated; with \(\alpha < 1\) it is damped.

The comparison:

Seventeen sweeps against 43 — a factor of 2.5 for one extra multiplication per bus per sweep.

The first three sweeps are worse. Accelerating a bad estimate overshoots, so the early error is larger. The benefit appears from sweep 4 onward and grows: by sweep 10 the accelerated error is twenty times smaller.

Why it works. Gauss–Seidel's error decays as \(\mu^{k}\) with \(\mu\) real and positive, which means every correction points the same way. The iteration is consistently under-shooting:

Each sweep moves only 22% of the remaining distance. Multiplying the step by \(\alpha\) moves \(0.22\alpha\) of it instead, and the ideal \(\alpha\) would be \(1/0.22 = 4.5\) — if the error decayed uniformly in every component, which it does not.

The reason 4.5 fails. The error has several components decaying at different rates, and the acceleration factor that suits the slowest over-corrects the fastest, driving them unstable. The usable \(\alpha\) is set by the fastest-decaying component rather than the slowest, which is why the practical range is 1.4 to 1.7 rather than 4.5.

The method's proper name is successive over-relaxation, and it applies to any linear stationary iteration. For symmetric positive-definite systems there is a formula for the optimal \(\alpha\) in terms of the Jacobi spectral radius; for the nonlinear load flow there is not, and \(\alpha\) is found by trial — the subject of the next problem.

Scanning \(\alpha\):

A minimum near 1.4, with a sharp penalty beyond 1.6.

The shape of the curve is characteristic of over-relaxation and asymmetric:

The penalty for undershooting is mild and the penalty for overshooting is severe, so the practical choice is deliberately conservative — 1.4 to 1.5 rather than the true optimum.

Why the optimum is not universal. It depends on the network's electrical structure through the spectral radius:

The last case is the surprise: a negative reactance in the network can make the iteration oscillatory, and then damping rather than acceleration is required.

Practical strategies, in order of sophistication:

The third is the most useful in practice: \(\mu\) is measurable from the ratio of consecutive errors, and \(\alpha\) can be reset every few sweeps.

And the diminishing returns. Acceleration takes 43 sweeps to 17 — a factor of 2.5. Newton–Raphson takes the same problem to four iterations, a factor of ten. Tuning \(\alpha\) is optimising the wrong thing once a Newton solver is available, and the historical interest in acceleration factors dates from the years when it was not.

The difficulty. The iteration formula needs \(S_i\), and at a PV bus \(Q_i\) is unknown. It must be computed afresh each sweep from the current voltage estimate:

Using the newest available value of every voltage, exactly as the main iteration does.

The three-step sweep at a PV bus:

Step 3 is what makes it a PV bus. The update gives a voltage of the wrong magnitude but the right angle, and the angle is the quantity actually being solved for.

The solution, after 34 sweeps at \(\alpha = 1.4\):

Comparison with the PQ case of Problem 7:

Every bus is slightly lower, because holding bus 2 at 1.045 rather than letting it rise to 1.047 requires 0.151 pu of reactive support instead of 0.200 — and the 0.05 pu difference is felt across the network.

The convergence is slower — 34 sweeps against 17 — because the rescaling of step 3 throws away part of each correction. A PV bus converges in angle only, and the magnitude information computed at every sweep is discarded.

Which classification is right depends on the machine's control. A generator on automatic voltage regulation holds its terminal voltage and is PV; a machine on constant excitation, or an equivalent representing a fixed reactive injection, is PQ. The two models of the same 40 MW machine give voltages differing by 0.2% — small here, and much larger in a system where the machine is a substantial part of the reactive supply.

The test, applied every sweep:

A bus that hits a limit stops controlling its voltage and becomes an ordinary load bus with a fixed reactive injection.

Here the limit binds. Problem 11 found that holding 1.045 pu requires \(Q_2 = 0.1507\), which exceeds the 0.10 available. So bus 2 becomes PQ with \(S_2 = 0.20+j0.10\):

The three cases side by side:

Every 0.05 pu of reactive support withdrawn from bus 2 costs about 0.0025 pu on every bus in the system.

The behaviour to watch for is limit oscillation. If the bus is PQ at the limit its voltage falls; at the lower voltage the required \(Q\) may drop back inside the limit, so it reverts to PV; then the required \(Q\) rises again and it re-limits. The iteration flips between the two and does not converge.

Three standard remedies:

The third is the crudest and the most common. It can give a slightly wrong answer if a bus limits transiently on the way to the solution, which is why production programs use the first two together.

The engineering significance is larger than the arithmetic. A generator on its reactive limit has lost voltage control, and the network has lost a support point. In a stressed system this is the sequence that leads to collapse: one machine limits, voltages fall, the next machine limits to compensate, and so on. Every large voltage-collapse blackout has this cascade in its timeline, and the load flow shows the first step of it.

The calculation is one row of the nodal equation, done once at the end:

The slack's injection is never used during the iteration — it is an output, not an input.

Evaluating with the converged voltages:

The slack machine generates 129.6 MW and absorbs 7.4 MVAr — it is running under-excited, because the network's charging is producing more vars than the load requires.

The real power balance:

4.59 MW of loss on 145 MW of load — 3.2%. That is high for a transmission system, reflecting this test network's \(x/r = 3\); a real 220 kV network would lose 1–2%.

The reactive balance:

Negative, meaning the network itself generates 17.4 MVAr net — its charging exceeds its \(I^{2}X\) absorption, so the system is operating below its surge impedance loading. Set 15's language, appearing in a load flow.

Why the slack must be able to take it. The slack absorbs whatever the rest of the system does not supply, and here that is 129.6 MW — 89% of the total load. If bus 1's machine cannot produce that, the case is infeasible however well it converges. A converged load flow with an implausible slack output is a wrong case, not a wrong answer.

And the sensitivity to the choice of slack. Choosing bus 5 as slack instead would give the same voltages relative to each other but a different loss allocation and a different set of generator outputs — because the loss depends on the flow pattern, which depends on which bus is free. Loss allocation between market participants is contentious for exactly this reason.

The line current includes the local half-charging:

And \(S_{ji}\) is computed the same way from the other end. The two are not negatives of each other — their sum is the loss in that line.

All seven lines, in per unit:

Balance at bus 1:

Exactly the slack injection of Problem 13. The two calculations are independent — one from the nodal equation and one from the line flows — so the agreement confirms both.

Balance at bus 4, which has a load of \(0.40+j0.05\):

The load, exactly.

Reading the flows. The largest is 88.9 MW on line 1–2, and the smallest is 6.3 MW on 4–5 — the same ranking the DC flow of Set 17 gave. Power leaves the slack on both its circuits and distributes through bus 2, which handles 107.5 MW of throughput across its three outgoing lines.

The reactive flows tell a different story. Several are negative at the sending end and positive at the receiving — line 3–4 sends \(-0.052\) and receives \(+0.032\). That is charging: the line is generating more vars than it consumes, so vars flow out of it at both ends. Every lightly loaded line in this system does the same, which is why the network is a net var source.

The loss in a line is the sum of the two end flows:

Not the difference. Both are defined as flowing into the line from their respective buses, so what remains is what the line consumed.

Line by line:

Matching the slack calculation of Problem 13 exactly.

Three lines carry 81% of the loss. Lines 1–2, 1–3 and 2–5 are the heavily loaded ones, and loss goes as the square of the current:

Line 1–2 carries 0.889 pu and loses 0.0141; line 4–5 carries 0.063 and loses 0.0003. A factor of 14 in flow gives a factor of 45 in loss — close to the square, the discrepancy being the different resistances.

The reactive column is negative almost throughout, meaning those lines generate rather than absorb. The single exception is 2–5, the only line loaded above its own surge impedance loading:

Line 4–5's \(-0.051\) is the extreme case: it carries almost nothing and is nearly pure charging.

The 3.2% loss is high, and the cause is the network rather than the operation. Every line here has \(x/r = 3\), against 8–15 for real transmission conductors. Scaling all seven resistances down by a factor of three would bring the loss to about 1.1%, which is a realistic transmission figure.

What reduces losses, in order of effectiveness on this network:

The two mechanisms, computed from the converged voltages and currents:

Matching the sum of the line-by-line reactive losses in Problem 15 exactly, and confirming that both were computed correctly. The network generates 17.4 MVAr more than it consumes.

The same calculation for real power, as a check on the method:

Exactly the figure obtained from the slack injection in Problem 13, by a completely different route. Two independent calculations of the loss agreeing to five decimals is the strongest available check on a converged solution.

The whole reactive account:

All in per unit. The slack machine runs under-excited, taking 7.4 MVAr out of the network — and the account closes to the last digit.

The balance equation closes:

What the negative slack reactive means operationally. A machine absorbing vars is under-excited, and under-excitation is the region its capability curve most restricts — the limit is set by end-region heating in the stator core rather than by field or armature current. A network that is lightly loaded relative to its charging pushes its machines into that region, which is exactly the condition Set 14's shunt reactors were installed to relieve.

And it identifies the operating regime. The network as a whole is generating vars, so it is running below its surge impedance loading — the same conclusion Set 15 reached from the sign of \(Q_S\). At heavy load the sign would reverse, the machines would go over-excited, and shunt capacitors rather than reactors would be needed.

The angles:

The DC model overstates every angle by 0.13 to 0.48°, and the error grows with distance from the slack.

The branch flows:

Every DC flow is low, by 0.07% to 3.9%. The errors are systematic, not random.

Why every error has the same sign. The DC model has no losses, so the total generation it computes is 1.25 pu against the AC's 1.296. The missing 0.046 pu has to come from somewhere, and it comes out of the flows on the circuits that carry it — which are the ones nearest the slack. That is why 1–2 has the largest error and 3–4 the smallest.

The two error sources, separated:

The third approximation contributes almost nothing at these angles. Nearly all of the error is the loss.

The verdict. Under 4% on every flow, and correct ranking throughout — the DC model identified the same three heavily loaded lines and the same lightly loaded ones. For a contingency screen that is entirely adequate. For a study that must confirm a thermal rating at 95% loading, 4% is not.

And on a real network the DC model does better, not worse. This system's \(x/r = 3\) gives 3.2% losses; a 400 kV network with \(x/r = 12\) loses under 1%, and the DC errors fall in proportion. The DC power flow's accuracy improves with voltage class, which is why it is a transmission tool and not a distribution one.

Scaling the loads, with bus 2's generation held:

At \(\alpha = 1.4\) and a tolerance of \(10^{-8}\).

The critical question: is there a solution at \(k = 3.2\) that the method failed to find, or is there no solution at all?

Testing with a direct nonlinear solver — Newton from the previous \(k\)'s answer, continued in small steps:

The solution ceases to exist at \(k \approx 3.16\), a total load of 4.58 pu against the base 1.45. That is the nose of the P–V curve — the voltage-collapse point of Set 15, found on a five-bus network.

So the failure is the system's, not the method's — beyond \(k = 3.16\). But between \(k = 3.0\) and \(k = 3.16\) a solution exists and Gauss–Seidel does not find it:

A narrow band, but a real one — and it is the band nearest the collapse point, which is exactly where an operator most needs an answer.

Why the method fails there. Near the nose the two solution branches approach each other, the Jacobian becomes near-singular, and the iteration's contraction ratio \(\mu\) rises towards 1. At \(k = 3.0\) it takes 86 sweeps against 24 at \(k = 1.0\); a little further and \(\mu\) exceeds unity and the iteration diverges.

The diagnosis matters. A non-converging load flow admits two very different readings:

Only the second is an engineering finding, and distinguishing them requires either a more robust solver or a continuation method that traces the P–V curve to its nose. Confusing the two is the commonest error in interpreting a failed case.

Per bus per sweep, the operations are:

About \(d+3\) complex operations, and \(d\) averages 3 in a real network. Call it 6 complex operations per bus.

Per sweep, therefore, the cost is proportional to \(n\) — not \(n^{2}\), because the sum runs only over the actual neighbours. Gauss–Seidel exploits sparsity automatically, without any of the ordering machinery of Set 17.

The number of sweeps grows with \(n\) too. Information travels one bus per sweep, so the count is bounded below by the network's diameter, which grows roughly as \(\sqrt{n}\) for a planar network but is observed to give:

So the total is \(O(n^{1.5})\) to \(O(n^{2})\) — and beyond about a thousand buses the iteration count grows faster than the size, and reliability collapses.

Against Newton–Raphson:

The decisive row is the second. Newton's iteration count does not grow with the system, so the crossover comes early — at perhaps 20 to 50 buses, and certainly by 100.

Gauss–Seidel's one remaining advantage is memory: one voltage vector against a factorised Jacobian. That mattered greatly in 1960, when core memory was measured in kilowords, and it is why the method dominated for a decade. It does not matter now.

Where it survives. As an initialiser — two or three Gauss–Seidel sweeps from a flat start cost almost nothing and give Newton a better starting point than the flat profile, improving its reliability on difficult cases. That is its role in modern software, and it is a real one.

The bus solution:

The branch summary, ordered by loading:

Sending-end quantities. Total loss 4.587 MW, 3.2% of the 145 MW load.

The limit checks:

Only bus 1's voltage is at the 1.06 specified value, above the 1.05 operating band — but it is the slack and was set there deliberately.

Finding 1 — the system is lightly loaded relative to its charging. The network generates 17.4 MVAr net and the slack absorbs 7.4. Under-excited operation is the reactive constraint here, not over-excited, and the remedy if it becomes binding is a shunt reactor rather than a capacitor.

Finding 2 — line 1–2 is the critical circuit. It carries 89 MVA, 31% of the losses, and — from the contingency analysis of Set 17 — rises to 125 MW on the loss of 1–3. Its rating governs the whole case.

Finding 3 — buses 3 and 4 are one bus. Their voltages differ by 0.0006 pu and their angles by 0.33°. The \(0.01+j0.03\) tie makes them electrically indistinguishable, which the impedance matrix of Set 18 said and the load flow confirms. Any study not concerned with the flow on that particular line could merge them.

Finding 4 — the loadability margin is 3.16. From Problem 18, the system collapses at 3.16 times its present load. That is a comfortable margin by any standard, and it is the single number a voltage-stability assessment would report.