SECTION 5: POWER FLOW

ESE 470 – Energy Distribution Systems

2 Introduction

K. Webb ESE 470

 Nodal Analysis
3

 Consider the following circuit

 Three voltage sources

 𝑉𝑉𝑠𝑠𝑠 , 𝑉𝑉𝑠𝑠𝑠 , 𝑉𝑉𝑠𝑠𝑠
 Generic branch impedances
 Could be any combination of R, L, and C
 Three unknown node voltages
 𝑉𝑉1 , 𝑉𝑉2 , and 𝑉𝑉3
 Would like to analyze the circuit
 Determine unknown node voltages
 One possible analysis technique is nodal analysis
K. Webb ESE 470
 Nodal Analysis
4

 Nodal analysis
 Systematic application of KCL at each unknown node
 Apply Ohm’s law to express branch currents in terms of node
voltages
 Sum currents at each unknown node

 We’ll sum currents leaving each node

and set equal to zero
 At node 𝑉𝑉1 , we have
𝑉𝑉1 − 𝑉𝑉𝑠𝑠𝑠 𝑉𝑉1 − 𝑉𝑉2
+ =0
𝑍𝑍𝑠𝑠𝑠 𝑍𝑍1
 Every current term includes division
by an impedance
 Easier to work with admittances instead
K. Webb ESE 470
 Nodal Analysis
5

 Now our first nodal equation becomes

𝑉𝑉1 − 𝑉𝑉𝑠𝑠𝑠 𝑌𝑌𝑠𝑠𝑠 + 𝑉𝑉1 − 𝑉𝑉2 𝑌𝑌1 = 0
where
𝑌𝑌𝑠𝑠𝑠 = 1/𝑍𝑍𝑠𝑠𝑠 and 𝑌𝑌1 = 1/𝑍𝑍1
 Rearranging to place all unknown node voltages on the left and all
source terms on the right
𝑌𝑌𝑠𝑠𝑠 + 𝑌𝑌1 𝑉𝑉1 − 𝑌𝑌1 𝑉𝑉2 = 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉𝑠𝑠𝑠
 Applying KCL at node 𝑉𝑉2
𝑉𝑉2 − 𝑉𝑉1 𝑌𝑌1 + 𝑉𝑉2 𝑌𝑌2 + 𝑉𝑉2 − 𝑉𝑉𝑠𝑠𝑠 𝑌𝑌𝑠𝑠𝑠 + 𝑉𝑉2 − 𝑉𝑉3 𝑌𝑌3 = 0
K. Webb ESE 470
 Nodal Analysis
6

 Rearranging
−𝑌𝑌1 𝑉𝑉1 + 𝑌𝑌1 + 𝑌𝑌2 + 𝑌𝑌𝑠𝑠𝑠 + 𝑌𝑌3 𝑉𝑉2 − 𝑌𝑌3 𝑉𝑉3 = 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉𝑠𝑠𝑠
 Finally, applying KCL at node 𝑉𝑉3 , gives
𝑉𝑉3 − 𝑉𝑉2 𝑌𝑌3 + 𝑉𝑉3 − 𝑉𝑉𝑠𝑠𝑠 𝑌𝑌𝑠𝑠𝑠 = 0
−𝑌𝑌3 𝑉𝑉2 + 𝑌𝑌3 + 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉3 = 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉𝑠𝑠𝑠
 Note that the source terms are the Norton equivalent
current sources (short-circuit currents) associated with each
voltage source
K. Webb ESE 470
 Nodal Analysis
7

 Putting the nodal equations into matrix form

𝑌𝑌𝑠𝑠1 + 𝑌𝑌𝑌 −𝑌𝑌1 0 𝑉𝑉1 𝑌𝑌𝑠𝑠1 𝑉𝑉𝑠𝑠1
−𝑌𝑌1 𝑌𝑌1 + 𝑌𝑌2 + 𝑌𝑌𝑠𝑠𝑠 + 𝑌𝑌3 −𝑌𝑌3 𝑉𝑉2 = 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉𝑠𝑠𝑠
0 −𝑌𝑌3 𝑌𝑌3 + 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉3 𝑌𝑌𝑠𝑠𝑠 𝑉𝑉𝑠𝑠𝑠

or
𝒀𝒀𝒀𝒀 = 𝑰𝑰
where
 𝒀𝒀 is the 𝑁𝑁 × 𝑁𝑁 admittance matrix
 𝑰𝑰 is an 𝑁𝑁 × 1 vector of known source currents
 𝑽𝑽 is an 𝑁𝑁 × 1 vector of unknown node voltages

 This is a system of 𝑁𝑁 (here, three) linear equations with 𝑁𝑁

unknowns
 We can solve for the vector of unknown voltages as
𝑽𝑽 = 𝒀𝒀−1 𝑰𝑰

K. Webb ESE 470

 The Admittance Matrix, 𝒀𝒀
8

 Take a closer look at the form of the admittance matrix, 𝒀𝒀

𝑌𝑌𝑠𝑠1 + 𝑌𝑌𝑌 −𝑌𝑌1 0 𝑌𝑌11 𝑌𝑌12 𝑌𝑌13
−𝑌𝑌1 𝑌𝑌1 + 𝑌𝑌2 + 𝑌𝑌𝑠𝑠𝑠 + 𝑌𝑌3 −𝑌𝑌3 = 𝑌𝑌21 𝑌𝑌22 𝑦𝑦23
0 −𝑌𝑌3 𝑌𝑌3 + 𝑌𝑌𝑠𝑠𝑠 𝑌𝑌31 𝑌𝑌32 𝑌𝑌33

 The elements of 𝒀𝒀 are

 Diagonal elements, 𝑌𝑌𝑘𝑘𝑘𝑘 :
 𝑌𝑌𝑘𝑘𝑘𝑘 = sum of all admittances connected to node 𝑘𝑘
 Self admittance or driving-point admittance
 Off-diagonal elements, 𝑌𝑌𝑘𝑘𝑘𝑘 (𝑘𝑘 ≠ 𝑛𝑛):
 𝑌𝑌𝑘𝑘𝑘𝑘 = −(total admittance between nodes 𝑘𝑘 and 𝑛𝑛)
 Mutual admittance or transfer admittance

 Note that, because the network is reciprocal, 𝒀𝒀 is symmetric

K. Webb ESE 470
 Nodal Analysis
9

 Nodal analysis allows us to solve for unknown voltages

given circuit admittances and current (Norton
equivalent) inputs
 An application of Ohm’s law
𝒀𝒀𝒀𝒀 = 𝑰𝑰
 A linear equation
 Simple, algebraic solution

 For power-flow analysis, things get a bit more

complicated

K. Webb ESE 470

10 Power-Flow Analysis

K. Webb ESE 470

 The Power-Flow Problem
11

 A typical power system is not entirely unlike the simple

circuit we just looked at
 Sources are generators
 Nodes are the system buses
 Buses are interconnected by impedances of transmission
lines and transformers
 Inputs and outputs now include power (P and Q)
 System equations are now nonlinear
 Can’t simply solve 𝒀𝒀𝒀𝒀 = 𝑰𝑰
 Must employ numerical, iterative solution methods

 Power system analysis to determine bus voltages and

power flows is called power-flow analysis or load-flow
analysis
K. Webb ESE 470
 System One-Line Diagram
12

 Consider the one-line diagram for a simple power system

 System includes:
 Generators
 Buses
 Transformers
 Treated as equivalent circuit impedances in per-unit
 Transmission lines
 Equivalent circuit impedances
 Loads
K. Webb ESE 470
 Bus Variables
13

 The buses are the system nodes

 Four variables associated with each bus, 𝑘𝑘
 Voltage magnitude, 𝑉𝑉𝑘𝑘
 Voltage phase angle, 𝛿𝛿𝑘𝑘

 Real power delivered to the bus, 𝑃𝑃𝑘𝑘

 Reactive power delivered to the bus, 𝑄𝑄𝑘𝑘

K. Webb ESE 470

 Bus Power
14

 Net power delivered to bus 𝑘𝑘 is the difference between power

flowing from generators to bus 𝑘𝑘 and power flowing from bus 𝑘𝑘 to
loads
𝑃𝑃𝑘𝑘 = 𝑃𝑃𝐺𝐺𝐺𝐺 − 𝑃𝑃𝐿𝐿𝐿𝐿
𝑄𝑄𝑘𝑘 = 𝑄𝑄𝐺𝐺𝐺𝐺 − 𝑄𝑄𝐿𝐿𝐿𝐿

 Even though we’ve introduced power flow into the analysis, we can
still write nodal equations for the system
 Voltage and current related by the bus admittance matrix, 𝒀𝒀𝑏𝑏𝑏𝑏𝑏𝑏
𝐈𝐈 = 𝐘𝐘𝑏𝑏𝑏𝑏𝑏𝑏 𝐕𝐕
 𝐘𝐘𝑏𝑏𝑏𝑏𝑏𝑏 contains the bus mutual and self admittances associated with
transmission lines and transformers
 For an 𝑁𝑁 bus system, 𝐕𝐕 is an 𝑁𝑁 × 1 vector of bus voltages
 𝐈𝐈 is an 𝑁𝑁 × 1 vector of source currents flowing into each bus
 From generators and loads

K. Webb ESE 470

 Types of Buses
15

 There are four variables associated with each bus

 𝑉𝑉𝑘𝑘 = 𝑽𝑽𝑘𝑘
 𝛿𝛿𝑘𝑘 = ∠𝑽𝑽𝑘𝑘
 𝑃𝑃𝑘𝑘
 𝑄𝑄𝑘𝑘
 Two variables are inputs to the power-flow problem
 Known
 Two are outputs
 To be calculated
 Buses are categorized into three types depending on which
quantities are inputs and which are outputs
 Slack bus (swing bus)
 Load bus (PQ bus)
 Voltage-controlled bus (PV bus)

K. Webb ESE 470

 Bus Types
16

 Slack bus (swing bus):

 Reference bus
 Typically bus 1
 Inputs are voltage magnitude, 𝑉𝑉1 , and phase angle, 𝛿𝛿1
 Typically 1.0∠0°
 Power, 𝑃𝑃1 and 𝑄𝑄1 , is computed

 Load bus (𝑷𝑷𝑷𝑷 bus):

 Buses to which only loads are connected
 Real power, 𝑃𝑃𝑘𝑘 , and reactive power, 𝑄𝑄𝑘𝑘 , are the knowns
 𝑉𝑉𝑘𝑘 and 𝛿𝛿𝑘𝑘 are calculated
 Majority of power system buses are load buses

K. Webb ESE 470

 Bus Types
17

 Voltage-controlled bus (𝑷𝑷𝑷𝑷 bus):

 Buses connected to generators
 Buses with shunt reactive compensation
 Real power, 𝑃𝑃𝑘𝑘 , and voltage magnitude, 𝑉𝑉𝑘𝑘 , are known
inputs
 Reactive power, 𝑄𝑄𝑘𝑘 , and voltage phase angle, 𝛿𝛿𝑘𝑘 , are
