Transcript PSA NOTES

EE6501
POWER SYSTEM ANALYSIS
1
UNIT I
INTRODUCTION
2
Power system network
3
SINGLE LINE DIAGRAM
It is a diagrammatic representation of a power system in
which the components are represented by their symbols.
4
COMPONENTS OF A POWER SYSTEM
1.Alternator
2.Power transformer
3.Transmission lines
4.Substation transformer
5.Distribution transformer
6.Loads
5
MODELLING OF GENERATOR AND
SYNCHRONOUS MOTOR
1Φ equivalent circuit of generator
1Φ equivalent circuit of synchronous motor
6
MODELLING OF TRANSFORMER
K
E2 N 2
I

 1
E1
N1 I 2
R2
K2
X
 X 1  X 2 '  X 1  22
K
R01  R1  R2 '  R1 
X 01
=Equivalent resistance referred to 1o
=Equivalent reactance referred to 1o
7
MODELLING OF TRANSMISSION LINE
Π type
T type
8
MODELLING OF INDUCTION MOTOR
1
Rr (  1) =Resistance representing load
s
R  RS  Rr ' =Equivalent resistance referred to stator
'
X  X S  X r'
=Equivalent reactance referred to stator
9
per unit=actual value/base value
Let KVAb=Base KVA
kVb=Base voltage
Zb=Base impedance in Ω
Zb
kVb 


2
MVAb
kVb 


2
KVAb
1000
10
Changing the base of per unit quantities
Let z = actual impedance(Ω)
Z b = base impedance (Ω)
Z p.u
Let
Z * MVAb
Z
Z



2
2
Zb  kVb 
 kVb 
MVAb
kVb,old & MVBb,old
kVb,new & MVBb,new
represent old base values
represent new base values
11
Z p.u ,old 
Z * MVAb ,old
 kV
b , old
Z

 (1)
2
Z p.u ,old * MVAb ,old
 kV
b , old
Z p.u , new 

 (2)
2
Z * MVAb ,new
 kV

kV

*
 kV
2
 (3)
b , new
Z p.u , new  Z p.u ,old
b , old
b , new


2
2
*
MVAb ,new
MVAb ,old
12
ADVANTAGES OF PER UNIT CALCULATIONS
 The p.u impedance referred to either side of a 1Φ
transformer is same
 The manufacturers provide the impedance value in p.u
 The p.u impedance referred to either side of a 3Φ
transformer is same regardless of the 3Φ connections YY,Δ-Y
 p.u value always less than unity.
13
IMPEDANCE DIAGRAM
• This diagram obtained by replacing each component by
their 1Φ equivalent circuit.
Following approximations are made to draw impedance
diagram
1. The impedance b/w neutral and ground omitted.
2. Shunt branches of the transformer equivalent circuit
neglected.
14
REACTANCE DIAGRAM
 It is the equivalent circuit of the power system in which
the various components are represented by their
respective equivalent circuit.
 Reactance diagram can be obtained after omitting all
resistances & capacitances of the transmission line from
impedance diagram
15
REACTANCE DIAGRAM FOR THE GIVEN POWER SYSTEM
NETWORK
16
PROCEDURE TO FORM REACTANCE
DIAGRAM FROM SINGLE DIAGRAM
1.Select a base power kVAb or MVAb
2.Select a base voltage kVb
3. The voltage conversion is achieved by means of transformer kVb on LT
section= kVb on HT section x LT voltage rating/HT voltage rating
4. When specified reactance of a component is in ohms
p.u reactance=actual reactance/base reactance
specified reactance of a component is in p.u
kV 

*
 kV 
2
X p.u ,new  X p.u ,old
b , old
b , new
2
*
MVAb ,new
MVAb ,old
17
p.u. calculation of 3 winding transformer
Zp=Impedance of primary winding
Zs’=Impedance of secondary winding
Zt’=Impedance of tertiary winding
Short circuit test conducted to find out the above 3
impedances
18
1
Z p   Z ps  Z pt  Z st ' 
2
1
'
Z s   Z ps  Z st '  Z pt 
2
1
'
Zt    Z ps  Z pt  Z st ' 
2
Z ps
= Leakage impedance measured in 1o with 2o short circuited
and tertiary open.
Z pt = Leakage impedance measured in 1o with tertiary short
circuited and 2o open.
Z st '
= Leakage impedance measured in 2o with tertiary short
circuited and 1o open and referred to primary
19
PRIMITIVE NETWORK
It is a set of unconnected elements which provides information
regarding the characteristics of individual elements. it can be
represented both in impedance & admittance form
20
BUS ADMITTANCE(Y BUS) MATRIX
Y BUS can be formed by 2 methods
1.Inspection method
2.Singular transformation
Y BUS =
 Y11 Y12   Y1n 


