Analysis Techniques for Large-Scale Electrical Systems: Newton-Raphson Power Flow and Matrix Algebra

ECE 530 Analysis Techniques for Large-Scale Electrical Systems Lecture 4: Newton-Raphson Power Flow Prof. Hao Zhu Dept. of Electrical and Computer Engineering University of Illinois at Urbana-Champaign haozhu@illinois.edu 9/2/2015 1

Announcements No class next Monday (labor day) Homework 1 due next Wednesday, in class A reference book on linear algebra: The Matrix Cookbook http://www.math.uwaterloo.ca/~hwolkowi/matrixcookb ook.pdf 2

Newton-Raphson Method To solve a set of nonlinear equations denoted by f(x)= 0 N-R proceeds by iterative linearization Derivation of N-R method is similar to the scalar case ( ) ( ) ( ) ( ) 0 ( ) x J x f x = = + + + f x f x f x J x f x x higher order terms ( ) J x x 1 ( ) ( ) + ( 1) ( ) v ( ) v v = + x x x + ( 1) ( ) v ( ) v 1 ( ) v v = x x J x f x ( ) ( ) ( ) v f x Iterate until ( ) 3

Power Flow Variables Assume the slack bus is the first bus (with a fixed voltage angle/magnitude). We then need to determine the voltage angle/magnitude at the other buses. = + ( ) x P P P 2 2 2 2 G D + + x ( ) ( ) x P Gn P Q P Q n n Dn = x ( ) f x V Q 2 2 2 2 G D + ( ) x V Q Q Q n n Gn Dn 4

N-R Power Flow Solution The power flow is solved using the same procedure discussed with the general Newton-Raphson: ( ) v = x x Set 0; make an initial guess of , v ( ) v f x While ( ) Do + ( 1) ( ) v ( ) v 1 ( ) v v = = x x J x f x ( ) ( ) + 1 v v End While 5

Power Flow Formulation First consider the case that all buses are PQ-bus; ?? and ?? are specified; |??| and ?? unknown This specification of the two bus variables leads to two equality constraints (recall the real-valued power balance equations) ? ???= ?? ?? ???cos???+ ???sin??? ???+ ??? ?=1 ? ?= ?? ?? ?? ???sin??? ???cos??? ???+ ??? ?=1 6

Jacobian Matrix The most difficult part of the algorithm is determining and inverting the n by n Jacobian matrix, ( ) ( ) ( ) f f x x f f x x = J x J x x x x ( ) x f 1 1 1 1 x 2 x n x ( ) ( ) ( ) x f 2 2 2 ( ) 1 2 n x x x ( ) x ( ) x ( ) x f f f n n n 1 2 n 7

Power Flow Jacobian Elements of Jacobian matrix: differentiating each function with respect to (wrt) every variable Consider the real power equation at bus i ? ???= ?? ?? ???cos???+ ???sin??? ???+ ??? ?=1 The corresponding Jacobian entries are: ???? ??? ? ? = ?? ?? ???sin???+ ???cos??? ? ??? ???= ??|??|(???sin??? ???cos???) 8

Power Flow Jacobian, contd By differentiating it wrt voltage magnitude ? ???= ?? ?? ???cos???+ ???sin??? ???+ ??? ?=1 We have ? ??? ?|??|= 2 ?????+ ? ??? ???cos???+ ???sin??? ? ??? ? ??= ??(???cos???+ ???sin???) Similar derivations for reactive power equations 9

Two Bus Newton-Raphson Example For the two-bus power system shown below, use the Newton-Raphson power flow to determine the voltage magnitude and angle at bus two. Assume Slack bus 1, SBase = 100 MVA, and ? = 1? 3. Line Z = 0.1j One 1.000 pu Two 1.000 pu 0 MW 0 MVR 200 MW 100 MVR 10 10 j 10 10 j j j 2 = = x Y bus V 2 10

Two Bus Example, contd General power balance equations n = + = P ( cos sin ) V V G B P P i i k ik ik ik ik Gi Di = 1 k n = = Q ( sin cos ) V V G B Q Q i i k ik ik ik ik Gi Di = 1 k Bus two power balance equations (10sin V V ) 2.0 + = 0 2 1 2 2 ( 10cos + (10) 1.0 + = ) 0 V V V 2 1 2 2 11