calculated

K. Webb ESE 470

 Solving the Power-Flow Problem
18

 The power-flow solution involves determining:

 𝑉𝑉𝑘𝑘 , 𝛿𝛿𝑘𝑘 , 𝑃𝑃𝑘𝑘 , and 𝑄𝑄𝑘𝑘
 There are 𝑁𝑁 buses
 Each with two unknown quantities
 There are 2𝑁𝑁 unknown quantities in total
 Need 2𝑁𝑁 equations
 𝑁𝑁 of these equations are the nodal equations
𝑰𝑰 = 𝒀𝒀𝒀𝒀 (1)
 The other 𝑁𝑁 equations are the power-balance equations
𝑺𝑺𝑘𝑘 = 𝑃𝑃𝑘𝑘 + 𝑗𝑗𝑄𝑄𝑘𝑘 = 𝑽𝑽𝑘𝑘 𝑰𝑰𝑘𝑘∗ (2)
 From (1), the nodal equation for the 𝑘𝑘𝑡𝑡𝑡 bus is
𝑰𝑰𝑘𝑘 = ∑𝑁𝑁
𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑽𝑽𝑛𝑛 (3)

K. Webb ESE 470

 Solving the Power-Flow Problem
19

 Substituting (3) into (2) gives

𝑃𝑃𝑘𝑘 + 𝑗𝑗𝑄𝑄𝑘𝑘 = 𝑽𝑽𝑘𝑘 ∑𝑁𝑁
𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑽𝑽𝑛𝑛
∗
(4)
 The bus voltages in (3) and (4) are phasors, which we can
represent as
𝑽𝑽𝑛𝑛 = 𝑉𝑉𝑛𝑛 𝑒𝑒 𝑗𝑗𝛿𝛿𝑛𝑛 and 𝑽𝑽𝑘𝑘 = 𝑉𝑉𝑘𝑘 𝑒𝑒 𝑗𝑗𝛿𝛿𝑘𝑘 (5)
 The admittances can also be written in polar form
𝑌𝑌𝑘𝑘𝑘𝑘 = 𝑌𝑌𝑘𝑘𝑘𝑘 𝑒𝑒 𝑗𝑗𝜃𝜃𝑘𝑘𝑘𝑘 (6)
 Using (5) and (6) in (4) gives
∗
𝑃𝑃𝑘𝑘 + 𝑗𝑗𝑄𝑄𝑘𝑘 = 𝑉𝑉𝑘𝑘 𝑒𝑒 𝑗𝑗𝛿𝛿𝑘𝑘 ∑𝑁𝑁 𝑌𝑌
𝑛𝑛=1 𝑘𝑘𝑘𝑘 𝑒𝑒 𝑗𝑗𝜃𝜃𝑘𝑘𝑘𝑘 𝑉𝑉 𝑒𝑒 𝑗𝑗𝛿𝛿𝑛𝑛
𝑛𝑛

𝑃𝑃𝑘𝑘 + 𝑗𝑗𝑄𝑄𝑘𝑘 = 𝑉𝑉𝑘𝑘 ∑𝑁𝑁

𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 𝑒𝑒
𝑗𝑗 𝛿𝛿𝑘𝑘 −𝛿𝛿𝑛𝑛 −𝜃𝜃𝑘𝑘𝑘𝑘
(7)
K. Webb ESE 470
 Solving the Power-Flow Problem
20

 In Cartesian form, (7) becomes

𝑃𝑃𝑘𝑘 + 𝑗𝑗𝑄𝑄𝑘𝑘 =
𝑉𝑉𝑘𝑘 ∑𝑁𝑁
𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 [cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘
+𝑗𝑗 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 ] (8)

 From (8), active power is

𝑃𝑃𝑘𝑘 = 𝑉𝑉𝑘𝑘 ∑𝑁𝑁
𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (9)

 And, reactive power is

𝑄𝑄𝑘𝑘 = 𝑉𝑉𝑘𝑘 ∑𝑁𝑁
𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (10)

K. Webb ESE 470

 Solving the Power-Flow Problem
21

𝑁𝑁
𝑃𝑃𝑘𝑘 = 𝑉𝑉𝑘𝑘 ∑𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (9)
𝑁𝑁
𝑄𝑄𝑘𝑘 = 𝑉𝑉𝑘𝑘 ∑𝑛𝑛=1 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (10)

 Solving the power-flow problem amounts to finding a

solution to a system of nonlinear equations, (9) and (10)
 Must be solved using numerical, iterative algorithms
 Typically Newton-Raphson
 In practice, commercial software packages are available for
power-flow analysis
 E.g. PowerWorld, CYME, ETAP

 We’ll now learn to solve the power-flow problem

 Numerical, iterative algorithm
 Newton-Raphson

K. Webb ESE 470

 Solving the Power-Flow Problem
22

 First, we’ll introduce a variety of numerical algorithms

for solving equations and systems of equations
 Linear system of equations
 Direct solution
 Gaussian elimination
 Iterative solution
 Jacobi
 Gauss-Seidel
 Nonlinear equations
 Iterative solution
 Newton-Raphson
 Nonlinear system of equations
 Iterative solution
 Newton-Raphson
K. Webb ESE 470
 Linear Systems of Equations –
23
Direct Solution

K. Webb ESE 470

 Solving Linear Systems of Equations
24

 Gaussian elimination
 Direct
(i.e. non-iterative) solution
 Two parts to the algorithm:
 Forward elimination
 Back substitution

K. Webb ESE 470

 Gaussian Elimination
25

 Consider a system of equations

−4𝑥𝑥1 + 7𝑥𝑥3 = −5
2𝑥𝑥1 − 3𝑥𝑥2 + 5𝑥𝑥3 = −12
𝑥𝑥2 − 3𝑥𝑥3 = 3

 This can be expressed in matrix form:

−4 0 7 𝑥𝑥1 −5
2 −3 5 𝑥𝑥2 = −12
0 1 −3 𝑥𝑥3 3
 In general
𝐀𝐀 ⋅ 𝐱𝐱 = 𝐲𝐲

 For a system of three equations with three unknowns:

𝐴𝐴11 𝐴𝐴12 𝐴𝐴13 𝑥𝑥1 𝑦𝑦1
𝐴𝐴21 𝐴𝐴22 𝐴𝐴23 𝑥𝑥2 = 𝑦𝑦2
𝐴𝐴31 𝐴𝐴32 𝐴𝐴33 𝑥𝑥3 𝑦𝑦3

K. Webb ESE 470

 Gaussian Elimination
26

 We’ll use a 3×3 system as an example to develop the Gaussian elimination

algorithm
𝐴𝐴11 𝐴𝐴12 𝐴𝐴13 𝑥𝑥1 𝑦𝑦1
𝐴𝐴21 𝐴𝐴22 𝐴𝐴23 𝑥𝑥2 = 𝑦𝑦2
𝐴𝐴31 𝐴𝐴32 𝐴𝐴33 𝑥𝑥3 𝑦𝑦3

 First, create the augmented system matrix

𝐴𝐴11 𝐴𝐴12 𝐴𝐴13 ⋮ 𝑦𝑦1
𝐴𝐴21 𝐴𝐴22 𝐴𝐴23 ⋮ 𝑦𝑦2
𝐴𝐴31 𝐴𝐴32 𝐴𝐴33 ⋮ 𝑦𝑦3

 Each row represents and equation

 𝑁𝑁 rows for 𝑁𝑁 equations
 Row operations do not affect the system
 Multiply a row by a constant
 Add or subtract rows from one another and replace row with the result

K. Webb ESE 470

 Gaussian Elimination – Forward Elimination
27

 Perform row operations to reduce the augmented

matrix to upper triangular
 Only zeros below the main diagonal
 Eliminate 𝑥𝑥𝑖𝑖 from the 𝑖𝑖 + 1 st through the 𝑁𝑁 th equations for
𝑖𝑖 = 1 … 𝑁𝑁
 Forward elimination

 After forward elimination, we have

𝐴𝐴11 𝐴𝐴12 𝐴𝐴13 ⋮ 𝑦𝑦1
0 𝐴𝐴′22 𝐴𝐴′23 ⋮ 𝑦𝑦2′
0 0 𝐴𝐴′33 ⋮ 𝑦𝑦3′
 Where the prime notation (e.g. 𝐴𝐴′22 ) indicates that the value
has been changed from its original value
K. Webb ESE 470
 Gaussian Elimination – Back Substitution
28

𝐴𝐴11 𝐴𝐴12 𝐴𝐴13 ⋮ 𝑦𝑦1

0 𝐴𝐴′22 𝐴𝐴′23 ⋮ 𝑦𝑦2′
0 0 𝐴𝐴′33 ⋮ 𝑦𝑦3′
 The last row represents an equation with only a single unknown
𝐴𝐴′33 ⋅ 𝑥𝑥3 = 𝑦𝑦3′
 Solve for 𝑥𝑥3
𝑦𝑦3′
𝑥𝑥3 = ′
𝐴𝐴33
 The second-to-last row represents an equation with two unknowns
𝐴𝐴′22 ⋅ 𝑥𝑥2 + 𝐴𝐴′23 ⋅ 𝑥𝑥3 = 𝑦𝑦2′
 Substitute in newly-found value of 𝑥𝑥3
 Solve for 𝑥𝑥2
 Substitute values for 𝑥𝑥2 and 𝑥𝑥3 into the first-row equation
 Solve for 𝑥𝑥1
 This process is back substitution
K. Webb ESE 470
 Gaussian elimination
29

 Gaussian elimination summary

 Create the augmented system matrix
 Forward elimination
 Reduce to an upper-triangular matrix
 Back substitution
 Starting with 𝑥𝑥𝑁𝑁 , solve for 𝑥𝑥𝑖𝑖 for 𝑖𝑖 = 𝑁𝑁 … 1

 A direct solution algorithm

 Exact value for each 𝑥𝑥𝑖𝑖 arrived at with a single execution of
the algorithm
 Alternatively, we can use an iterative algorithm
 The Jacobi method

K. Webb ESE 470

 Linear Systems of Equations –
30
Iterative Solution – Jacobi Method

K. Webb ESE 470

 Jacobi Method
31

 Consider a system of 𝑁𝑁 linear equations

𝐀𝐀 ⋅ 𝐱𝐱 = 𝐲𝐲
𝐴𝐴1,1 ⋯ 𝐴𝐴1,𝑁𝑁 𝑥𝑥1 𝑦𝑦1
⋮ ⋱ ⋮ ⋮ = ⋮
𝐴𝐴𝑁𝑁,1 ⋯ 𝐴𝐴𝑁𝑁,𝑁𝑁 𝑥𝑥𝑁𝑁 𝑦𝑦𝑁𝑁

 The 𝑘𝑘th equation (𝑘𝑘th row) is