Y
Y


Y
21
22
2n 

Y Y   Y 
nn 
 n1 n 2
21
INSPECTION METHOD
For n bus system
Diagonal element of Y BUS
n
Yii   yij
j 1
Off Diagonal element of Y BUS
Yij   yij
22
SINGULAR TRANSFORMATION METHOD
A  y A
T
Y BUS =
Where [y]=primitive admittance
A=bus incidence matrix
23
ALGORITHM FOR FORMATION OF THE
BUS IMPEDANCE MATRIX
•
Modification of Zbus matrix involves any one of the following 4 cases
Case 1:adding a branch impedance zb from a new bus p to
the reference bus
Addition of new bus will increase the order the Zbus matrix by 1
 zoriginal 0 
Zbus ,new  

zb 
 0
(n+1)th column and row elements are zero except the diagonal
diagonal element is zb
24
Case 2: adding a branch impedance zb from a new bus p
to the existing bus q
Addition of new bus will increase the order the Zbus matrix by 1
The elements of (n+1)th column and row are the elements of
qth column and row and the diagonal element is Zqq+Zb
Case 3:adding a branch impedance zb from an existing bus p to
the reference bus
The elements of (n+1)th column and row are the elements of
qth column and row and the diagonal element is Zqq+Zb and
(n+1)th row and column should be eliminated using the following
formula
Z jk ,act  Z jk 
Z j ( n1) Z ( n1) k
Z ( n1)( n1)
j  1,2...n; k  1,2..n
25
Case 4:adding a branch impedance zb between existing buses h and q
elements of (n+1)th column are elements of bus h column –
bus q column and elements of (n+1)th row are elements of
bus h row – bus q row the diagonal element= Z b  Z hh  Z qq  2 Z hq
and (n+1)th row and column should be eliminated using the following
formula
Z jk ,act  Z jk 
Z j ( n1) Z ( n1) k
Z ( n1)( n1)
j  1,2...n; k  1,2..n
26
UNIT II
POWER FLOW ANALYSIS
27
BUS CLASSIFICATION
1.Slack bus or Reference bus or Swing bus:
|V| and δ are specified. P and Q are un specified, and to be
calculated.
2.Generator bus or PV bus or Voltage controlled bus:
P and |V| are specified. Q and δ are un specified, and to be
calculated
3.Load bus or PQ bus:
P and Q are specified. |V| and δ are un specified, and to be
calculated
28
ITERATIVE METHOD
n
I p   YpqVq
q 1
S p  Pp  jQ p  VP I p
*
Pp  jQ p
VP
*
n
  YpqVq
q 1
The above Load flow equations are non linear and
can be solved by following iterative methods.
1.Gauss seidal method
2.Newton Raphson method
3.Fast Decoupled method
29
GAUSS SEIDAL METHOD
For load bus calculate |V| and δ from Vpk+1 equation
p 1
n
P

jQ

1
p
p
k 1
k 1
k
Vp  
 YpqVq   YpqVq 
k *
Ypp  (VP )
q 1
q  p 1

For generator bus calculate Q from QPK+1 equation
p 1
n




k 1
k *
k 1
k 
Q p  1*Im (VP )   YpqVq   YpqVq  
q p

 q 1
 
30
• Check Qp,calk+1 with the limits of Qp
• If Qp,calk+1 lies within the limits bus p remains as PV bus
otherwise it will change to load bus
• Calculate δ for PV bus from Vpk+1 equation
• Acceleration factor α can be used for faster convergence
• Calculate change in bus-p voltage
k 1
k 1
Vp  Vp  Vp
k
• If |ΔVmax |< ε, find slack bus power otherwise increase
the iteration count (k)
• Slack bus power=
 S  S
G
L
31
NEWTON RAPHSON METHOD
n
Pi  Qi   Vi V j Yij ij   i   j
j 1
n
Pi   Vi V j Yij cos(ij   i   j )
j 1
n
Qi   Vi V j Yij sin(ij   i   j )
j 1
J 2    


J 4   V 
 Pi k
 P   J1
 Q    J

  3
