N-R load flow corrections and algorithm

The Newton–Raphson (N-R) method linearises the nonlinear power equations about the current estimate by a first-order Taylor series and solves for the corrections in bus voltage angles and magnitudes, repeating until the power mismatches vanish.

Power equations

With $Y_{kn} = |Y_{kn}|\angle\theta_{kn}$ and $V_k = |V_k|\angle\delta_k$:

$$\begin{aligned} P_k &= \sum_{n=1}^{N}|V_k||V_n||Y_{kn}|\cos(\theta_{kn} + \delta_n - \delta_k) \\ Q_k &= -\sum_{n=1}^{N}|V_k||V_n||Y_{kn}|\sin(\theta_{kn} + \delta_n - \delta_k) \end{aligned}$$

So $P_k$ and $Q_k$ are functions of all $\delta$ and $|V|$.

Unknowns and equations

Let bus 1 be the slack bus, buses $2 \ldots m+1$ PV buses ($m$ of them), and the rest PQ buses.

• Unknowns: $\delta_k$ at all $(N-1)$ non-slack buses; $|V_k|$ at the $(N-1-m)$ PQ buses.
• Equations: $P_k$ at the $(N-1)$ non-slack buses; $Q_k$ at the $(N-1-m)$ PQ buses.
• Total: $2N - 2 - m$ equations and unknowns.

Linearisation (Taylor series)

Let $\delta_k^{(0)}$, $|V_k|^{(0)}$ be initial estimates, and $\Delta\delta_k$, $\Delta|V_k|$ the corrections. Expanding $P_k$ about the estimate and neglecting higher-order terms:

$$P_k^{sp} = P_k^{(0)} + \sum_{n=2}^{N}\frac{\partial P_k}{\partial \delta_n}\Delta\delta_n + \sum_{n \in PQ}\frac{\partial P_k}{\partial |V_n|}\Delta|V_n|$$

and similarly for $Q_k$. Writing the mismatches $\Delta P_k = P_k^{sp} - P_k^{(0)}$ and $\Delta Q_k = Q_k^{sp} - Q_k^{(0)}$:

$$\begin{bmatrix} \Delta P \\ \Delta Q \end{bmatrix} = \begin{bmatrix} J_1 & J_2 \\ J_3 & J_4 \end{bmatrix} \begin{bmatrix} \Delta\delta \\ \Delta|V| \end{bmatrix}, \quad J_1 = \frac{\partial P}{\partial \delta},\ J_2 = \frac{\partial P}{\partial |V|},\ J_3 = \frac{\partial Q}{\partial \delta},\ J_4 = \frac{\partial Q}{\partial |V|}$$

General expression for corrections

$$\begin{bmatrix} \Delta\delta \\ \Delta|V| \end{bmatrix}^{(r)} = [J^{(r)}]^{-1} \begin{bmatrix} \Delta P \\ \Delta Q \end{bmatrix}^{(r)}, \qquad \begin{aligned} \delta_k^{(r+1)} &= \delta_k^{(r)} + \Delta\delta_k^{(r)} \\ |V_k|^{(r+1)} &= |V_k|^{(r)} + \Delta|V_k|^{(r)} \end{aligned}$$

In practice the linear system is solved by triangular factorisation (Gaussian elimination), not by explicit inversion.

Elements of the Jacobian

Differentiating the given expressions (write $\phi_{kn} = \theta_{kn} + \delta_n - \delta_k$):

$J_1$ ($\partial P/\partial\delta$), order $(N-1)\times(N-1)$:

$$\frac{\partial P_k}{\partial \delta_k} = \sum_{n \ne k}|V_k||V_n||Y_{kn}|\sin\phi_{kn},\qquad \frac{\partial P_k}{\partial \delta_n} = -|V_k||V_n||Y_{kn}|\sin\phi_{kn}\ (n \ne k)$$

$J_2$ ($\partial P/\partial|V|$), order $(N-1)\times(N-1-m)$:

$$\frac{\partial P_k}{\partial |V_k|} = 2|V_k||Y_{kk}|\cos\theta_{kk} + \sum_{n \ne k}|V_n||Y_{kn}|\cos\phi_{kn},\qquad \frac{\partial P_k}{\partial |V_n|} = |V_k||Y_{kn}|\cos\phi_{kn}$$

$J_3$ ($\partial Q/\partial\delta$), order $(N-1-m)\times(N-1)$:

$$\frac{\partial Q_k}{\partial \delta_k} = \sum_{n \ne k}|V_k||V_n||Y_{kn}|\cos\phi_{kn},\qquad \frac{\partial Q_k}{\partial \delta_n} = -|V_k||V_n||Y_{kn}|\cos\phi_{kn}$$

$J_4$ ($\partial Q/\partial|V|$), order $(N-1-m)\times(N-1-m)$:

$$\frac{\partial Q_k}{\partial |V_k|} = -2|V_k||Y_{kk}|\sin\theta_{kk} - \sum_{n \ne k}|V_n||Y_{kn}|\sin\phi_{kn},\qquad \frac{\partial Q_k}{\partial |V_n|} = -|V_k||Y_{kn}|\sin\phi_{kn}$$

Off-diagonal elements are zero where $Y_{kn} = 0$, so the Jacobian is as sparse as $Y_{bus}$. It is not symmetric and must be recomputed every iteration.

N-R algorithm

1. Read bus and line data; form $Y_{bus}$ and express its elements in polar form.
2. Flat start: $\delta_k^{(0)} = 0$ at all non-slack buses; $|V_k|^{(0)} = 1.0$ at PQ buses; $|V_k| = |V_k|^{sp}$ at PV buses. Set $r = 0$, tolerance $\varepsilon$.
3. Compute $P_k^{(r)}$ (non-slack buses) and $Q_k^{(r)}$ (PQ buses) from the power equations.
4. For PV buses, compute $Q_k^{(r)}$ and check limits. If violated, fix $Q_k$ at the limit and treat the bus as a PQ bus (its $|V_k|$ becomes an unknown, adding one row/column).
5. Compute mismatches $\Delta P_k^{(r)} = P_k^{sp} - P_k^{(r)}$, $\Delta Q_k^{(r)} = Q_k^{sp} - Q_k^{(r)}$.
6. If $\max(|\Delta P|, |\Delta Q|) < \varepsilon$, go to step 10.
7. Evaluate the Jacobian elements $J_1$–$J_4$ with the present $\delta$ and $|V|$.
8. Solve $[J][\Delta\delta;\ \Delta|V|] = [\Delta P;\ \Delta Q]$ for the corrections.
9. Update $\delta_k^{(r+1)} = \delta_k^{(r)} + \Delta\delta_k$, $|V_k|^{(r+1)} = |V_k|^{(r)} + \Delta|V_k|$; set $r = r + 1$ and go to step 3.
10. Compute slack bus power, PV-bus reactive powers, line flows $S_{ik}$ and losses; print results.

 Ybus -> flat start -> P,Q calc -> dP, dQ
                         ^             |
                         |      converged? --yes--> flows
                         |             | no
   update d,|V| <- solve J*dx <- form Jacobian

N-R converges quadratically, usually in 3–5 iterations regardless of system size.