𝐴𝐴𝑘𝑘,1 𝑥𝑥1 + 𝐴𝐴𝑘𝑘,2 𝑥𝑥2 + ⋯ + 𝐴𝐴𝑘𝑘,𝑘𝑘 𝑥𝑥𝑘𝑘 + ⋯ + 𝐴𝐴𝑘𝑘,𝑁𝑁 𝑥𝑥𝑁𝑁 = 𝑦𝑦𝑘𝑘 (11)

 Solve (11) for 𝑥𝑥𝑘𝑘

1
𝑥𝑥𝑘𝑘 = [𝑦𝑦𝑘𝑘 − (𝐴𝐴𝑘𝑘,1 𝑥𝑥1 + 𝐴𝐴𝑘𝑘,2 𝑥𝑥2 + ⋯ + 𝐴𝐴𝑘𝑘,𝑘𝑘−1 𝑥𝑥𝑘𝑘−1 + (12)
𝐴𝐴𝑘𝑘,𝑘𝑘
+𝐴𝐴𝑘𝑘,𝑘𝑘+1 𝑥𝑥𝑘𝑘+1 + ⋯ + 𝐴𝐴𝑘𝑘,𝑁𝑁 𝑥𝑥𝑁𝑁 )]
K. Webb ESE 470
 Jacobi Method
32

 Simplify (12) using summing notation

𝑘𝑘−1 𝑁𝑁
1
𝑥𝑥𝑘𝑘 = 𝑦𝑦 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛 , 𝑘𝑘 = 1 … 𝑁𝑁 (13)
𝐴𝐴𝑘𝑘,𝑘𝑘 𝑘𝑘
𝑛𝑛=1 𝑛𝑛=𝑘𝑘+1

 An equation for 𝑥𝑥𝑘𝑘

 But, of course, we don’t yet know all other 𝑥𝑥𝑛𝑛 values
 Use (13) as an iterative expression
𝑘𝑘−1 𝑁𝑁
1
𝑥𝑥𝑘𝑘,𝑖𝑖+1 = 𝑦𝑦 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 , 𝑘𝑘 = 1 … 𝑁𝑁 (14)
𝐴𝐴𝑘𝑘,𝑘𝑘 𝑘𝑘
𝑛𝑛=1 𝑛𝑛=𝑘𝑘+1

 The 𝑖𝑖 subscript indicates iteration number

 𝑥𝑥𝑘𝑘,𝑖𝑖+1 is the updated value from the current iteration
 𝑥𝑥𝑛𝑛,𝑖𝑖 is a value from the previous iteration
K. Webb ESE 470
 Jacobi Method
33

𝑘𝑘−1 𝑁𝑁
1
𝑥𝑥𝑘𝑘,𝑖𝑖+1 = 𝑦𝑦𝑘𝑘 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 , 𝑘𝑘 = 1 … 𝑁𝑁 (14)
𝐴𝐴𝑘𝑘,𝑘𝑘
𝑛𝑛=1 𝑛𝑛=𝑘𝑘+1

 Old values of 𝑥𝑥𝑛𝑛 , on the right-hand side, are used to

update 𝑥𝑥𝑘𝑘 on the left-hand side
 Start with an initial guess for all unknowns, 𝐱𝐱 0
 Iterate until adequate convergence is achieved
a specified stopping criterion is satisfied
 Until

 Convergence is not guaranteed

K. Webb ESE 470

 Convergence
34

 An approximation of 𝐱𝐱 is refined on each iteration

 Continue to iterate until we’re close to the right answer
for the vector of unknowns, 𝐱𝐱
 Assume we’ve converged to the right answer when 𝐱𝐱
changes very little from iteration to iteration
 On each iteration, calculate a relative error quantity
𝑥𝑥𝑘𝑘,𝑖𝑖+1 − 𝑥𝑥𝑘𝑘,𝑖𝑖
𝜀𝜀𝑖𝑖 = max , 𝑘𝑘 = 1 … 𝑁𝑁
𝑥𝑥𝑘𝑘,𝑖𝑖
 Iterate until
𝜀𝜀𝑖𝑖 ≤ 𝜀𝜀𝑠𝑠
where 𝜀𝜀𝑠𝑠 is a chosen stopping criterion
K. Webb ESE 470
 Jacobi Method – Matrix Form
35

 The Jacobi method iterative formula, (14), can be rewritten in matrix form:
𝐱𝐱 𝑖𝑖+1 = 𝐌𝐌𝐱𝐱 𝑖𝑖 + 𝐃𝐃−1 𝐲𝐲 (15)
where 𝐃𝐃 is the diagonal elements of A
𝐴𝐴1,1 0 ⋯ 0
0 𝐴𝐴2,2 0 ⋮
𝐃𝐃 =
⋮ 0 ⋱ 0
0 ⋯ 0 𝐴𝐴𝑁𝑁,𝑁𝑁
and
𝐌𝐌 = 𝐃𝐃−1 𝐃𝐃 − 𝐀𝐀 (16)
 Recall that the inverse of a diagonal matrix is given by inverting each diagonal
element
1/𝐴𝐴1,1 0 ⋯ 0
0 1/𝐴𝐴2,2 0 ⋮
𝐃𝐃−𝟏𝟏 =
⋮ 0 ⋱ 0
0 ⋯ 0 1/𝐴𝐴𝑁𝑁,𝑁𝑁

K. Webb ESE 470

 Jacobi Method – Example
36

 Consider the following system of equations

−4𝑥𝑥1 + 7𝑥𝑥3 = −5
2𝑥𝑥1 − 3𝑥𝑥2 + 5𝑥𝑥3 = −12
𝑥𝑥2 − 3𝑥𝑥3 = 3

 In matrix form:
−4 0 7 𝑥𝑥1 −5
2 −3 5 𝑥𝑥2 = −12
0 1 −3 𝑥𝑥3 3
 Solve using the Jacobi method

K. Webb ESE 470

 Jacobi Method – Example
37

 The iteration formula is

𝐱𝐱 𝑖𝑖+1 = 𝐌𝐌𝐱𝐱 𝑖𝑖 + 𝐃𝐃−1 𝐲𝐲
where
−4 0 0 −0.25 0 0
𝐃𝐃 = 0 −3 0 𝐃𝐃−1 = 0 −0.333 0
0 0 −3 0 0 −0.333
0 0 1.75
𝐌𝐌 = 𝐃𝐃−1 𝐃𝐃 − 𝐀𝐀 = 0.667 0 1.667
0 0.333 0
 To begin iteration, we need a starting point
 Initial guess for unknown values, 𝐱𝐱
 Often, we have some idea of the answer
 Here, arbitrarily choose
𝑇𝑇
𝐱𝐱 0 = 10 25 10
K. Webb ESE 470
 Jacobi Method – Example
38

 At each iteration, calculate

𝐱𝐱 𝑖𝑖+1 = 𝐌𝐌𝐱𝐱 𝑖𝑖 + 𝐃𝐃−1 𝐲𝐲
𝑥𝑥1,𝑖𝑖+1 0 0 1.75 𝑥𝑥1,𝑖𝑖 1.25
𝑥𝑥2,𝑖𝑖+1 = 0.667 0 1.667 𝑥𝑥2,𝑖𝑖 + 4
𝑥𝑥3,𝑖𝑖+1 0 0.333 0 𝑥𝑥3,𝑖𝑖 −1
 For 𝑖𝑖 = 1:
𝑥𝑥1,1 0 0 1.75 10 1.25
𝐱𝐱1 = 𝑥𝑥2,1 = 0.667 0 1.667 25 + 4
𝑥𝑥3,1 0 0.333 0 10 −1
𝑇𝑇
𝐱𝐱1 = 18.75 27.33 7.33

 The relative error is

𝑥𝑥𝑘𝑘,1 − 𝑥𝑥𝑘𝑘,0
𝜀𝜀1 = max = 0.875
𝑥𝑥𝑘𝑘,0

K. Webb ESE 470

 Jacobi Method – Example
39

 For 𝑖𝑖 = 2:
𝑥𝑥1,2 0 0 1.75 18.75 1.25
𝐱𝐱 2 = 𝑥𝑥2,2 = 0.667 0 1.667 27.33 + 4
𝑥𝑥3,2 0 0.333 0 7.33 −1
𝑇𝑇
𝐱𝐱 2 = 14.08 28.72 8.11

 The relative error is

𝑥𝑥𝑘𝑘,2 − 𝑥𝑥𝑘𝑘,1
𝜀𝜀2 = max = 0.249
𝑥𝑥𝑘𝑘,1

 Continue to iterate until relative error falls below a specified

stopping condition

K. Webb ESE 470

 Jacobi Method – Example
40

 Automate with computer code, e.g. MATLAB

 Setup the system of equations

 Initialize matrices and parameters for iteration

K. Webb ESE 470

 Jacobi Method – Example
41

 Loop to continue iteration as long as:

 Stopping criterion is not satisfied
 Maximum number of iterations is not exceeded

 On each iteration
 Use previous 𝐱𝐱 values to update 𝐱𝐱
 Calculate relative error
 Increment the number of iterations

K. Webb ESE 470

 Jacobi Method – Example
42

 Set 𝜀𝜀𝑠𝑠 = 1 × 10−6 and iterate:

𝒊𝒊 𝐱𝐱 𝒊𝒊 𝜺𝜺𝒊𝒊
0 𝑇𝑇 -
10 25 10
1 𝑇𝑇 0.875
18.75 27.33 7.33
2 𝑇𝑇 0.249
14.08 28.72 8.11
3 𝑇𝑇 0.097
15.44 26.91 8.57
4 𝑇𝑇 0.071
16.25 28.59 7.97
5 𝑇𝑇 0.070
15.20 28.12 8.53
6 𝑇𝑇 0.065
16.18 28.35 8.37
⋮ ⋮ ⋮
371 𝑇𝑇 0.995×10-6
20.50 36.00 11.00

 Convergence achieved in 371 iterations

K. Webb ESE 470
 Linear Systems of Equations –
43
Iterative Solution – Gauss-Seidel

K. Webb ESE 470

 Gauss-Seidel Method
44

 The iterative formula for the Jacobi method is

𝑘𝑘−1 𝑁𝑁
1
𝑥𝑥𝑘𝑘,𝑖𝑖+1 = 𝑦𝑦 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 , 𝑘𝑘 = 1 … 𝑁𝑁 (14)
𝐴𝐴𝑘𝑘,𝑘𝑘 𝑘𝑘
𝑛𝑛=1 𝑛𝑛=𝑘𝑘+1

 Note that only old values of 𝑥𝑥𝑛𝑛 (i.e. 𝑥𝑥𝑛𝑛,𝑖𝑖 ) are used to
update the value of 𝑥𝑥𝑘𝑘
 Assume the 𝑥𝑥𝑘𝑘,𝑖𝑖+1 values are determined in order of
increasing 𝑘𝑘
 When updating 𝑥𝑥𝑘𝑘,𝑖𝑖+1 , all 𝑥𝑥𝑛𝑛,𝑖𝑖+1 values are already known
for 𝑛𝑛 < 𝑘𝑘
 We can use those updated values to calculate 𝑥𝑥𝑘𝑘,𝑖𝑖+1
 The Gauss-Seidel method

K. Webb ESE 470

 Gauss-Seidel Method