Two Bus Example, contd = ) 2.0 + = x P ( ) (10sin 0 V 2 2 2 2 = ( 10cos + (10) 1.0 + = ( ) x ) 0 Q V V 2 2 2 2 Now calculate the power flow Jacobian P ( ) = = P ( ) x x 2 2 V 2 2 x 2 ( ) x J x Q ( ) Q ( ) 2 2 V 2 2 2 + 10 10 cos sin 10sin V V 2 2 2 20 10cos V 2 2 2 2 12

Two Bus Example, First Iteration 0 1 (0) = = x Set 0, guess v Calculate ) 2.0 + (10sin V 2.0 1.0 2 2 V (0) = = x f( ) 2 ( 10cos + + ) (10) 1.0 + V 2 2 2 10sin 10 10 cos sin V V 10 0 0 2 2 2 20 (0) = = J x ( ) 10cos V 10 2 2 2 2 12.0 1.0 0 1 10 0 0 0.2 0.9 (1) = = x Solve 10 13

Two Bus Example, Convergence 0.9(10sin( 0.2)) 2.0 + 0.212 0.279 (1) = = x f( ) 2 + 10 1.0 + 0.9( 10cos( 0.2)) 8.82 1.788 0.9 1.986 8.199 (1) = J x ( ) 1 0.2 0.9 = = 8.82 1.788 1.986 8.199 0.212 0.279 0.236 0.8554 0.233 0.8586 (2) = = x 0.0145 0.0190 0.0000906 0.0001175 (2) (3) = x x f( ) (3) = x f( ) Done! V 0.8554 13.52 2 14

Two Bus Solved Values Once the voltage angle and magnitude at bus 2 are known we can calculate all the other system values, such as the line flows and the generator reactive power output 200.0 MW 168.3 MVR -200.0 MW -100.0 MVR Line Z = 0.1j One 1.000 pu Two 0.855 pu -13.522 Deg 200.0 MW 168.3 MVR 200 MW 100 MVR PowerWorld Case Name: Bus2_Intro 15

Two Bus Case, Low Voltage Solution This case actually has two solutions! The second "low voltage" is found by using a low initial guess. = = 0 (0) x Set 0, guess v 0.25 Calculate ) 2.0 + (10sin V 2 2 2 V (0) = = x f( ) 2 0.875 ( 10cos + (10) 1.0 + + ) V 2 2 2 10sin 10 10 cos sin V V 2.5 0 0 2 2 2 20 (0) = = J x ( ) 10cos V 5 2 2 2 2 16

Low Voltage Solution, cont'd 1 0.8 0 2.5 0 0 2 (1) = = x Solve 0.875 0.25 1.462 0.534 5 0.075 0.921 0.220 1.42 0.2336 (2) (2) (3) = = = f x x x ( ) Low voltage solution 200.0 MW 831.7 MVR -200.0 MW -100.0 MVR Line Z = 0.1j One 1.000 pu Two 0.261 pu -49.914 Deg 200.0 MW 831.7 MVR 200 MW 100 MVR 17

Practical Power Flow Software Most commercial software packages have built in defaults to prevent convergence to low voltage solutions. One approach is to automatically change the load model from constant power to constant current or constant impedance when the load bus voltage gets too low In PowerWorld these defaults can be modified on the Tools, Simulator Options, Advanced Options page; note you also need to disable the Initialize from Flat Start Values option The PowerWorld case Bus2_Intro_Low is set solved to the low voltage solution Initial bus voltages can be set using the Bus Information Dialog 18

NR Initialization A textbook starting solution for the NR is to use flat start values in which all the angles are set to the slack bus angle, all the PQ bus voltages are set to 1.0, and all the PV bus voltages are set to their PV setpoint values This approach usually works for small systems It seldom works for large systems. The usual approach for solving large systems is to start with an existing solution, modify it, then hope it converges If not make the modification smaller More robust methods are possible, but convergence is certainly not guaranteed!! 19

PV Buses At PV bus i (generator bus), ??and |??| are fixed Similar to slack bus, no need to include these |??| s in x or write the reactive power balance equations the reactive power output of the generator varies to maintain the fixed terminal voltage (within limits) Otherwise, it can be included in x using an additional equality for PV-bus i ? ???= ?? ?? ???cos???+ ???sin??? ???+ ??? ?=1 setpoint ?= ?? ?? ?? 20

Three Bus PV Case Example For this three bus case we have = + + x x ( ) ( ) Q P P P P P P 2 2 2 2 G D = = x ( ) f x 0 3 3 3 Q 3 G + D ( ) x V 2 2 2 D Line Z = 0.1j 0.941 pu -7.469 Deg One 1.000 pu Two 200 MW 100 MVR 170.0 MW 68.2 MVR Line Z = 0.1j Line Z = 0.1j Three 1.000 pu 30 MW 63 MVR 21

Modeling Voltage Dependent Load So far we've assumed that the load is independent of the bus voltage (i.e., constant power). However, the power flow can be easily extended to include voltage depedence with both the real and reactive l is done by making P and Q a function of oad. This V : Di Di i n + + = ( cos sin ) ( ) 0 V V G B P P V i k ik ik ik ik Gi Di i = 1 k n + = ( sin cos ) ( ) 0 V V G B Q Q V i k ik ik ik ik Gi Di i = 1 k 22

Voltage Dependent Load Example In previous two bus example now assume the load is constant impedance, so 2 = ) 2.0 + = x P ( ) (10sin 0 V V 2 2 2 2 2 2 = ( 10cos + (10) 1.0 + = ( ) x ) 0 Q V V V 2 2 2 2 2 Now calculate the power flow Jacobian 10 ( ) 10 2 V + + cos sin 10sin 4.0 V V V + 2 2 2 20 2 2.0 = x J 10cos V 2 2 2 2 23

Voltage Dependent Load, cont'd 0 1 (0) = = x Again set 0, guess v Calculate 2 ) 2.0 + (10sin V V 2.0 1.0 2 2 V 2 + (0) = = x f( ) 2 2 ( 10cos + ) (10) 1.0 V V 2 2 2 2 10 0 4 (0) = J x ( ) 12 1 0 1 10 0 4 2.0 1.0 0.1667 0.9167 (1) = = x Solve 12 24

Voltage Dependent Load, cont'd With constant impedance load the MW/Mvar load at bus 2 varies with the square of the bus 2 voltage magnitude. This if the voltage level is less than 1.0, the load is lower than 200/100 MW/Mvar 160.0 MW 120.0 MVR -160.0 MW -80.0 MVR Line Z = 0.1j 0.894 pu -10.304 Deg One 1.000 pu Two 160.0 MW 120.0 MVR 160 MW 80 MVR PowerWorld Case Name: Bus2_Intro_Z 25

Generator Reactive Power Limits The reactive power output of generators varies to maintain the terminal voltage; on an actual generator this is done by the exciter To maintain higher voltages requires more reactive power Generators have reactive power limits, which are dependent upon the generator's MW output So for PV-bus i, we have ?? These limits must be satisfied by the power flow solution. max min ?? ?? 26

Generator Reactive Limits, cont'd During power flow once a solution is obtained, check to make sure every generator reactive power output ??is within its limits If the reactive power is outside of the limits, fix Q at the max or min value, and resolve treating the generator as a PQ bus this is know as "type-switching" also need to check if a PQ generator can again regulate Rule of thumb: to raise system voltage we need to supply more Vars 27

Switching Bus Status At iteration ?, we use ?(? 1) to compute ?(?) Suppose that after N-R update ?(?) violates the reactive power limits at a PV-bus k Then, in iteration ?, we need to change bus k to a PQ-bus ??= ?? ?? ?? ???? And ?? should be included in x After computing ?(?+1), we may need to switch bus k back to a PV-one by comparing ?? max max min ?? ???? > ?? < ?? min setpoint (?+1)with ?? How do we want to code up the N-R algorithm? 28

The N-R Power Flow: 5-bus Example T2 800 MVA 345/15 kV T1 1 5 4 3 520 MVA Line 3 345 kV 50 mi 400 MVA 15 kV 800 MVA 15 kV 400 MVA 15/345 kV 345 kV 100 mi 40 Mvar 80 MW Line 2 Line 1 345 kV 200 mi 2 280 Mvar 800 MW Single-line diagram This five bus example is taken from Chapter 6 of Power System Analysis and Design by Glover, Sarma, and Overbye, 5th Edition, 2011 29

The N-R Power Flow: 5-bus Example V per unit PG per unit QG per unit PL per unit 0 QL per unit 0 QGmax per unit QGmin per unit Bus Type degrees 1 Swing 1.0 0 Table 1. Bus input data 2 Load 0 0 8.0 2.8 3 Constant voltage 1.05 5.2 0.8 0.4 4.0 -2.8 4 5 Load Load 0 0 0 0 0 0 0 0 Maximum MVA per unit R X G B Bus-to- Bus 2-4 2-5 4-5 per unit per unit per unit per unit Table 2. Line input data 0.0090 0.0045 0.00225 0.100 0.050 0.025 0 0 0 1.72 0.88 0.44 12.0 12.0 12.0 30

The N-R Power Flow: 5-bus Example Maximum TAP Setting per unit R per unit X Gc per unit Bm per unit Maximum MVA per unit per unit Table 3. Transformer input data Bus-to- Bus 1-5 0.00150 0.02 0 0 6.0 3-4 0.00075 0.01 0 0 10.0 Bus Input Data Unknowns V1 = 1.0, 1 = 0 1 P1, Q1 V2, 2 2 P2 = PG2-PL2 = -8 Q2 = QG2-QL2 = -2.8 V3 = 1.05 P3 = PG3-PL3 = 4.4 P4 = 0, Q4 = 0 Table 4. Input data and unknowns Q3, 3 3 V4, 4 V5, 5 4 5 P5 = 0, Q5 = 0 31

Five Bus Case Ybus PowerWorld Case Name: Bus5_GSO 32

Ybus Calculation Details Elements of Ybus connected to bus 2 =Y = 0 Y 21 23 1 1 + = = = 89276 . 0 + j . 9 91964 Y per unit 24 + ' ' . 0 009 j 1 . 0 R jX 24 24 1 1 + = = = 78552 . 1 + 83932 . 19 Y j per unit 25 + ' ' 0045 . 0 j . 0 05 R jX 25 25 ' ' 1 1 B B = + + + 25 24 Y j j 22 + + ' ' ' ' 2 2 R jX R jX 24 24 25 25 + . 1 j 72 . 0 88 . 0 ( = = . 2 + + 89276 91964 . 9 j 28 ) . 1 ( = 78552 19 j 83932 . 84 ) j 2 2 67828 4590 . 5847 . 28 624 . j per unit 33

Here are the Initial Bus Mismatches 34

And the Initial Power Flow Jacobian 35

And the Hand Calculation Details = Y Y Y 0 . 8 ) 0 ( 2 ) 0 ( 2 ] 24 ) 0 ( 5 84 cos( ) 0 . 1 ( 5847 . 28 { 0 . 1 95 cos( ) 0 . 1 ( 95972 . 9 + cos( ) 0 . 1 ( 9159 . 19 + 10 89 . 2 ( 0 . 8 = ) 0 ( 1 ) 0 ( 3 ) 0 ( 2 P ( 0 ( 2 V ] 22 ) 0 ( 2 ) 0 ( cos[ 2 ) ){ Y cos[ cos[ ] ] P + + + P V V V x cos[ cos[ P Y V ) 0 ( 4 V 2 2 2 21 1 21 + 22 2 23 3 23 24 4 ]} 25 5 25 624 . = ) 143 . 143 . ) )} 99972 . 7 95 ) = = 4 per unit = ) 0 ( 2 sin[ ) 0 . 1 )( unit per ) 0 ( 4 . 95 ] 1 ) 0 ( 24 ) 0 ( 4 V 95972 . 9 )( 0 . 1 ( = 91964 . 9 = ) 0 ( 2 V sin[ ] J Y 24 24 143 36

Five Bus Power System Solved One Five Four Three A A MVA MVA 395 MW 520 MW A MVA 114 Mvar 337 Mvar slack 1.000 pu 0.000 Deg 0.974 pu -4.548 Deg 1.019 pu -2.834 Deg 80 MW 40 Mvar A A MVA MVA 1.050 pu -0.597 Deg 0.834 pu -22.406 Deg Two 800 MW 280 Mvar 37

Quasi-Newton Power Flow Methods We will introduced some modified versions of the Newton-Raphson (N-R) power flow Since most of the computation in the N-R PF is associated with building and solving the Jacobian matrix, J, the focus is on trying to reduce this computation Quasi-Newton methods of a more general update ?(?+1)= ?(?) ? ?? 1?(??) Even if ? ?, the solution is still valid for if the quasi-Newton iterations do converge f(x) = 0 38

Dishonest Newton Method The simplest modification results when J is kept constant for a number of iterations, say ? iterations Sometimes known as the Dishonest Newton method ?(?+1)= ?(?) ? ?( ?/? ?) 1?(??) The approach balances increased speed per iteration, with potentially more iterations to perform There is also an increased possibility for divergence 39

Dishonest N-R Example, contd 1 x + ( 1) ( ) v ( ) 2 v v = (( ) -2) x x x (0) 2 (0) = Guess x 1. Iteratively solving we get We pay a price in increased iterations, but with decreased computation per iteration ( ) v ( ) v v 0 1 2 3 4 (honest) 1 1.5 1.41667 1.41422 1.41422 (dishonest) 1 1.5 1.375 1.429 1.408 x x 40

Numerical differentiation Another approach to avoid directly computing the Jacobian matrix J at every iteration Based on numerical approximation ??? ??? with small perturbation ?? (of about 1%) ????+ ?? ???? ?? 2 ?? 41

Secant Method Consider a linear approximation of ? ? ? = ?0+ ?1? with ? ?(? 1)= ? ?(? 1) and ? ?(?)= ? ?(?) The iterative solution to ? ? = 0 becomes 1 ? ?? ? ?? 1 ?? ?? 1 ?(?+1)= ?(?) ? ?(?) Fig. 3.8 of Crow s book, 3rd edition 42