Pi k  Pi sch
Qi k  Qi sch  Qi k
32
• Calculate |V| and δ from the following equation
 i k 1   i k   k
Vi k 1  Vi k   Vi k
• If
Pi k  
Qi k  
• stop the iteration otherwise increase the iteration count (k)
33
FAST DECOUPLED METHOD
 J2 & J3 of Jacobian matrix are zero
0    


J 4   V 
 P 
 P  J 1   


  
 P   J 1
 Q    0

 
 Q 
Q  J 4  V  
 V

V


P
  B ' 
 Vi
Q
  B ''  V
 Vi
    B ' 
1
 V    B '' 
P
V
1
Q
V
34
 i k 1   i k   k
Vi k 1  Vi k   Vi k
This method requires more iterations than NR
method but less time per iteration
It is useful for in contingency analysis
35
COMPARISION BETWEEN ITERATIVE
METHODS
Gauss – Seidal Method
1. Computer memory requirement is less.
2. Computation time per iteration is less.
3. It requires less number of arithmetic operations to
complete an iteration and ease in programming.
4. No. of iterations are more for convergence and rate of
convergence is slow (linear convergence
characteristic.
5. No. of iterations increases with the increase of no. of
buses.
36
NEWTON – RAPHSON METHOD






Superior convergence because of quadratic convergence.
It has an 1:8 iteration ratio compared to GS method.
More accurate.
Smaller no. of iterations and used for large size systems.
It is faster and no. of iterations is independent of the no. of buses.
Technique is difficult and calculations involved in each iteration are
more and thus computation time per iteration is large.
 Computer memory requirement is large, as the elements of jacobian
matrix are to be computed in each iteration.
 Programming logic is more complex.
37
UNIT III
FAULT ANALYSIS-BALANCED FAULT
38
Need for fault analysis
 To determine the magnitude of fault current throughout
the power system after fault occurs.
 To select the ratings for fuses, breakers and switchgear.
 To check the MVA ratings of the existing circuit breakers
when new generators are added into a system.
39
BALANCED THREE PHASE FAULT
 All the three phases are short circuited to each other and to earth.
 Voltages and currents of the system balanced after the symmetrical fault
occurred. It is enough to consider any one phase for analysis.
SHORT CIRCUIT CAPACITY


It is the product of magnitudes of the prefault voltage
and the post fault current.
It is used to determine the dimension of a bus bar and the
interrupting capacity of a circuit breaker.
40
Short Circuit Capacity (SCC)
SCC  V 0 I F
IF
VT

ZT
2
SCC 1
Sb ,1
VT


MVA / 
ZT
ZT p.u
SCC 3 
If 
Sb ,3
ZT
MVA
p .u
SCC 3 *106
3 *VL ,b *106
41
Procedure for calculating short circuit capacity
and fault current
 Draw a single line diagram and select common base Sb
MVA and kV
 Draw the reactance diagram and calculate the total p.u
impedance from the fault point to source (Thevenin
impedance ZT)
 Determine SCC and If
42
ALGORITHM FOR SHORT CIRCUIT ANALYSIS
USING BUS IMPEDANCE MATRIX
• Consider a n bus network. Assume that three phase fault
is applied at bus k through a fault impedance zf
• Prefault voltages at all the buses are
V1 (0) 
V (0) 
 2 
.

Vbus (0)  

V
(0)
 k 
.



Vn (0) 
• Draw the Thevenin equivalent circuit i.e Zeroing all voltage sources
and add voltage source V (0) at faulted bus k and draw the reactance
k
diagram
43
• The change in bus voltage due to fault is
Vbus
 V1 
.



.




V
 k
.



 Vn 
• The bus voltages during the fault is
Vbus ( F )  Vbus (0)  Vbus
• The current entering into all the buses is zero.the current entering
into faulted bus k is –ve of the current leaving the bus k
44
Vbus  Z bus I bus
 Z11 . Z1k .

 . . . .
Vbus   Z k1 . Z kk .

 . . . .
Z . Z
.
nk
 n1
Vk ( F )  Vk (0)  Z kk I k ( F )
Z1n  0



.  .

Z kn    I k ( F ) 


.  .

Z nn  0

Vk ( F )  Z f I k ( F )
Vk (0)
Ik (F ) 
Z kk  Z f
Vi ( F )  Vi (0)  Z ik I k ( F )
45
UNIT IV
FAULT ANALYSIS – UNBALANCED
FAULTS
46
INTRODUCTION
UNSYMMETRICAL FAULTS
o One or two phases are involved
o Voltages and currents become unbalanced and each phase is to be
treated individually
o The various types of faults are
Shunt type faults
1.Line to Ground fault (LG)
2. Line to Line fault (LL)
3. Line to Line to Ground fault (LLG)
Series type faults
Open conductor fault (one or two conductor open fault)
47
FUNDAMENTALS OF SYMMETRICAL
COMPONENTS
 Symmetrical components can be used to transform
three phase unbalanced voltages and currents to
balanced voltages and currents
 Three phase unbalanced phasors can be resolved into
following three sequences
1.Positive sequence components
2. Negative sequence components
3. Zero sequence components
48
Positive sequence components
Three phasors with equal magnitudes, equally displaced from one
another by 120o and phase sequence is same as that of original
phasors. Va1 ,Vb1 ,Vc1
Negative sequence components
Three phasors with equal magnitudes, equally displaced from one
another by 120o and phase sequence is opposite to that of original
phasors. V ,V ,V
a2
b2
c2
Zero sequence components
Three phasors with equal magnitudes and displaced from one another
by 0o Va 0 ,Vb 0 ,Vc 0
49
RELATIONSHIP BETWEEN UNBALANCED VECTORS
AND SYMMETRICAL COMPONENTS
Va  Va 0  Va1  Va 2
Vb  Vb 0  Vb1  Vb 2
Vc  Vc 0  Vc1  Vc 2
Va   1 1
V    1 a 2
 b 
Vc   1 a
1 1

A  1 a 2
1 a

1  Va 0 

a  Va1 
a 2  Va 2 
1

a
a 2 
Similarly we can obtain for currents also
50
SEQUENCE IMPEDANCE
Impedances offered by power system components to positive,
negative and zero sequence currents.
Positive sequence impedance
The impedance of a component when positive sequence
currents alone are flowing.
Negative sequence impedance
The impedance of a component when negative sequence
currents alone are flowing.
Zero sequence impedance
The impedance of a component when zero sequence currents
alone are flowing.
51
SEQUENCE NETWORK
SEQUENCE NETWORK FOR GENERATOR
positive sequence network
negative sequence network
Zero sequence network
52
SEQUENCE NETWORK FOR
TRANSMISSION LINE
positive sequence network
negative sequence network
Zero sequence network
53
SEQUENCE NETWORK FOR
TRANSFORMER
positive sequence network
negative sequence network
Zero sequence network
54
SEQUENCE NETWORK FOR LOAD
positive sequence network
negative sequence network
Zero sequence network
55
SINGLE LINE TO GROUND FAULT
Ib  0
Ic  0
Va  Zf Ia
Ia1  Ia2  Ia0  Ia / 3
Consider a fault between phase a and
ground through an impedance zf
Ea
Ia1 
Z1  Z 2  Z 0  3Z f
56
LINE TO LINE (LL) FAULT
Ia  0
I c   Ib
Vb  Vc  Ib Z f
Ia2   Ia1
Ia 0  0
Consider a fault between phase b and c
through an impedance zf
Va1  Va2  Zf Ia1
Ia1 
Ea
Z1  Z 2  3Z f
Ib   I c 
 jEa
Z1  Z 2  3Z f
57
DOUBLE LINE TO GROUND (LLG) FAULT
Ia0  0
Ia1  Ia2  Ia0  0
Vb =Vc  Z f (Ib  Ic )  3Zf Ia0
Va0  Va1  Vb  3Zf Ia0
Ea
Ia1 
Z1  Z 2 ( Z 0  3Z f ) / ( Z 2  Z 0  3Z f )
Consider a fault between phase b and c
through an impedance zf to ground
58
UNBALANCED FAULT ANALYSIS USING
BUS IMPEDANCE MATRIX
SINGLE LINE TO GROUND FAULT USING Zbus
 Consider a fault between phase a and ground through
an impedance zf at bus k
For a fault at bus k the symmetrical
components of fault current
Ik 0  Ik1  Ik 2 
Vk (0)
Z kk1  Z kk 2  Z kk 0  3Z f
59
LINE TO LINE (LL) FAULT
Consider a fault between phase b and c through an impedance zf
Ik 0  0
I k1  I k 2 
Vk (0)
Z kk1  Z kk 2  Z f
60
DOUBLE LINE TO GROUND (LLG) FAULT
Consider a fault between phase b and c through an impedance zf to ground
Ik1 
Ik 2
Ik
0
Vk (0)
2
0
f
Z
(
Z

3
Z
)
Z kk 1  kk2 kk 0
Z kk  Z kk  3Z f
Vk (0)  Z kk 1I k 1

Z kk 2
Vk (0)  Z kk 1I k 1

Z kk 0  3Z f
Ik ( F )  Ik b  Ik c
61
BUS VOLTAGES AND LINE CURRENTS
DURING FAULT
Vi 0 ( F )  0  Z ik 0 I k 0
Vi1 ( F )  Vi 0 (0)  Z ik 1I k 1
Vi 2 ( F )  0  Z ik 2 I k 2
Iij 0 
Iij1 
Iij 2 
Vi 0 ( F )  V j 0 ( F )
Z ij 0
Vi1 ( F )  V j1 ( F )
Z ij1
Vi 2 ( F )  V j 2 ( F )
Z ij 2
62
UNIT V
STABILITY ANALYSIS
63
STABILITY
 The tendency of a power system to develop restoring
forces equal to or greater than the disturbing forces to
maintain the state of equilibrium.
 Ability to keep the machines in synchronism with another
machine
64
CLASSIFICATION OF STABILITY
 Steady state stability
Ability of the power system to regain synchronism after small and
slow disturbances (like gradual power changes)
 Dynamic stability
Ability of the power system to regain synchronism after small
disturbances occurring for a long time (like changes in turbine
speed, change in load)
 Transient stability
This concern with sudden and large changes in the network
conditions i.e. . sudden changes in application or removal of loads,
line switching operating operations, line faults, or loss of excitation.
65
 Steady state limit is the maximum power that can be
transferred without the system become unstable when
the load in increased gradually under steady state
conditions.
 Transient limit is the maximum power that can be
transferred without the system becoming unstable when
a sudden or large disturbance occurs.
66
Swing Equation for Single Machine Infinite
Bus System
• The equation governing the motion of the rotor of a
synchronous machine
d 2m
J 2  Ta  Tm  Te
dt
where
J=The total moment of inertia of the rotor(kg-m2)
 m =Singular displacement of the rotor
Tm=Mechanical torque (N-m)
Te=Net electrical torque (N-m)
Ta=Net accelerating torque (N-m)
67
 m  smt   m
d m
d m
 sm 
dt
dt
d 2 m d 2 m
 2
2
dt
dt
d 2 m
J m 2  pa  pm  pe
dt
• Where pm is the shaft power input to the machine
pe is the electrical power
pa is the accelerating power
68
J m  M
d 2 m
M
 pa  pm  pe
dt 2
2H
M 
S machine
sm
pa
pm  pe
2 H d 2 m


sm dt 2
S machine
S machine
2 H d 2
 pa  pm  pe
s dt 2
H=machine inertia constant
s  2 f
H d 2
 pa   pm   pe
 f 0 dt 2
f
d 2  f 0

 pm  p2 max sin    0 pa
2
dt
H
H
d
 
dt
d   f 0
d 2

pa 
p.u
dt
H
dt 2
p.u
δ and ωs are in electrical
radian
69
Swing Equation for Multimachine System
Smachine =machine rating(base)
S system
=system base
H system d 2
 pa  pm  pe p.u
2
 f dt
Smachine
H system  H machine
S system
70
Rotor Angle Stability
• It is the ability of interconnected synchronous machines
of a power system to maintain in synchronism. The
stability problem involves the study of the electro
mechanical oscillations inherent in power system.
• Types of Rotor Angle Stability
1. Small Signal Stability (or) Steady State Stability
2. Transient stability
71
Voltage Stability
 It is the ability of a power system to maintain steady
acceptable voltages at all buses in the system under
normal operating conditions and after being subjected to
a disturbance.
 The major factor for instability is the inability of the
power system to meet the demand for reactive power.
72
• Mid Term Stability
It represents transition between short term and long
term responses.
Typical ranges of time periods.
1. Short term : 0 to 10s
2. Mid Term : 10 to few minutes
3. Long Term : a few minutes to 10’s of minutes
• Long Term Stability
Usually these problem be associated with
1. Inadequacies in equipment responses.
2. Poor co-ordination of control and protection equipment.
3. Insufficient active/reactive power reserves.
73
Equal Area Criterion
• This is a simple graphical method to predict the transient
stability of two machine system or a single machine
against infinite bus. This criterion does not require swing
equation or solution of swing equation to determine the
stability condition.
• The stability conditions are determined by equating the
areas of segments on power angle diagram.
74
Power-angle curve for equal area criterion
multiplying swing equation by
on both sides
Multiplying both sides of the above equation by dt and then integrating between two arbitrary angles δ0 and δc
75
Once a fault occurs, the machine starts accelerating. Once the fault is cleared, the
machine keeps on accelerating before it reaches its peak at δc ,
The area of accelerating A1
The area of deceleration is given by A2
If the two areas are equal, i.e., A1 = A2, then the power system will be stable
76
Critical Clearing Angle (δcr) maximum allowable value of the clearing time and angle for
the system to remain stable are known as critical clearing time and angle.
δcr expression can be obtained by substituting δc = δcr in the equation A1 = A2
Substituting Pe = 0 in swing equation
Integrating the above equation
77
Replacing δ by δcr and t by tcr in the above equation, we get
the critical clearing time as
78
Factors Affecting Transient Stability
• Strength of the transmission network within the system
and of the tie lines to adjacent systems.
• The characteristics of generating units including inertia of
rotating parts and electrical properties such as transient
reactance and magnetic saturation characteristics of the
stator and rotor.
• Speed with which the faulted lines or equipments can be
disconnected.
79
Numerical Integration methods
 Modified Euler’s method
 Runge-Kutta method
80
MODIFIED EULER’S METHOD
• Using first derivative of the initial point next point is obtained
dX
•
the step t1  t0  t
X 1 p  X 0  t
dt
• Using this x1p dx/dt at x1p=f(t1, x1p)
• Corrected value is
  dx 
  dt 

 X0
X 1P  X 0  



  dx 
  dt 

 Xi
X ic1  X i  





 t




 dx 

 p 
 dt  X i1 
t

2


 dx 


 dt  X1p
2
81
Numerical Solution of the swing equation
• Input power pm=constant
• At steady state pe=pm,
 pm 

 p1max 
 0  sin 1 
p1max 
E' V
X1
• At synchronous speed
 0  0
p2 max 
E' V
X2
82
The swing equation
H d 2
 pa   pm   pe
2
 f 0 dt
 f0
d 2  f 0

pa
 pm  p2max sin   
2
dt
H
H
d
 
dt
d   f 0
d 2

pa  2
dt
H
dt
Applying Modified Eulers method to above equation
t1  t0  t
 d 
 t
dt

i
 i p1   i  
 d  
ip1  i  
 t
 dt i
83
•
The derivatives at the end of interval
 d 
p

 p  i 1
 dt  i1
  f0

 d  
pa 

 p 
 dt  i1  H
 p
i 1
The corrected value
  d 

 d 


 p 
  dt 
dt



 i1 


c
i
 i 1   i  
t


2




  d  
 d  




 p
 dt
dt



 i1

i
ic1  i  

2




 t



84
Runge-Kutta Method
•
•
•
•
•
Obtain a load flow solution for pretransient conditions
Calculate the generator internal voltages behind transient reactance.
Assume the occurrence of a fault and calculate the reduced admittance matrix
Initialize time count K=0,J=0
Determine the eight constants
K1k  f1 ( k ,  k ) t
l1k  f 2 ( k ,  k ) t
K1k
l1k
k
K  f1 ( 
,   ) t
2
2
Kk
lk
l2k  f 2 ( k  1 ,  k  1 ) t
2
2
k
K2
l2k
k
k
k
K 3  f1 ( 
,   ) t
2
2
k
k
K
l
l3k  f 2 ( k  2 ,  k  2 ) t
2
2
k
K3
l3k
k
k
k
K 4  f1 ( 
,   ) t
2
2
Kk
lk
l4k  f 2 ( k  3 ,  k  3 ) t
2
2
K1k  2 K 2k  2 K 3k  K 4k 

k
 
6
k
k
l1  2l2  2l3k  l4k 

k
 
6
k
2
k
85
•
Compute the change in state vector

k
K

 k 
•
k
1
 2 K 2k  2 K3k  K 4k 
6
 l1k  2l2k  2l3k  l4k 
6
Evaluate the new state vector
 k 1   k   k
 k 1   k   k
•
Evaluate the internal voltage behind transient reactance using the relation
E pk 1  E pk cos  pk 1  j E pk sin  pk 1
•
•
•
•
Check if t<tc yes K=K+1
Check if j=0,yes modify the network data and obtain the new reduced
admittance matrix and set j=j+1
set K=K+1
Check if K<Kmax, yes start from finding 8 constants