45

 Now use the 𝑥𝑥𝑛𝑛 values already updated on the

current iteration to update 𝑥𝑥𝑘𝑘
 That is, 𝑥𝑥𝑛𝑛,𝑖𝑖+1 for 𝑛𝑛 < 𝑘𝑘
 Gauss-Seidel iterative formula
𝑘𝑘−1 𝑁𝑁
1
𝑥𝑥𝑘𝑘,𝑖𝑖+1 = 𝑦𝑦 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖+1 − � 𝐴𝐴𝑘𝑘,𝑛𝑛 𝑥𝑥𝑛𝑛,𝑖𝑖 , 𝑘𝑘 = 1 … 𝑁𝑁 (17)
𝐴𝐴𝑘𝑘,𝑘𝑘 𝑘𝑘
𝑛𝑛=1 𝑛𝑛=𝑘𝑘+1

 Note that only the first summation has changed

 For already updated 𝑥𝑥 values
 𝑥𝑥𝑛𝑛 for 𝑛𝑛 < 𝑘𝑘
 Number of already-updated values used depends on 𝑘𝑘

K. Webb ESE 470

 Gauss-Seidel – Matrix Form
46

 In matrix form the iterative formula is the same as for the Jacobi
method
𝐱𝐱 𝑖𝑖+1 = 𝐌𝐌𝐱𝐱 𝑖𝑖 + 𝐃𝐃−1 𝐲𝐲 (15)
where, again
𝐌𝐌 = 𝐃𝐃−1 𝐃𝐃 − 𝐀𝐀 (16)

but now 𝐃𝐃 is the lower triangular part of 𝐀𝐀

𝐴𝐴1,1 0 ⋯ 0
𝐴𝐴 𝐴𝐴2,2 0 ⋮
𝐃𝐃 = 2,1
⋮ ⋮ ⋱ 0
𝐴𝐴𝑁𝑁,1 𝐴𝐴𝑁𝑁,2 ⋯ 𝐴𝐴𝑁𝑁,𝑁𝑁

 Otherwise, the algorithm and computer code is identical to that of

the Jacobi method

K. Webb ESE 470

 Gauss-Seidel – Example
47

 Apply Gauss-Seidel to our previous example

𝑇𝑇
 𝑥𝑥0 = 10 25 10
 𝜀𝜀𝑠𝑠 = 1 × 10−6
𝒊𝒊 𝐱𝐱 𝒊𝒊 𝜺𝜺𝒊𝒊
0 𝑇𝑇 -
10 25 10
1 𝑇𝑇 0.875
18.75 33.17 10.06
2 𝑇𝑇 0.005
18.85 33.32 10.11
3 𝑇𝑇 0.005
18.94 33.47 10.16
4 𝑇𝑇 0.005
19.03 33.61 10.20
⋮ ⋮ ⋮
151 𝑇𝑇 0.995×10-6
20.50 36.00 11.00

 Convergence achieved in 151 iterations

 Compared to 371 for the Jacobi method
K. Webb ESE 470
48 Nonlinear Equations

K. Webb ESE 470

 Nonlinear Equations
49

 Solution methods we’ve seen so far work only for

linear equations
 Now, we introduce an iterative method for solving a
single nonlinear equation
 Newton-Raphson method

 Next, we’ll apply the Newton-Raphson method to a

system of nonlinear equations
 Finally, we’ll use Newton-Raphson to solve the
power-flow problem
K. Webb ESE 470
 Newton-Raphson Method
50

 Want to solve
𝑦𝑦 = 𝑓𝑓(𝑥𝑥)
where 𝑓𝑓 𝑥𝑥 is a nonlinear function
 That is, we want to find 𝑥𝑥, given a known nonlinear
function, 𝑓𝑓, and a known output, 𝑦𝑦
 Newton-Raphson method
 Based on a first-order Taylor series approximation to
𝑓𝑓 𝑥𝑥
 The nonlinear 𝑓𝑓 𝑥𝑥 is approximated as linear to update
our approximation to the solution, 𝑥𝑥, on each iteration
K. Webb ESE 470
 Taylor Series Approximation
51

 Taylor series approximation

 Given:
A function, 𝑓𝑓 𝑥𝑥
 Value of the function at some value of 𝑥𝑥, 𝑓𝑓 𝑥𝑥0
 Approximate:
 Value of the function at some other value of 𝑥𝑥

 First-order Taylor series approximation

 Approximate 𝑓𝑓 𝑥𝑥 using only its first derivative
 𝑓𝑓 𝑥𝑥 approximated as linear – constant slope

𝑑𝑑𝑑𝑑
𝑦𝑦 = 𝑓𝑓 𝑥𝑥 ≈ 𝑓𝑓 𝑥𝑥0 + � 𝑥𝑥 − 𝑥𝑥0 = 𝑦𝑦�
𝑑𝑑𝑑𝑑 𝑥𝑥=𝑥𝑥0
K. Webb ESE 470
 First-Order Taylor Series Approximation
52

 Approximate value of the function at 𝑥𝑥

𝑓𝑓 𝑥𝑥 ≈ 𝑦𝑦� = 𝑓𝑓 𝑥𝑥0 + 𝑓𝑓 ′ 𝑥𝑥0 𝑥𝑥 − 𝑥𝑥0

K. Webb ESE 470

 Newton-Raphson Method
53

 First order Taylor series approximation is

𝑦𝑦 ≈ 𝑓𝑓 𝑥𝑥0 + 𝑓𝑓 ′ 𝑥𝑥0 𝑥𝑥 − 𝑥𝑥0

 Letting this be an equality and rearranging gives an iterative formula

for updating an approximation to 𝑥𝑥
𝑦𝑦 = 𝑓𝑓 𝑥𝑥𝑖𝑖 + 𝑓𝑓 ′ 𝑥𝑥𝑖𝑖 𝑥𝑥𝑖𝑖+1 − 𝑥𝑥𝑖𝑖

𝑓𝑓 ′ 𝑥𝑥𝑖𝑖 𝑥𝑥𝑖𝑖+1 − 𝑥𝑥𝑖𝑖 = 𝑦𝑦 − 𝑓𝑓 𝑥𝑥𝑖𝑖

1
𝑥𝑥𝑖𝑖+1 = 𝑥𝑥𝑖𝑖 + ′
𝑦𝑦 − 𝑓𝑓 𝑥𝑥𝑖𝑖 (18)
𝑓𝑓 𝑥𝑥𝑖𝑖
 Initialize with a best guess at 𝑥𝑥, 𝑥𝑥0
 Iterate (18) until
 A stopping criterion is satisfied, or
 The maximum number of iterations is reached

K. Webb ESE 470

 First-Order Taylor Series Approximation
54

1
𝑥𝑥𝑖𝑖+1 = 𝑥𝑥𝑖𝑖 + ′ 𝑦𝑦 − 𝑓𝑓 𝑥𝑥𝑖𝑖
𝑓𝑓 𝑥𝑥𝑖𝑖

 Now using the Taylor

series approximation
in a different way
 Not approximating
the value of y = f(x)
at x, but, instead
 Approximating the
value of x where
f(x) = y

K. Webb ESE 470

 Newton-Raphson – Example
55

 Consider the following nonlinear equation

𝑦𝑦 = 𝑓𝑓 𝑥𝑥 = 𝑥𝑥 3 + 10 = 20
 Apply Newton-Raphson to solve
 Find 𝑥𝑥, such that 𝑦𝑦 = 𝑓𝑓 𝑥𝑥 = 20
 The derivative function is
𝑓𝑓 ′ 𝑥𝑥 = 3𝑥𝑥 2
 Initial guess for 𝑥𝑥
𝑥𝑥0 = 1
 Iterate using the formula given by (18)
K. Webb ESE 470
 Newton-Raphson – Example
56

 𝑖𝑖 = 1:
𝑥𝑥1 = 𝑥𝑥0 + 𝑓𝑓 ′ 𝑥𝑥0 −1
𝑦𝑦 − 𝑓𝑓 𝑥𝑥0

𝑥𝑥1 = 1 + 3 ⋅ 12 −1
20 − 13 + 10

𝑥𝑥1 = 4

𝑥𝑥1 − 𝑥𝑥0
𝜀𝜀1 =
𝑥𝑥0

4−1
𝜀𝜀1 = =3
1
𝑥𝑥1 = 4, 𝜀𝜀1 = 3

K. Webb ESE 470

 Newton-Raphson – Example
57

 𝑖𝑖 = 2:
𝑥𝑥2 = 𝑥𝑥1 + 𝑓𝑓 ′ 𝑥𝑥1 −1
𝑦𝑦 − 𝑓𝑓 𝑥𝑥1
𝑥𝑥2 = 4 + 3 ⋅ 42 −1
20 − 43 + 10
𝑥𝑥2 = 2.875
𝑥𝑥2 − 𝑥𝑥1
𝜀𝜀2 =
𝑥𝑥1
2.875 − 4
𝜀𝜀2 =
4
𝜀𝜀2 = 0.281
𝑥𝑥2 = 2.875, 𝜀𝜀2 = 0.281

K. Webb ESE 470

 Newton-Raphson – Example
58

 𝑖𝑖 = 3:
𝑥𝑥3 = 𝑥𝑥2 + 𝑓𝑓 ′ 𝑥𝑥2 −1
𝑦𝑦 − 𝑓𝑓 𝑥𝑥2

𝑥𝑥3 = 2.875 + 3 ⋅ 2.8752 −1

20 − 2.8753 + 10

𝑥𝑥3 = 2.32

𝑥𝑥3 − 𝑥𝑥2
𝜀𝜀3 =
𝑥𝑥2

2.32 − 2.875
𝜀𝜀3 =
2.875

𝜀𝜀3 = 0.193

𝑥𝑥3 = 2.32, 𝜀𝜀3 = 0.193

K. Webb ESE 470

 Newton-Raphson – Example
59

 𝑖𝑖 = 4:
𝑥𝑥4 = 2.166, 𝜀𝜀4 = 0.066
 𝑖𝑖 = 5:
𝑥𝑥5 = 2.155, 𝜀𝜀5 = 0.005
 𝑖𝑖 = 6:
𝑥𝑥6 = 2.154, 𝜀𝜀6 = 28.4 × 10−6
 𝑖𝑖 = 7:
𝑥𝑥7 = 2.154, 𝜀𝜀7 = 0.808 × 10−9

 Convergence achieved very quickly

 Next, we’ll see how to apply Newton-Raphson to a system
of nonlinear equations
K. Webb ESE 470
60 Example Problems

K. Webb ESE 470

 Perform three iterations toward the solution of the following system
of equations using the Jacobi method. Let 𝐱𝐱 0 = 0, 0 𝑇𝑇 .
2𝑥𝑥1 + 𝑥𝑥2 = 12
2𝑥𝑥1 + 3𝑥𝑥2 = 5

61 K. Webb ESE 470

62 K. Webb ESE 470
63 K. Webb ESE 470
64 K. Webb ESE 470
 Perform three iterations toward the solution of the following system
of equations using the Gauss-Seidel method. Let 𝐱𝐱 0 = 0, 0 𝑇𝑇 .
2𝑥𝑥1 + 𝑥𝑥2 = 12
2𝑥𝑥1 + 3𝑥𝑥2 = 5

65 K. Webb ESE 470

66 K. Webb ESE 470
67 K. Webb ESE 470
68 K. Webb ESE 470
 Perform three iterations toward the solution of the following
equation using the Newton-Raphson method. Let 𝐱𝐱 0 = 0.
𝑓𝑓 𝑥𝑥 = cos 𝑥𝑥 + 3𝑥𝑥 = 10

69 K. Webb ESE 470

70 K. Webb ESE 470
71 Nonlinear Systems of Equations

K. Webb ESE 470

 Nonlinear Systems of Equations
72

 Now, consider a system of nonlinear equations

 Can be represented as a vector of 𝑁𝑁 functions
 Each is a function of an 𝑁𝑁-vector of unknown variables
𝑦𝑦1 𝑓𝑓1 𝑥𝑥1 , 𝑥𝑥2 , ⋯ , 𝑥𝑥𝑁𝑁
𝑦𝑦2 𝑓𝑓 𝑥𝑥 , 𝑥𝑥 , ⋯ , 𝑥𝑥𝑁𝑁
𝐲𝐲 = ⋮ = 𝐟𝐟 𝐱𝐱 = 2 1 2
⋮
𝑦𝑦𝑁𝑁 𝑓𝑓𝑁𝑁 𝑥𝑥1 , 𝑥𝑥2 , ⋯ , 𝑥𝑥𝑁𝑁

 We can again approximate this function using a first-order Taylor series

𝐲𝐲 = 𝐟𝐟 𝐱𝐱 ≈ 𝐟𝐟 𝐱𝐱 0 + 𝐟𝐟 ′ 𝐱𝐱 0 𝐱𝐱 − 𝐱𝐱 0 (19)
 Note that all variables are 𝑁𝑁-vectors
 𝐟𝐟 is an 𝑁𝑁-vector of known, nonlinear functions
 𝐱𝐱 is an 𝑁𝑁-vector of unknown values – this is what we want to solve for
 𝐲𝐲 is an 𝑁𝑁-vector of known values
 𝐱𝐱𝟎𝟎 is an 𝑁𝑁-vector of 𝐱𝐱 values for which 𝐟𝐟 𝐱𝐱0 is known

K. Webb ESE 470

 Newton-Raphson Method
73

 Equation (19) is the basis for our Newton-Raphson iterative formula

 Again, let it be an equality and solve for 𝐱𝐱
𝐲𝐲 − 𝐟𝐟 𝐱𝐱 0 = 𝐟𝐟 ′ 𝐱𝐱 0 𝐱𝐱 − 𝐱𝐱 0
𝐟𝐟 ′ 𝐱𝐱 0 −𝟏𝟏
𝐲𝐲 − 𝐟𝐟 𝐱𝐱 0 = 𝐱𝐱 − 𝐱𝐱 0
𝐱𝐱 = 𝐱𝐱 0 + 𝐟𝐟 ′ 𝐱𝐱 0 −𝟏𝟏
𝐲𝐲 − 𝐟𝐟 𝐱𝐱 0

 This last expression can be used as an iterative formula

𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + 𝐟𝐟 ′ 𝐱𝐱 𝑖𝑖 −𝟏𝟏
𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖
 The derivative term on the right-hand side of (20) is an 𝑁𝑁 × 𝑁𝑁 matrix
 The Jacobian matrix, 𝐉𝐉
𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + 𝐉𝐉𝑖𝑖−𝟏𝟏 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖 (20)

K. Webb ESE 470

 The Jacobian Matrix
74

𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + 𝐉𝐉𝑖𝑖−𝟏𝟏 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖 (20)

 Jacobian matrix
 𝑁𝑁 × 𝑁𝑁 matrix of partial derivatives for 𝐟𝐟 𝐱𝐱
 Evaluated at the current value of 𝐱𝐱, 𝐱𝐱 𝑖𝑖

𝜕𝜕𝑓𝑓1 𝜕𝜕𝑓𝑓1 𝜕𝜕𝑓𝑓1

⋯
𝜕𝜕𝑥𝑥1 𝜕𝜕𝑥𝑥2 𝜕𝜕𝑥𝑥𝑁𝑁
𝜕𝜕𝑓𝑓2 𝜕𝜕𝑓𝑓2 𝜕𝜕𝑓𝑓2
𝐉𝐉𝑖𝑖 = 𝜕𝜕𝑥𝑥1 ⋯
𝜕𝜕𝑥𝑥2 𝜕𝜕𝑥𝑥𝑁𝑁
⋮ ⋮ ⋱ ⋮
𝜕𝜕𝑓𝑓𝑁𝑁 𝜕𝜕𝑓𝑓𝑁𝑁 𝜕𝜕𝑓𝑓𝑁𝑁
⋯
𝜕𝜕𝑥𝑥1 𝜕𝜕𝑥𝑥2 𝜕𝜕𝑥𝑥𝑁𝑁 𝐱𝐱=𝐱𝐱 𝑖𝑖
K. Webb ESE 470
 Newton-Raphson Method
75

𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + 𝐉𝐉𝑖𝑖−𝟏𝟏 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖 (20)

 We could iterate (20) until convergence or a maximum

number of iterations is reached
 Requires inversion of the Jacobian matrix
 Computationally expensive and error prone
 Instead, go back to the Taylor series approximation
𝐲𝐲 = 𝐟𝐟 𝐱𝐱 𝑖𝑖 + 𝐉𝐉𝑖𝑖 𝐱𝐱 𝑖𝑖+1 − 𝐱𝐱 𝑖𝑖
𝐲𝐲 − 𝐟𝐟 𝐱𝐱𝑖𝑖 = 𝐉𝐉𝑖𝑖 𝐱𝐱 𝑖𝑖+1 − 𝐱𝐱 𝑖𝑖 (21)
 Left side of (21) represents a difference between the known and
approximated outputs
 Right side represents an increment of the approximation for 𝐱𝐱
Δ𝐲𝐲𝑖𝑖 = 𝐉𝐉𝑖𝑖 Δ𝐱𝐱 𝑖𝑖 (22)
K. Webb ESE 470
 Newton-Raphson Method
76

Δ𝐲𝐲𝑖𝑖 = 𝐉𝐉𝑖𝑖 Δ𝐱𝐱 𝑖𝑖 (22)

 On each iteration:
 Compute Δ𝐲𝐲𝑖𝑖 and 𝐉𝐉𝑖𝑖
 Solve for Δ𝐱𝐱 𝑖𝑖 using Gaussian elimination
 Matrix
inversion not required
 Computationally robust

 Update 𝐱𝐱
(23)
𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + Δ𝐱𝐱 𝑖𝑖

K. Webb ESE 470

 Newton-Raphson – Example
77

 Apply Newton-Raphson to solve the following system of

nonlinear equations
𝐟𝐟 𝐱𝐱 = 𝐲𝐲
𝑥𝑥12 + 3𝑥𝑥2 21
=
𝑥𝑥1 𝑥𝑥2 12
 Initial condition: 𝐱𝐱 0 = 1 2 𝑇𝑇
 Stopping criterion: 𝜀𝜀𝑠𝑠 = 1 × 10−6
 Jacobian matrix
𝜕𝜕𝑓𝑓1 𝜕𝜕𝑓𝑓1
𝜕𝜕𝑥𝑥1 𝜕𝜕𝑥𝑥2 2𝑥𝑥 3
𝐉𝐉𝑖𝑖 = = 𝑥𝑥 1,𝑖𝑖 𝑥𝑥1,𝑖𝑖
𝜕𝜕𝑓𝑓2 𝜕𝜕𝑓𝑓2 2,𝑖𝑖
𝜕𝜕𝑥𝑥1 𝜕𝜕𝑥𝑥2 𝐱𝐱=𝐱𝐱 𝑖𝑖
K. Webb ESE 470
 Newton-Raphson – Example
78

Δ𝐲𝐲𝑖𝑖 = 𝐉𝐉𝑖𝑖 Δ𝐱𝐱 𝑖𝑖 (22)

𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + Δ𝐱𝐱 𝑖𝑖 (23)

 Adjusting the indexing, we can equivalently write (22)

and (23) as:
Δ𝐲𝐲𝑖𝑖−1 = 𝐉𝐉𝑖𝑖−1 Δ𝐱𝐱 𝑖𝑖−1 (22)
𝐱𝐱 𝑖𝑖 = 𝐱𝐱 𝑖𝑖−1 + Δ𝐱𝐱 𝑖𝑖−1 (23)

 For iteration 𝑖𝑖:

 Compute Δ𝐲𝐲𝑖𝑖−1 and 𝐉𝐉𝑖𝑖−1
 Solve (22) for Δ𝐱𝐱 𝑖𝑖−1
 Update 𝐱𝐱 using (23)

K. Webb ESE 470

 Newton-Raphson – Example
79

 𝑖𝑖 = 1:
21 7 14
Δ𝑦𝑦0 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 0 = − =
12 2 10
2𝑥𝑥1,0 3 2 3
𝐉𝐉0 = 𝑥𝑥 𝑥𝑥1,0 =
2,0 2 1
4
Δ𝐱𝐱 0 =
2
1 4 5
𝐱𝐱1 = 𝐱𝐱 0 + Δ𝐱𝐱 0 = + =
2 2 4
𝑥𝑥𝑘𝑘,1 − 𝑥𝑥𝑘𝑘,0
𝜀𝜀1 = max , 𝑘𝑘 = 1 … 𝑁𝑁
𝑥𝑥𝑘𝑘,0

5
𝑥𝑥1 = , 𝜀𝜀1 = 4
4

K. Webb ESE 470

 Newton-Raphson – Example
80

 𝑖𝑖 = 2:
21 37 −16
Δ𝑦𝑦1 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱1 = − =
12 20 −8
2𝑥𝑥1,1 3 10 3
𝐉𝐉1 = 𝑥𝑥 𝑥𝑥1,1 =
2,1 4 5
−1.474
Δ𝐱𝐱1 =
−0.421
5 −1.474 3.526
𝐱𝐱 2 = 𝐱𝐱1 + Δ𝐱𝐱1 = + =
4 −0.421 3.579
𝑥𝑥𝑘𝑘,2 − 𝑥𝑥𝑘𝑘,1
𝜀𝜀2 = max , 𝑘𝑘 = 1 … 𝑁𝑁
𝑥𝑥𝑘𝑘,1

3.526
𝑥𝑥2 = , 𝜀𝜀2 = 0.295
3.579

K. Webb ESE 470

 Newton-Raphson – Example
81

 𝑖𝑖 = 3:
21 23.172 −2.172
Δ𝑦𝑦2 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 2 = − =
12 12.621 −0.621
2𝑥𝑥1,2 3 7.053 3
𝐉𝐉2 = 𝑥𝑥 𝑥𝑥1,2 =
2,2 3.579 3.526
−0.410
Δ𝐱𝐱 2 =
0.240
3.526 −0.410 3.116
𝐱𝐱 3 = 𝐱𝐱 2 + Δ𝐱𝐱 2 = + =
3.579 0.240 3.819
𝑥𝑥𝑘𝑘,3 − 𝑥𝑥𝑘𝑘,2
𝜀𝜀3 = max , 𝑘𝑘 = 1 … 𝑁𝑁
𝑥𝑥𝑘𝑘,2

3.116
𝑥𝑥3 = , 𝜀𝜀3 = 0.116
3.819

K. Webb ESE 470

 Newton-Raphson – Example
82

 𝑖𝑖 = 7:
21 21.000 −7
Δ𝑦𝑦6 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 6 = − = −0.527 × 10
12 12.000 0.926 × 10−7
2𝑥𝑥 3 6.000 3
𝐉𝐉6 = 𝑥𝑥 1,6 𝑥𝑥1,6 =
2,6 4.000 3.000
−6
Δ𝐱𝐱 6 = −0.073 × 10
0.128 × 10−6
3.000 −6 3.000
𝐱𝐱 7 = 𝐱𝐱 6 + Δ𝐱𝐱 6 = + −0.073 × 10 =
4.000 0.128 × 10−6 4.000
𝑥𝑥𝑘𝑘,7 − 𝑥𝑥𝑘𝑘,6
𝜀𝜀7 = max , 𝑘𝑘 = 1 … 𝑁𝑁
𝑥𝑥𝑘𝑘,6

3.000
𝑥𝑥7 = , 𝜀𝜀7 = 31.9 × 10−9
4.000
K. Webb ESE 470
 Newton-Raphson – MATLAB Code
83

 Define the system of equations

 Initialize 𝐱𝐱

 Set up solution parameters

K. Webb ESE 470

 Newton-Raphson – MATLAB Code
84

 Iterate:
 Compute Δ𝐲𝐲𝑖𝑖−1 and 𝐉𝐉𝑖𝑖−1
 Solve for Δ𝐱𝐱 𝑖𝑖−1
 Update 𝐱𝐱

K. Webb ESE 470

85 Example Problems

K. Webb ESE 470

 Perform three iterations toward the solution of the following system
of equations using the Newton-Raphson method. Let 𝐱𝐱 0 = 1, 1 𝑇𝑇 .
10𝑥𝑥12 + 𝑥𝑥2 = 20
𝑒𝑒 𝑥𝑥1 − 𝑥𝑥2 = 10

86 K. Webb ESE 470

87 K. Webb ESE 470
88 K. Webb ESE 470
89 K. Webb ESE 470
90 Power-Flow Solution – Overview

K. Webb ESE 470

 Solving the Power-Flow Problem - Overview
91

 Consider an 𝑁𝑁-bus power-flow problem

 1 slack bus
 𝑛𝑛𝑃𝑃𝑃𝑃 PV buses
 𝑛𝑛𝑃𝑃𝑃𝑃 PQ buses
𝑁𝑁 = 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃 + 1

 Each bus has two unknown quantities

 Two of 𝑉𝑉𝑘𝑘 , 𝛿𝛿𝑘𝑘 , 𝑃𝑃𝑘𝑘 , and 𝑄𝑄𝑘𝑘
 For the N-R power-flow problem, 𝑉𝑉𝑘𝑘 and 𝛿𝛿𝑘𝑘 are the
unknown quantities
 These are the inputs to the nonlinear system of equations – the
𝑃𝑃𝑘𝑘 and 𝑄𝑄𝑘𝑘 equations – of (9) and (10)
 Finding unknown 𝑉𝑉𝑘𝑘 and 𝛿𝛿𝑘𝑘 values allows us to determine
unknown 𝑃𝑃𝑘𝑘 and 𝑄𝑄𝑘𝑘 values

K. Webb ESE 470

 Solving the Power-Flow Problem - Overview
92

 The nonlinear system of equations is

𝐲𝐲 = 𝐟𝐟(𝐱𝐱)
 The unknowns , 𝐱𝐱, are bus voltages
 Unknown phase angles from PV and PQ buses
 Unknown magnitudes from PQ bus

𝛿𝛿2
⋮ 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃
𝛿𝛿𝑛𝑛𝑃𝑃𝑃𝑃 +𝑛𝑛𝑃𝑃𝑃𝑃+1
𝛅𝛅
𝐱𝐱 = = (24)
𝐕𝐕 𝑉𝑉𝑛𝑛𝑃𝑃𝑃𝑃 +2
𝑛𝑛𝑃𝑃𝑃𝑃
⋮
𝑉𝑉𝑛𝑛𝑃𝑃𝑃𝑃 +𝑛𝑛𝑃𝑃𝑃𝑃+1
K. Webb ESE 470
 Solving the Power-Flow Problem - Overview
93

𝐲𝐲 = 𝐟𝐟(𝐱𝐱)
 The knowns , 𝐲𝐲, are bus powers
 Known real power from PV and PQ buses
 Known reactive power from PQ bus

𝑃𝑃2
⋮ 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃
𝑃𝑃𝑛𝑛𝑃𝑃𝑃𝑃 +𝑛𝑛𝑃𝑃𝑃𝑃 +1
𝐏𝐏
𝐲𝐲 = = (25)
𝐐𝐐 𝑄𝑄𝑛𝑛𝑃𝑃𝑃𝑃 +2
𝑛𝑛𝑃𝑃𝑃𝑃
⋮
𝑄𝑄𝑛𝑛𝑃𝑃𝑃𝑃 +𝑛𝑛𝑃𝑃𝑃𝑃 +1
K. Webb ESE 470
 Solving the Power-Flow Problem - Overview
94

𝐲𝐲 = 𝐟𝐟(𝐱𝐱)
 The system of equations , 𝐟𝐟, consists of the
nonlinear functions for 𝐏𝐏 and 𝐐𝐐
 Nonlinear functions of 𝐕𝐕 and 𝛅𝛅
𝑃𝑃2 𝐱𝐱
⋮ 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃
𝐏𝐏(𝐱𝐱) ⋮
𝐟𝐟(𝐱𝐱) = = (26)
𝐐𝐐(𝐱𝐱)
𝑄𝑄𝑛𝑛𝑃𝑃𝑃𝑃 +2 𝐱𝐱
𝑛𝑛𝑃𝑃𝑃𝑃
⋮
⋮
K. Webb ESE 470
 Solving the Power-Flow Problem - Overview
95

 𝑃𝑃𝑘𝑘 𝐱𝐱 and 𝑄𝑄𝑘𝑘 𝐱𝐱 are given by

𝑁𝑁

𝑃𝑃𝑘𝑘 = 𝑉𝑉𝑘𝑘 � 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (9)

𝑛𝑛=1

𝑁𝑁

𝑄𝑄𝑘𝑘 = 𝑉𝑉𝑘𝑘 � 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (10)

𝑛𝑛=1

 Admittance matrix terms are

𝑌𝑌𝑘𝑘𝑘𝑘 = 𝑌𝑌𝑘𝑘𝑘𝑘 ∠𝜃𝜃𝑘𝑘𝑘𝑘

K. Webb ESE 470

 Solving the Power-Flow Problem - Overview
96

 The iterative N-R formula is

𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + Δ𝐱𝐱 𝑖𝑖
 The increment term, Δ𝐱𝐱 𝑖𝑖 , is computed through
Gaussian elimination of
Δ𝐲𝐲𝑖𝑖 = 𝐉𝐉𝑖𝑖 Δ𝐱𝐱 𝑖𝑖
 The Jacobian, 𝐽𝐽𝑖𝑖 , is computed on each iteration
 The power mismatch vector is

Δ𝐲𝐲𝑖𝑖 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖
 𝐲𝐲 is the vector of known powers, as given in (25)
 𝐟𝐟 𝐱𝐱 𝑖𝑖 are the 𝑃𝑃 and 𝑄𝑄 equations given by (9) and (10)
K. Webb ESE 470
97 Power-Flow Solution – Procedure

K. Webb ESE 470

 Solving the Power-Flow Problem - Procedure
98

 The following procedure shows how to set up and

solve the power-flow problem using the N-R
algorithm
1. Order and number buses
 Slack bus is #1

 Group all PV buses together next

 Group all PQ buses together last

2. Generate the bus admittance matrix, 𝐘𝐘

 And magnitude, Y = 𝐘𝐘 , and angle, θ = ∠𝐘𝐘, matrices

K. Webb ESE 470

 Solving the Power-Flow Problem - Procedure
99

3. Initialize known quantities

 Slack bus: 𝑉𝑉1 and 𝛿𝛿1
 PV buses: 𝑉𝑉𝑘𝑘 and 𝑃𝑃𝑘𝑘
 PQ buses: 𝑃𝑃𝑘𝑘 and 𝑄𝑄𝑘𝑘
 Output vector:

𝐏𝐏
𝐲𝐲 =
𝐐𝐐
4. Initialize unknown quantities
0
⋮ 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃
𝛅𝛅𝟎𝟎 0
𝐱𝐱 𝒐𝒐 = = (24)
𝐕𝐕𝟎𝟎
1.0
𝑛𝑛𝑃𝑃𝑃𝑃
⋮
1.0

K. Webb ESE 470

 Solving the Power-Flow Problem - Procedure
100

5. Set up Newton-Raphson parameters

 Tolerance for convergence, reltol
 Maximum # of iterations, max_iter
 Initialize relative error: 𝜀𝜀0 > reltol, e.g. 𝜀𝜀0 = 10
 Initialize iteration counter: 𝑖𝑖 = 0

6. while (𝜀𝜀 > reltol) && (𝑖𝑖 < max_iter)

 Update bus voltage phasor vector, 𝐕𝐕𝑖𝑖 , using magnitude
and phase values from 𝐱𝐱 𝑖𝑖 and from knowns
 Calculate the current injected into each bus, a vector of
phasors
𝐈𝐈𝑖𝑖 = 𝐘𝐘 ⋅ 𝐕𝐕𝑖𝑖
K. Webb ESE 470
 Solving the Power-Flow Problem - Procedure
101

6. while (𝜀𝜀 > reltol) && (𝑖𝑖 < max_iter) – cont’d

 Calculate complex, real, and reactive power injected into each bus
 This can be done using 𝐕𝐕𝑖𝑖 and 𝐈𝐈𝑖𝑖 vectors and element-by-element
multiplication (the .* operator in MATLAB)
∗
𝐒𝐒𝑘𝑘,𝑖𝑖 = 𝐕𝐕𝑘𝑘,𝑖𝑖 ⋅ 𝐈𝐈𝑘𝑘,𝑖𝑖
𝑃𝑃𝑘𝑘,𝑖𝑖 = 𝑅𝑅𝑅𝑅 𝐒𝐒𝑘𝑘,𝑖𝑖
𝑄𝑄𝑘𝑘,𝑖𝑖 = 𝐼𝐼𝐼𝐼 𝐒𝐒𝑘𝑘,𝑖𝑖
 Create 𝑓𝑓 𝑥𝑥𝑖𝑖 from 𝐏𝐏𝑖𝑖 and 𝐐𝐐𝑖𝑖 vectors
 Calculate power mismatch, Δ𝐲𝐲𝑖𝑖
Δ𝐲𝐲𝑖𝑖 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 𝑖𝑖
 Compute the Jacobian, 𝐉𝐉𝑖𝑖 , using voltage magnitudes and phase
angles from 𝐕𝐕𝑖𝑖
K. Webb ESE 470
 Solving the Power-Flow Problem - Procedure
102

6. while (𝜀𝜀 > reltol) && (𝑖𝑖 < max_iter) – cont’d

 Solve for Δ𝐱𝐱 𝑖𝑖 using Gaussian elimination
Δ𝐲𝐲𝑖𝑖 = 𝐉𝐉𝑖𝑖 Δ𝐱𝐱 𝑖𝑖
 Use the mldivide (\, backslash) operator in MATLAB: Δ𝐱𝐱 𝑖𝑖 = 𝐉𝐉𝑖𝑖 \Δ𝐲𝐲𝑖𝑖
 Update 𝐱𝐱
𝐱𝐱 𝑖𝑖+1 = 𝐱𝐱 𝑖𝑖 + Δ𝐱𝐱 𝑖𝑖
 Check for convergence using power mismatch
𝑦𝑦𝑘𝑘 − 𝑓𝑓𝑘𝑘 𝐱𝐱
𝜀𝜀𝑖𝑖+1 = max
𝑦𝑦𝑘𝑘
 Update the number of iterations
𝑖𝑖 = 𝑖𝑖 + 1
K. Webb ESE 470
 The Jacobian Matrix
103

 The Jacobian matrix has four quadrants of varying

dimension depending on the number of different
types of buses:
𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃 𝑛𝑛𝑃𝑃𝑃𝑃

𝜕𝜕𝐏𝐏 𝜕𝜕𝐏𝐏
𝜕𝜕𝛅𝛅 𝜕𝜕𝐕𝐕 𝑛𝑛𝑃𝑃𝑃𝑃 + 𝑛𝑛𝑃𝑃𝑃𝑃
𝐉𝐉1 𝐉𝐉2
𝐉𝐉 =
𝜕𝜕𝐐𝐐 𝜕𝜕𝐐𝐐
𝑛𝑛𝑃𝑃𝑃𝑃
𝜕𝜕𝛅𝛅 𝜕𝜕𝐕𝐕
𝐉𝐉3 𝐉𝐉4
K. Webb ESE 470
 The Jacobian Matrix
104

 Jacobian elements are partial derivatives of (9) and (10) with

respect to 𝛿𝛿 or 𝑉𝑉
 Formulas for the Jacobian elements:
 𝑛𝑛 ≠ 𝑘𝑘
𝜕𝜕𝑃𝑃𝑘𝑘
𝐉𝐉1𝑘𝑘𝑘𝑘 = = 𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (27)
𝜕𝜕𝛿𝛿𝑛𝑛

𝜕𝜕𝑃𝑃𝑘𝑘
𝐉𝐉2𝑘𝑘𝑘𝑘 = = 𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (28)
𝜕𝜕𝑉𝑉𝑛𝑛

𝜕𝜕𝑄𝑄𝑘𝑘
𝐉𝐉3𝑘𝑘𝑘𝑘 = = −𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (29)
𝜕𝜕𝛿𝛿𝑛𝑛

𝜕𝜕𝑄𝑄𝑘𝑘
𝐉𝐉4𝑘𝑘𝑘𝑘 = = 𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (30)
𝜕𝜕𝑉𝑉𝑛𝑛
K. Webb ESE 470
 The Jacobian Matrix
105

 Formulas for the Jacobian elements, cont’d:

 𝑛𝑛 = 𝑘𝑘
𝑁𝑁
𝜕𝜕𝑃𝑃𝑘𝑘
𝐉𝐉1𝑘𝑘𝑘𝑘 = = −𝑉𝑉𝑘𝑘 � 𝑌𝑌𝑘𝑘𝑛𝑛 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (31)
𝜕𝜕𝛿𝛿𝑘𝑘
𝑛𝑛=1
𝑛𝑛≠𝑘𝑘

𝑁𝑁
𝜕𝜕𝑃𝑃𝑘𝑘
𝐉𝐉2𝑘𝑘𝑘𝑘 = = 𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 cos 𝜃𝜃𝑘𝑘𝑘𝑘 + � 𝑌𝑌𝑘𝑘𝑛𝑛 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (32)
𝜕𝜕𝑉𝑉𝑘𝑘
𝑛𝑛=1

𝑁𝑁
𝜕𝜕𝑄𝑄𝑘𝑘
𝐉𝐉3𝑘𝑘𝑘𝑘 = = 𝑉𝑉𝑘𝑘 � 𝑌𝑌𝑘𝑘𝑛𝑛 𝑉𝑉𝑛𝑛 cos 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (33)
𝜕𝜕𝛿𝛿𝑘𝑘
𝑛𝑛=1
𝑛𝑛≠𝑘𝑘

𝑁𝑁
𝜕𝜕𝑄𝑄𝑘𝑘
𝐉𝐉4𝑘𝑘𝑘𝑘 = = −𝑉𝑉𝑘𝑘 𝑌𝑌𝑘𝑘𝑘𝑘 sin 𝜃𝜃𝑘𝑘𝑘𝑘 + � 𝑌𝑌𝑘𝑘𝑘𝑘 𝑉𝑉𝑛𝑛 sin 𝛿𝛿𝑘𝑘 − 𝛿𝛿𝑛𝑛 − 𝜃𝜃𝑘𝑘𝑘𝑘 (34)
𝜕𝜕𝑉𝑉𝑘𝑘
𝑛𝑛=1

K. Webb ESE 470

106 Power-Flow Solution – Example

K. Webb ESE 470

 Power-Flow Solution – Buses
107

 Determine all bus voltages and power flows for the following three-
bus power system

 Three buses, 𝑛𝑛𝑃𝑃𝑃𝑃 = 1, 𝑛𝑛𝑃𝑃𝑃𝑃 = 1, ordered PV first, then PQ:

 Bus 1: slack bus
 𝑉𝑉1 and 𝛿𝛿1 are known, find 𝑃𝑃1 and 𝑄𝑄1

 Bus 2: PV bus
 𝑃𝑃2 and 𝑉𝑉2 are known, find 𝛿𝛿2 and 𝑄𝑄2

 Bus 3: PQ bus
 𝑃𝑃3 and 𝑄𝑄3 are known, find 𝑉𝑉3 and 𝛿𝛿3
K. Webb ESE 470
 Power-Flow Solution – Admittance Matrix
108

 Per-unit, per-length impedance of all transmission lines:

𝑧𝑧 = 31.1 + 𝑗𝑗𝑗𝑗𝑗 × 10−6 𝑝𝑝𝑝𝑝/𝑘𝑘𝑘𝑘
 Admittance of each line:
1
𝑌𝑌12 = 𝑌𝑌23 = = 2.06 − 𝑗𝑗𝑗𝑗.9 𝑝𝑝𝑝𝑝
𝑧𝑧 ⋅ 150 𝑘𝑘𝑘𝑘
1
𝑌𝑌13 = = 1.54 − 𝑗𝑗15.7 𝑝𝑝𝑝𝑝
𝑧𝑧 ⋅ 200 𝑘𝑘𝑘𝑘
K. Webb ESE 470
 Power-Flow Solution – Admittance Matrix
109

 The admittance matrix (see p. 8):

𝑌𝑌11 𝑌𝑌12 𝑌𝑌13 3.6 − 𝑗𝑗𝑗𝑗.6 −2.06 + 𝑗𝑗𝑗𝑗.9 −1.5 + 𝑗𝑗𝑗𝑗.7
𝐘𝐘 = 𝑌𝑌21 𝑌𝑌22 𝑌𝑌23 = −2.06 + 𝑗𝑗𝑗𝑗.9 4.1 − 𝑗𝑗𝑗𝑗.8 −2.06 + 𝑗𝑗𝑗𝑗.9
𝑌𝑌31 𝑌𝑌32 𝑌𝑌33 −1.5 + 𝑗𝑗𝑗𝑗.7 −2.06 + 𝑗𝑗𝑗𝑗.9 3.6 − 𝑗𝑗𝑗𝑗.6

 Admittance magnitude and angle matrices:

36.8 21.0 15.8 −84.4° 95.6° 95.6°
Y = 𝐘𝐘 = 21.0 42.0 21.0 , 𝛉𝛉 = 95.6° −84.4° 95.6°
15.8 21.0 36.8 95.6° 95.6° −84.4°

K. Webb ESE 470

 Power-Flow Solution – Initialize Knowns
110

 Known quantities
 Slack bus: 𝑉𝑉1 = 1.0 𝑝𝑝𝑝𝑝, 𝛿𝛿1 = 0°
 PV bus: 𝑉𝑉2 = 1.05 𝑝𝑝𝑝𝑝, 𝑃𝑃2 = 2.0 𝑝𝑝𝑝𝑝
 PQ bus: 𝑃𝑃3 = −5.0 𝑝𝑝𝑝𝑝, 𝑄𝑄3 = −1.0 𝑝𝑝𝑝𝑝
 Output vector
𝑃𝑃2 2.0
𝐏𝐏
𝐲𝐲 = = 𝑃𝑃3 = −5.0
𝐐𝐐
𝑄𝑄3 −1.0

K. Webb ESE 470

 Power-Flow Solution – Initialize Unknowns
111

 The vector of unknown quantities to be solved for is

𝛿𝛿2
𝛅𝛅
𝐱𝐱 = = 𝛿𝛿3
𝐕𝐕
𝑉𝑉3

 Initialize all unknown bus voltage phasors to 𝐕𝐕𝑘𝑘 = 1.0∠0° 𝑝𝑝𝑝𝑝

𝛿𝛿2,0 0
𝛅𝛅0
𝐱𝐱 0 = = 𝛿𝛿3,0 = 0
𝐕𝐕0
𝑉𝑉3,0 1.0

 The complete vector of bus voltage phasors – partly known, partly

unknown – is
𝑉𝑉1 ∠𝛿𝛿1 1.0∠0° 1.0∠0°
𝐕𝐕 = 𝑉𝑉2 ∠𝛿𝛿2 = 1.05∠𝛿𝛿2,0 = 1.05∠0°
𝑉𝑉3 ∠𝛿𝛿3 𝑉𝑉3,0 ∠𝛿𝛿3,0 1.0∠0°
K. Webb ESE 470
 Power-Flow Solution – Jacobian Matrix
112

 The Jacobian matrix for this system is

𝜕𝜕𝑃𝑃2 𝜕𝜕𝑃𝑃2 𝜕𝜕𝑃𝑃2
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3
𝜕𝜕𝑃𝑃3 𝜕𝜕𝑃𝑃3 𝜕𝜕𝑃𝑃3
𝐉𝐉 =
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3
𝜕𝜕𝑄𝑄3 𝜕𝜕𝑄𝑄3 𝜕𝜕𝑄𝑄3
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3

 This matrix will be computed on each iteration using

the current approximation to the vector of
unknowns, 𝐱𝐱 𝑖𝑖
K. Webb ESE 470
 Power-Flow Solution – Set Up and Iterate
113

 Set up N-R iteration parameters

 reltol = 1e-6
 max_iter = 1e3
 𝜀𝜀0 = 10
 𝑖𝑖 = 0

 Iteratively update the approximation to the vector of

unknowns as long as
 Stopping criterion is not satisfied
𝜀𝜀𝑖𝑖 > 𝜀𝜀𝑠𝑠

 Maximum number of iterations is not exceeded

𝑖𝑖 ≤ 𝑚𝑚𝑚𝑚𝑚𝑚 _𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖

K. Webb ESE 470

 Power-Flow Solution – Iterate
114

 𝑖𝑖 = 0:
 Vector of bus voltage phasors

𝑉𝑉1 ∠𝛿𝛿1 1.0∠0°

𝐕𝐕0 = 𝑉𝑉2 ∠𝛿𝛿2,0 = 1.05∠0°
𝑉𝑉3,0 ∠𝛿𝛿3,0 1.0∠0°

 Current injected into each bus

𝐈𝐈0 = 𝐘𝐘 ⋅ 𝐕𝐕0

3.6 − 𝑗𝑗𝑗𝑗.6 −2.1 + 𝑗𝑗𝑗𝑗.9 −1.5 + 𝑗𝑗𝑗𝑗.7 1.0∠0°

𝐈𝐈0 = −2.1 + 𝑗𝑗𝑗𝑗.9 4.1 − 𝑗𝑗𝑗𝑗.8 −2.1 + 𝑗𝑗𝑗𝑗.9 1.05∠0°
−1.5 + 𝑗𝑗𝑗𝑗.7 −2.1 + 𝑗𝑗𝑗𝑗.9 3.6 − 𝑗𝑗𝑗𝑗.6 1.0∠0°

1.05∠95.6°
𝐈𝐈0 = 2.10∠ − 84.4°
1.05∠95.6°

K. Webb ESE 470

 Power-Flow Solution – Iterate
115

 𝑖𝑖 = 0:
 Complex power injected into each bus

𝐒𝐒0 = 𝐕𝐕0 .∗ 𝐈𝐈0∗

∗
1.0∠0° 1.05∠95.6°
𝐒𝐒0 = 1.05∠0° .∗ 2.10∠ − 84.4°
1.0∠0° 1.05∠95.6°
−0.103 − 𝑗𝑗𝑗.045
𝐒𝐒0 = 0.216 + 𝑗𝑗𝑗.195
−0.103 − 𝑗𝑗𝑗.045
 Real and reactive power
−0.103 −1.045
𝐏𝐏0 = 0.216 , 𝐐𝐐0 = 2.195
−0.103 −1.045

K. Webb ESE 470

 Power-Flow Solution – Iterate
116

 𝑖𝑖 = 0:
 Power mismatch

Δ𝐲𝐲0 = 𝐲𝐲 − 𝐟𝐟 𝐱𝐱 0
2.0 0.216 1.784
Δ𝐲𝐲0 = −5.0 − −0.103 = −4.897
−1.0 −1.045 0.045
 Next, compute the Jacobian matrix
𝜕𝜕𝑃𝑃2 𝜕𝜕𝑃𝑃2 𝜕𝜕𝑃𝑃2
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3
𝜕𝜕𝑃𝑃3 𝜕𝜕𝑃𝑃3 𝜕𝜕𝑃𝑃3
𝐉𝐉0 =
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3
𝜕𝜕𝑄𝑄3 𝜕𝜕𝑄𝑄3 𝜕𝜕𝑄𝑄3
𝜕𝜕𝛿𝛿2 𝜕𝜕𝛿𝛿3 𝜕𝜕𝑉𝑉3 𝐱𝐱=𝐱𝐱 0

K. Webb ESE 470

 Power-Flow Solution – Iterate
117

 𝑖𝑖 = 0:
 Elements of the Jacobian matrix are computed using 𝑉𝑉
and 𝛿𝛿 values from 𝐕𝐕0 and 𝑌𝑌 and 𝜃𝜃 values from 𝐘𝐘:
1.0
𝑉𝑉0 = 1.05
1.0
0°
𝛿𝛿0 = 0°
0°
36.8 21.0 15.8
𝑌𝑌 = 21.0 42.0 21.0
15.8 21.0 36.8
−84.4° 95.6° 95.6°
𝜃𝜃 = 95.6° −84.4° 95.6°
95.6° 95.6° −84.4°
K. Webb ESE 470
 Power-Flow Solution – Iterate
118

 𝑖𝑖 = 0:
 Jacobian, 𝐉𝐉1
𝜕𝜕𝑃𝑃2
= −𝑉𝑉2 𝑌𝑌21 𝑉𝑉1 sin 𝛿𝛿2 − 𝛿𝛿1 − 𝜃𝜃21 + 𝑌𝑌23 𝑉𝑉3 sin 𝛿𝛿2 − 𝛿𝛿3 − 𝜃𝜃23
𝜕𝜕𝛿𝛿2
𝜕𝜕𝑃𝑃3
= −𝑉𝑉3 𝑌𝑌31 𝑉𝑉1 sin 𝛿𝛿3 − 𝛿𝛿1 − 𝜃𝜃31 + 𝑌𝑌32 𝑉𝑉2 sin 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32
𝜕𝜕𝛿𝛿3
𝜕𝜕𝑃𝑃2
= 𝑉𝑉2 𝑌𝑌23 𝑉𝑉3 sin 𝛿𝛿2 − 𝛿𝛿3 − 𝜃𝜃23
𝜕𝜕𝛿𝛿3
𝜕𝜕𝑃𝑃3
= 𝑉𝑉3 𝑌𝑌32 𝑉𝑉2 sin 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32
𝜕𝜕𝛿𝛿2

K. Webb ESE 470

 Power-Flow Solution – Iterate
119

 𝑖𝑖 = 0:
 Jacobian, 𝐉𝐉2
𝜕𝜕𝑃𝑃2
= 𝑉𝑉2 𝑌𝑌23 cos 𝛿𝛿2 − 𝛿𝛿3 − 𝜃𝜃23
𝜕𝜕𝑉𝑉3
𝜕𝜕𝑃𝑃3
= 2 ⋅ 𝑉𝑉3 𝑌𝑌33 cos 𝜃𝜃33 +
𝜕𝜕𝑉𝑉3
𝑌𝑌31 𝑉𝑉1 cos 𝛿𝛿3 − 𝛿𝛿1 − 𝜃𝜃31 + 𝑌𝑌32 𝑉𝑉2 cos 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32

 Jacobian, 𝐉𝐉3
𝜕𝜕𝑄𝑄3
= −𝑉𝑉3 𝑌𝑌32 𝑉𝑉2 cos 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32
𝜕𝜕𝛿𝛿2
𝜕𝜕𝑄𝑄3
= 𝑉𝑉3 𝑌𝑌31 𝑉𝑉1 cos 𝛿𝛿3 − 𝛿𝛿1 − 𝜃𝜃31 + 𝑌𝑌32 𝑉𝑉2 cos 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32
𝜕𝜕𝛿𝛿3

K. Webb ESE 470

 Power-Flow Solution – Iterate
120

 𝑖𝑖 = 0:
 Jacobian, 𝐉𝐉4
𝜕𝜕𝑄𝑄3
= 𝑉𝑉3 𝑌𝑌33 cos 𝜃𝜃33 +
𝜕𝜕𝑉𝑉3
𝑌𝑌31 𝑉𝑉1 cos 𝛿𝛿3 − 𝛿𝛿1 − 𝜃𝜃31 + 𝑌𝑌32 𝑉𝑉2 cos 𝛿𝛿3 − 𝛿𝛿2 − 𝜃𝜃32

 Evaluating the Jacobian expressions using 𝑉𝑉 and 𝛿𝛿

values from 𝐕𝐕0 and 𝑌𝑌 and 𝜃𝜃 values from 𝐘𝐘, gives
43.89 −21.95 −2.160
𝐉𝐉0 = −21.95 37.62 3.497
2.160 −3.702 35.53

K. Webb ESE 470

 Power-Flow Solution – Iterate
121

 𝑖𝑖 = 0:
 Use Gaussian elimination to solve for Δ𝐱𝐱 0

43.89 −21.95 −2.160 Δ𝑥𝑥1,0 1.784

Δ𝐲𝐲0 = 𝐉𝐉0 Δ𝐱𝐱 0 = −21.95 37.62 3.497 Δ𝑥𝑥2,0 = −4.897
2.160 −3.702 35.53 Δ𝑥𝑥3,0 0.045
−0.0345
Δ𝐱𝐱 0 = −0.1492
−0.0122

 Update the vector of unknowns, 𝐱𝐱

0 −0.0345 −0.0345
𝐱𝐱1 = 𝐱𝐱 0 + Δ𝐱𝐱 0 = 0 + −0.1492 = −0.1492
1.0 −0.0122 0.9878

K. Webb ESE 470

 Power-Flow Solution – Iterate
122

 𝑖𝑖 = 0:
 Use power mismatch to check for convergence
𝑦𝑦𝑘𝑘 − 𝑓𝑓𝑘𝑘 𝑥𝑥
𝜀𝜀0 = max = 0.9794
𝑦𝑦𝑘𝑘

 Move on to the next iteration, 𝑖𝑖 = 1

 Create 𝐕𝐕1 using 𝐱𝐱1 values
 Calculate 𝐈𝐈1
 Calculate 𝐒𝐒1 , 𝐏𝐏1 , 𝐐𝐐1
 Create 𝐟𝐟 𝐱𝐱1 from 𝐏𝐏1 and 𝐐𝐐1
 Calculate Δ𝐲𝐲1 , 𝐉𝐉1 , Δ𝐱𝐱1
 Update 𝐱𝐱 to 𝐱𝐱 2
 Check for convergence
…

K. Webb ESE 470

 Power-Flow Solution
123

 Convergence is achieved after four iterations

1.0∠0° 3.08 − 𝑗𝑗𝑗.82
𝐕𝐕4 = 1.1∠ − 2.1° , 𝐒𝐒4 = 2.0 + 𝑗𝑗𝑗.67
0.97∠ − 8.8° −5.0 − 𝑗𝑗𝑗.0

𝜀𝜀4 = 0.41 × 10−6

K. Webb ESE 470

124 Example Problems

K. Webb ESE 470

 For the power system shown, determine
a) The type of each bus
b) The first row of the admittance matrix, 𝐘𝐘
c) The vector of unknowns, 𝐱𝐱
d) The vector of knowns, 𝐲𝐲
e) The Jacobian matrix, 𝐉𝐉, in symbolic form

125 K. Webb ESE 470

126 K. Webb ESE 470
127 K. Webb ESE 470
128 K. Webb ESE 470