Power System Stability-Chapter 3
Power System Stability-Chapter 3
Power System Stability-Chapter 3
5
Power System Stability
5.1 Introduction
The present day electric power system network forms a complete, nonlinear
dynamical system. Several controls are provided in order to perform the power
system into a proper way. During the normal operation, all these controls try
to bring the system to an operating equilibrium ensuring the balance of real
and reactive powers in the system. Following a disturbance, the balance of
real and reactive powers gets disturbed. The dynamical power system network
undergoes a transition period and may settle down to an operating equilibrium
with the help of the above controls, which may or may not be the same as
the predisturbance equilibrium point. The capability of the system to achieve
an operating equilibrium, after disturbances, depends on its inherent strength,
nature and amount of disturbances. The system becomes unstable if it is not
capable of regaining the operating equilibrium. Thus, a general definition of
stability is given as follows:
“The stability of a dynamical system is its property or ability to remain
in a state of operating equilibrium under normal operating conditions
and to regain an acceptable state of equilibrium after being subjected to
a disturbance”.
During the early part of the 20th century, the concern of the power system
engineers was to maximize the real power transfer from the remotely located
generating stations to the load centre. The problem of maintaining synchronous
operation started when two or more generators were connected in the network to
share the system power demand. The difficulty in maintaining the synchronous
operations was experienced specifically in the case of severe disturbances
such as network faults and outage of large generating plants. This was called
the transient stability problem. Several practical measures were suggested to
improve the transient stability including the fast exciters and protection system.
Although these measures helped in improving the synchronizing capability and
the transient stability limit, a few of these resulted into deterioration in system
damping. With poor damping, the system becomes oscillatory unstable, called
217
218 Electrical Power Systems: Analysis, Security and Deregulation
Voltage stability
The voltage stability also called load stability refers to the “ability of the
system to maintain load bus voltages within acceptable limit, following
some disturbance or change in power demand”.
The voltage stability can be further classified as follows.
(i) Large disturbance voltage stability, which is the ability of the system
to regain voltages at all the buses within the acceptable steady state
220 Electrical Power Systems: Analysis, Security and Deregulation
EG - EM
I = (5.3)
jX
Since the resistances of the machines and the transmission lines are neglected,
the power output of the generator is also the power input to the motor and
is given by
*
P = Real part of ( EG ¥ I )
= Re ; EG c G mE
* E - EM
(5.4)
jX
Power System Stability 223
d 2q
Ta = Ja = J (5.8)
dt 2
where q is the angular position of the rotor in radians at any instant of time,
and t is the time in seconds.
It is convenient to measure q with respect to a reference axis that is rotating
at the synchronous speed. If d is the angular displacement of the rotor in
electrical degrees from the synchronously rotating reference axis and ws is the
synchronous speed in electrical radians, then q can be expressed as the sum
of: (1) time varying angle wst on the rotating reference axis, and (2) the torque
angle d of the rotor with respect to the rotating reference axis. In other words,
q = wst + d electrical radians (5.9)
Differentiating with respect to t, we get
dq dd
= ws + (5.10)
dt dt
Differentiating once again with respect to t, we get
d 2q d 2d
= (5.11)
dt 2 dt 2
226 Electrical Power Systems: Analysis, Security and Deregulation
d 2q
Angular acceleration of rotor, a =
dt 2
d 2d
a= electrical radians (5.12)
dt 2
In a synchronous generator, the accelerating torque Ta is equal to the difference
of input shaft torque Tm and the output electromagnetic torque Te.
Ta = Tm – Te
d 2q d 2d
J◊ =J◊ (5.13)
dt 2 dt 2
d 2d
Ta = Tm – Te = J ◊
dt 2
Multiplying both sides by angular velocity, w, we get
d 2d
wTa = wTm - wTe = J w (5.14)
dt 2
If Pa, Pm and Pe denote the accelerating power, mechanical power input and
electrical power output respectively, we get
d 2d
Pa = Pm - Pe = J w (5.15)
dt 2
Also, since angular momentum M = Jw, therefore Eq. (5.15) can be written as
d 2d
Pa = Pm - Pe = M (5.16)
dt 2
Equation (5.16) is called the swing equation. It is a nonlinear differential
equation of the second order.
N
= where S is the MVA rating
S
N = SH
Equivalent H constant
Consider a system in which ‘n’ number of generators are connected in parallel
to the same bus bar.
Let S1, S2, S3, …, Sn be the MVA rating of individual machines
H1, H2, H3, …, Hn be the inertia constants of individual machines
N1, N2, N3, …, Nn be the kinetic energy stored in individual machines
Se be the MVA rating of equivalent machine
He be the inertia constant of equivalent machine
Ne be the kinetic energy stored in equivalent machine and Sb be the
base MVA.
The energy stored by the equivalent machine is given by the sum of energies
stored by individual machines.
N = N1 + N2 + + Nn
Se He = S1H1 + S2H2 + + SnHn
where Se = S1 + S2 + + Sn
If the base MVA, Sb, is equal to the combined MVA rating of individual
machines Se , i.e. Sb = Se, we get
ÊS ˆ ÊS ˆ ÊS ˆ
H e = H1 Á 1 ˜ + H 2 Á 2 ˜ + + H n Á n ˜
Ë b¯
S Ë b¯
S Ë Sb ¯
If the machines are identical, we have
S1 = S2 = = Sn = S
H 1 = H 2 = = Hn = H
then
Ê HS ˆ
He = n ¥ ÁË ˜
nS ¯
With identical machines
Sb = Se = n S
Substituting Sb, we get
Ê HS ˆ
He = n ¥ Á
Ë nS ˜¯
He = H
Thus the equivalent H constant of several identical machines operating in
parallel is the same as that of any one of the machines.
Power System Stability 229
Equivalent M constant of two machines
Two synchronous machines connected by a reactance can be replaced by one
equivalent machine connected through a reactance to an infinite bus as follows.
The swing equation of machine 1 is given by
d 2d1
M1 = Pm1 - Pe1 (5.17)
dt 2
The swing equation of machine 2 is given by
d 2d 2
M2 = Pm 2 - Pe 2 (5.18)
dt 2
From Eq. (5.17)
d 2d1 Pm1 - Pe1
2
= (5.19)
dt M1
From Eq. (5.18)
d 2d 2 Pm 2 - Pe 2
= (5.20)
dt 2
M2
Subtracting Eq. (5.20) from Eq. (5.19), we get
M1M 2 d 2d M P - M1 Pm 2 M 2 Pe1 - M1 Pe 2
= 2 m1 - (5.23)
M1 + M 2 dt 2 M1 + M 2 M1 + M 2
The swing equation of an equivalent machine is given by
d 2d
= Pm¢ - Pe¢
M¢ (5.24)
dt 2
From Eq. (5.24) we can say that the equivalent values of M, Pm and Pe are
given by
230 Electrical Power Systems: Analysis, Security and Deregulation
M 1M 2
M =
M1 + M 2
M 2 Pm1 - M1 Pm 2
Pm =
M1 + M 2
M 2 Pe1 - M1 Pe 2
and Pe =
M1 + M 2
Relationship between inertia constant M and inertia constant H
We have
2N 2N N
M= = =
w s 360 f 180 f
SH
= MJ · s/electrical degree (since N = SH)
180 f
If the angle is expressed in radians
SH
M= MJ · s/electrical radian
pf
EXAMPLE 5.2 The moment of inertia of a 4 pole, 100 MVA, 11 kV,
3-f, 0.8 power factor, 50 Hz turbo alternator is 10000 kg · m2. Calculate H
and M.
Solution:
J = 10000 kg · m2
120 f 120 ¥ 50
Ns = = = 1500 rpm
p 4
N s 1500
ns = = = 25 rps
60 60
ws = 2p ns = 50p
1 1
N= J w s2 = ¥ 10000 ¥ (50p ) 2 = 123.37 MJ
2 2
N 123.37
H= = = 1.2337 MJ/MVA
S 100
SH 100 ¥ 1.2337
M= = = 0.0137 MJ ◊ s /electrical degree
180 f 180 ¥ 50
EXAMPLE 5.3 A 50 Hz, 4 pole, turbo alternator rated 100 MVA, 11 kV
has an inertia constant of 8 MJ/MVA. Determine
1. the energy stored in the rotor at synchronous speed.
2. find the rotor acceleration if the mechanical input is suddenly raised
to 80 MW for an electric load of 50 MW. (Neglect mechanical and
electrical losses).
Power System Stability 231
Solution:
(i) H = 8 MJ/MVA; S = 100 MVA
We know that
N = HS = 800 MJ
d 2q
(ii) Swing equation is M = Pa = Pm - Pe
dt 2
Here for alternator, Pm = 80 MW
Pe = 50 MW
Pa = 30 MW
SH N 800
Also M = = = = 0.0889 MJ◊ s
180 f 180 f 180 ¥ 50
d 2d Pa 30
Acceleration= = = = 337.5 electrical degree/s
dt 2
M 0.0889
EXAMPLE 5.4 A 50 Hz, 4 pole turbo generator rated 20 MVA, 11 kV has
an inertia constant of H = 9 kW· s/kVA. Find the kinetic energy stored in the
rotor at synchronous speed. Find the acceleration, if the input less the rotational
losses is 26800 HP and the electrical power developed is 16 MW.
Solution:
S = 20 MVA, 11 kV,
H = 9 kW· s/kVA = 9 kJ/kVA
H = 9 MJ/MVA
Kinetic energy = N = HS = 9 20 = 180 MJ
Pa = Pm – Pe
Pm = 26800 746 = 19992800 W = 19.99 MW
Pa = 19.99 – 16 = 3.99 MW
N 180
M= = = 0.02 MW ◊ s 2 / electrical degree
180 f 180 ¥ 50
We know that
d 2d d 2d Pa 3.99
M 2
= Pa fi 2
= =
dt dt M 0.02
d d 2
Acceleration 2
= 199.5 electrical degree/s 2
dt
EXAMPLE 5.5 A power station with 4 generators each 80 MVA, 8 MJ/MVA
is in proximity with another power station having 3 generators each 200 MVA,
3.5 MJ/MVA. Determine the inertia constant of a single equivalent machine
for use in stability studies. Assume a base value of 100 MVA.
232 Electrical Power Systems: Analysis, Security and Deregulation
Solution:
4 generators are each of 80 MVA, 8 MJ/MVA
3 generators are each of 200 MVA, 3.5 MJ/MVA
We know that
HeSe = H1S1 + H2S2 +
4 3
= Âi =1
H1S1 + ÂH S
i =1
2 2
= (4 80 8) + (3 200 3.5)
= 2560 + 2100 = 4660 MJ
4660 4660
He = = = 46.6 MJ/MVA
Se 100
EXAMPLE 5.6 Two turbo alternators specified below are interconnected
using a short line:
Machine 1 : 4 poles, 50 Hz, 125 MVA, 0.8 lag, 25000 kg · m2
Machine 2 : 4 poles, 50 Hz, 150 MVA, 0.9 lag, 20000 kg · m2
Determine the inertia constant of the single equivalent machine on a base of
150 MVA.
Solution:
Machine 1
N1
M1 =
180 f
1
N1 = J w s2
2
120 ¥ 50
Ns = = 1500 rpm
4
1500
ns = = 25 rps
60
1
N1 = ¥ 25000 ¥ w s2
2
ws = 2p ns = 157.0796 rad/s = 308.425 MJ
N1
M1 = = 0.03426 MJ ◊ s/electrical degree
180 f
M1 in p.u. on a base of 150 MVA = 2.2846 10–4 p.u.
Machine 2
ns = 25 rps
1 2 1 2
N2 = J w s = ¥ 20000 ¥ (157.0796) = 246.74
2 2
Power System Stability 233
N2
M2 = = 0.274155 MJ ◊ s/electrical degree
180 f
M2 in p.u on a base of 150 MVA = 1.8277 10–4 p.u.
M 1M 2
M = = 1.01538 ¥ 10-4 p.u.
M1 + M 2
EXAMPLE 5.8 Two turbo alternators given below are interconnected using
a short line.
Machine 1 : 4 poles, 50 Hz, 75 MVA, 0.8 lag, 30000 kg · m2
Machine 2 : 2 poles, 50 Hz, 100 MVA, 0.85 lag, and 10000 kg · m2
Determine the inertia constant of the single equivalent machine on a base of
200 MVA.
Solution:
Machine 1
120 f
ws = 2p ns = 2p ¥ = 157.079
4 ¥ 60
234 Electrical Power Systems: Analysis, Security and Deregulation
1 1
Kinetic energy = J w s2 = ¥ 30000 ¥ (157.079)2
2 2
N1 = 370
N1
H1 = = 4.935 MJ/MVA
S1
Machine 2
ws = 314.159
1 1
N2 = J w s2 = ¥ 10000 ¥ (314.159) 2 = 493.48
2 2
N 2 493.48
H2 = = = 4.9348 MJ/MVA
S2 100
Equivalent inertia constant
75
¥ 4.935
SH1
M1 = = 200 = 2.056 ¥ 10-4 p.u.
180 f 180 ¥ f
SH 2
M2 = = 2.738 ¥ 10-4 p.u.
180 f
M 1M 2
M= = 1.174 ¥ 10-4 p.u.
M1 + M 2
Úd Ú
-1/ 2 2Pa
dd = ◊ dt
M
d -1/ 2 +1 2Pa
= ◊t
1 M
- +1
2
d 1/ 2 2Pa
= ◊t (5.26)
1 M
2
1 2Pa Pa
d 1/ 2 = 2 M ◊ t = 2M
◊t
dd 2 Pa Pa 2 Pa2
= ◊ ◊t = ◊t
dt M 2M 2M 2
dd ÊP ˆ
= Á a ˜ t is the change in rotor angle at anytime
dt ËM¯
d d Pa
Change in rotor angle = = ◊t
dt M
Time period for five cycles = no. of cycles time period for the given frequency
1
= 5 ¥ = 0.1 s
50
dd
Change in rotor angle = = 15 0.1 = 1.5 electrical rad/s
dt
dd
= 1.5 electrical rad/s
dt
Since pole pair p = 1, we have qm = qe
dd
= 1.5 mechanical rad/s = 2p ns
dt
1.5
ns = = 0.2387 rps
2p
Ns = 0.2387 60 = 14.32 rpm
Since the machine is subjected to constant acceleration, the speed will increase
by 14.32 rpm.
New speed = synchronous speed + Ns
120 f
Synchronous speed = = 3000 rpm
p
New speed = 3014.32 rpm
EXAMPLE 5.10 A 200 MVA, 11 kV, 50 Hz, 4 pole turbo alternator has an
inertia constant of 6 MJ/MVA.
(i) Determine the stored energy in the rotor at synchronous speed.
(ii) The machine is operating at a load of 120 MW when the load suddenly
increases to 160 MW. Determine the rotor retardation. Neglect losses.
(iii) The retardation calculated above is maintained for 5 cycles. Determine
the change in power angle and the rotor speed in rpm at the end of
this period.
Solution:
H = 6 MJ/MVA, S = 200 MVA
(i) N = HS = 6 200 = 1200 MJ
(ii) Pa = 40 (retardation)
N 1200
M= = = 0.1333 MJ ◊ s/ electrical degree.
180 f 180 ¥ 50
d 2d Pa 40
= = = 300 electrical degree / s 2
dt 2
M 0.1333
Power System Stability 237
We know that change in rotor angle
d d Pa
= ◊t
dt M
where t is the time period for acceleration or retardation.
1
(iii) Here for five cycles t = 5 ¥ = 0.1 s
50
dd
= 300 0.1 = 30 electrical degree/s
dt
p
= 30 ¥ electrical rad/s = 0.5236 electrical rad/s
180
0.5236
= mechanical rad/s
p
Pole pairs, p = 2
dd 0.5236
= = 0.2618 mechancial rad/s = 2pns
dt 2
0.2618
ns = = 0.04167 rps
2p
Ns = ns 60 = 2.5 rpm
Since the retardation is for 5 cycles,
New speed = synchronous speed – Ns
New speed = 1500 – 2.5 = 1497.5 rpm
dd
Multiplying both sides by 2
dt
d d d 2d P dd
2 =2 a
dt dt 2
M dt
d ÈÊ dd ˆ ˘
2
P dd
ÍÁ ˜ ˙ = 2 a
dt Î Ë dt ¯ ˚ M dt
ÈÊ dd ˆ 2 ˘ P
d ÍÁ ˜ ˙ = 2 a dd
Ë
Î dt ˚¯ M
Upon integration,
2 dt
Ê dd ˆ Pa
ÁË ˜¯ =
dt
d0
Ú 2 M dd
dt
2
=
M Ú P dd
d0
a
dt
dd 2
dt
=
M Ú P dd
d0
a
If dd/dt = w, the velocity of the displacement angle with respect to the infinite
bus, then
dt
dd 2
w=
dt
=
M Ú P dd
d0
a
If the machine is continuously swinging, then the above equation will be non-
zero. The stability is indicated by the zero value, i.e. if the integral vanishes.
i.e. Ú P dd
a =0
Ú (P m - Pe ) d d = 0
Ú P dd - Ú P dd = 0
m e
Ú P dd = Ú P dd
m e
A1 = Ú ( P - P ) dd
d0
s e
A2 = Ú ( P - P ) dd
ds
e s
ds dm
and A1 – A2 = Ú ( P - P ) dd - Ú ( P - P ) dd
d0
s e
ds
e s
dm
A1 – A2 = Ú ( P - P ) dd
d0
s e
d
2( Ps - Pe )
Equation Ú
d0
M
d d = 0 is satisfied and dd /dt = 0 when A1 = A2. The
Figure 5.12 Electric power input to a motor as a function of torque angle for the
maximum sudden increase of load without loss of stability.
EXAMPLE 5.11 A synchronous motor is receiving 30% of the power which
is capable of receiving from an infinite bus. If the load on the motor is doubled,
calculate the maximum value of d during the swinging of the motor around
its new equilibrium position.
Solution: The load in the motor is 30% and let the initial operating load
angle be d0, i.e. Pe = 0.3Pmax.
244 Electrical Power Systems: Analysis, Security and Deregulation
Figure 5.13
Area A1 = Ú (0.6P
d0
max - Pmax sin d ) d d
Area A2 = Ú (P
d1
max sin d - 0.6 Pmax ) d d
p
cos d m = cos d 0 + 0.6(d 0 - d m ) ¥
180
Ê 0.6 p d m ˆ
= 0.95393 + 0.1828 - Á
Ë 180 ˜¯
cos d m = 1.1368 - 0.01047d m
d m = 58.15∞
Ê X ¥ X TL 2 ˆ
where X I = X d¢ + Á TL1
Ë X TL1 + X TL 2 ˜¯
During fault condition: Upon the occurrence of three phase fault at the
sending of the transmission line 2, the generator gets isolated from the power
system for the purpose of power flow as shown in Figure 5.16. Thus during
the period of fault lasts,
Pe2 = 0
Obviously, Pmax2 < Pmax1. The rotor now starts to decelerate as shown in
Figure 5.18. The system will be stable if a decelerating area A2 can be found
equal to accelerating area A1 before d reaches the maximum allowable value
dmax. As area A1 depends upon the clearing time tc, the clearing time must be
less than a certain value (critical clearing time) for the system to be stable.
It is to be observed that the equal area criterion helps to determine the critical
clearing angle and not the critical clearing time. Critical clearing time can be
obtained by numerical solution of the swing equation.
Figure 5.18
Prefault condition: Before the occurrence of a fault, both the lines are
connected as shown in Figure 5.20.
During fault condition: The circuit model of the system during fault is
shown in Figure 5.21. This circuit reduces to that of Figure 5.22 through one
delta to star and one star to delta conversion.
248 Electrical Power Systems: Analysis, Security and Deregulation
Using delta to star conversion, the circuit becomes as shown in Figure 5.22,
X TL 2
X TL1 ¥
XA = 2 = XB
X TL 2 X TL 2
X TL1 + +
2 2
X TL 2 X TL 2
¥
XF = 2 2
X TL 2 X TL2
X TL1 + +
2 2
Convert star connection to delta connection, the circuit becomes as shown in
Figure 5.23,
( X d¢ + X A ) X F + X F X B + ( X d¢ + X A ) X B
X II =
XF
Ê P ˆ
where d max = p - sin -1 Á m ˜
Ë Pmax 3 ¯
Integrating, we get
d d
( Pmd + Pmax 2 cos d ) d cr + ( Pmax 3 cos d + Pmd ) d max = 0
0 cr
Solution:
Prefault condition Transfer reactance during prefault operation is
0.5 ¥ 0.5
X I = 0.25 + + 0.06 = 0.56
0.5 + 0.5
|E ¢| |V | 1.25 ¥ 1.0
Pe1 = sin d = sin d = 2.232 sin d
XI 0.56
Power System Stability 251
During fault condition The fault occurs at the end of the line 2 or near
bus 2. Therefore during the short circuit fault the circuit separates by the circuit
breaker for finding the transfer reactance as shown in the Figure below. During
the clearing of fault, no power is transferred from the circuit, i.e. Pe2 = 0.
Postfault condition With the opening of the faulted line, say by simultaneous
opening of the circuit breakers CB2 and CB4, the postfault transfer reactance is
X III = 0.25 + 0.5 + 0.06 = 0.81
|E ¢| |V | 1.25 ¥ 1.0
Pe3 = sin d = sin d = 1.543 sin d
X III 0.81
The initial power angle d0 is calculated as
Pm 0 = Pe0 = 1.0 = Pmax 1 sin d 0
Ê 1.0 ˆ
d 0 = sin -1 Á = 26.32∞
Ë 2.232 ˜¯
and
Ê P ˆ
d max = 180∞ - sin -1 Á m 0 ˜
Ë Pmax 3 ¯
Ê 1.0 ˆ
d max = 180∞ - sin -1 Á = 139.6∞
Ë 1.543 ˜¯
p
Pm (d max - d 0 ) - Pmax 2 cos d 0 + Pmax 3 cos d max
cos d cr = 180
Pmax 3 - Pmax 2
p
1.0(139.6 - 26.62) - 0 + 1.543 cos139.6∞
= 180 = 0.51640
1.543 - 0
d cr = 58.9∞
Solution:
Given: infinite bus voltage V = 1.0 p.u.
terminal voltage of generator Et =1.0 p.u.
Xd = j0.3
Pm0 = Pe0 = 0.8
The current flowing into the infinite bus is
S* 0.8 - j 0.074
I= = = 0.8 - j 0.074
V *
1.0
The transfer reactance between the internal voltage and the infinite bus before
fault is
0.3 ¥ 0.3
X I = 0.3 + 0.2 + = 0.65
0.3 + 0.3
E ¢ = V + jX I I = 1.0 + j (0.65) (0.8 - j 0.074) = 1.17–26.387∞
(i) Three-phase fault occurs at the sending end of one of the lines
Prefault condition
0.3 ¥ 0.3
X I = 0.3 + 0.2 += 0.65
0.3 + 0.3
|E ¢| |V | 1.17 ¥ 1.0
Pe1 = sin d = sin d = 1.8 sin d
XI 0.65
During the fault condition: Pe2 = 0
Postfault condition Since both the lines are intact when the fault is
cleared therefore the power angle equation of prefault and postfault
are the same.
|E ¢| |V | 1.17 ¥ 1.0
Pe3 = sin d = sin d = 1.8 sin d
XI 0.65
Power System Stability 253
The initial power angle d0 is calculated as
Pm 0 = Pe0 = 0.8 = Pmax 1 sin d 0
Ê 0.8 ˆ
d 0 = sin -1 Á = 26.388∞ = 0.46055 rad
Ë 1.8 ˜¯
and
Ê P ˆ
d max = 180∞ - sin -1 Á m 0 ˜
Ë Pmax 3 ¯
Ê 0.8 ˆ
d max = 180∞ - sin -1 Á = 153.612∞
Ë 1.8 ˜¯
p
Pm (d max - d 0 ) - Pmax 2 cos d 0 + Pmax 3 cos d max
cos d cr = 180
Pmax 3 - Pmax 2
p
0.8(153.612 - 26.388) - 0 + 1.8 cos153.612∞
= 180 = 1.48 rad
1.8 - 0
d cr = 84.775∞
The critical clearing time
2 H (d cr - d 0 )
tc =
p fPm
2 ¥ 5(1.48 - 0.46055)
tc = = 0.26 s
p ¥ 60 ¥ 0.8
(ii) Three-phase fault occurs at the middle of one of the lines
Prefault condition From the reactance diagram,
Ê 0.3 ¥ 0.3 ˆ
X T = X I = 0.3 + 0.2 + Á = 0.65
Ë 0.3 + 0.3 ˜¯
|E ¢| |V | 1.17 ¥ 1.0
Pe1 = sin d = sin d = 1.8 sin d
XI 0.65
During the fault condition (fault at the middle of the transmission
line 2)
254 Electrical Power Systems: Analysis, Security and Deregulation
Dt 2
Dd n = Dd n -1 + Dw n -1 ◊ Dt + ( Pa ( n -1) - Pa ( n - 2) ) (5.32)
2M
But from Eq. (5.30)
Dt
Dw n -1 = Pa ( n - 2)
M
Substituting in Eq. (5.32)
Dt Dt 2
Dd n = Dd n -1 + ◊ Pa ( n - 2) ◊ Dt + ( Pa ( n -1) - Pa ( n - 2) )
M 2M
Dt 2
= Dd n -1 + ( Pa ( n -1) - Pa ( n - 2) )
2M
This equation gives the increment in angle during any interval in terms of the
increment for the previous interval.
Time interval t should be short enough to give the required accuracy.
If it is too short, it will increase the number of calculations to plot a swing
curve, and thus it provides accuracy.
d 2d n -1 Pa ( n -1)
As an–1 = 2
=
dt M (5.35)
Pa ( n -1)
Thus, wn–1/2 = ◊ Dt
M
Again, the angular velocity, wn–1/2 remains constant for t = (n – 1)t, during
the nth interval. From Figure 5.28(c), the displacement angle dn–1 increases to
dn over this interval by an amount dn.
To begin with, the power output at t = 0– and 0+ are averaged out and
the average value of Pa is determined. Then d is calculated from the value
during the preceding interval.
EXAMPLE 5.15 Consider a 60 Hz machine for which H = 2.7 MJ/MVA
and it is initially operating in steady state with input and output of 1 p.u.
and an angular displacement of 45 electrical degree with respect to an
infinite bus bar. Upon occurrence of a fault, assume that the input remains
constant and the output is given by Pe = d/90°. Calculate and plot the swing
curve by the step-by-step method II. Using the time interval t = 0.05 s.
Up to t = 1 s. Step-by-step method-II using the time interval t = 0.05 s and
upto t = 1 s.
Solution:
H = 2.7 MJ/MVA
At time t = 0–, Pm = Pe = 1 p.u.
At time t = 0+, Pm = 1 p.u., Pe = d/90°
SH
M = = 2.5 ¥ 10-4 p.u.
180 f
d
Pa = 1 -
90∞
2
Dt 0.052
= = 10
M 2.5 ¥ 10 -4
Now the equation for dn is given by
Dt 2
Dd n = Dd n -1 + Pa ( n -1)
M
= Dd n -1 + 10 Pa ( n -1)
From the swing curve, we see that the value of d increases and then
decreases, which makes the system stable.
Note: In the method-II, the average of the accelerating powers was calculated
from the prefault (0–) and during the fault (0+) powers, and used in the further
iterations. Similarly, if the fault is cleared by some means at any intermediate
time t, the same method of calculation of average power is to be adopted for
during the fault (t –) and the postfault (t+) time periods, and these values are
to be used in further calculations.
(Contd...)
From the swing curve, the value of d increases and then decreases, which
makes the system stable.
EXAMPLE 5.16 A 20 MVA, 3-f, 50 Hz generator delivers rated power at
unity power factor via a double circuit transmission line to an infinite bus
bar. The generator unit has a kinetic energy of 2.5 MJ/MVA at rated speed.
Its X d¢ = 0.3 p.u. The transformer circuits have negligible resistances and
each has a reactance of 0.3 p.u. on a 20 MVA base. The voltage behind the
transient reactance is 1.05 p.u. and the voltage of the infinite bus is 1 p.u.
A 3-f short circuit occurs at the middle of one of the transformer circuits.
It involves ground (a) What is the initial displacement angle of the machine?
The fault is cleared in 0.4 s by simultaneous opening of CBs at both ends of
Power System Stability 269
the faulted transmission line. (b) Calculate and plot the swing curve for the
system and ascertain whether the system is stable or not.
Take t = 0.05 s and tmax = 1 s.
Solution:
SH 1 ¥ 2.5
M = = = 2.778 ¥ 10 -4
180 f 180 ¥ 50
MW 20 ¥ 1 20
Pi = = = = 1 p.u.
MVA 20 20
The power angle equation for the three conditions should be taken into account.
1. Prefault; 2. During fault; 3. Postfault
Prefault condition
0.3 ¥ 0.3
X1 = X AB = 0.3 + = 0.3 + 0.15
0.3 + 0.3
X1 = 0.45 p.u.
During fault
X2 = XAB is obtained by converting the star connection to delta connection.
0.3 ¥ 0.3
\ X 2 = X AB = 0.3 + 0.3 + = 0.6 + 0.6
0.15
X 2 = 1.2 p.u.
Postfault The faulty transmission line is made out of service by simultaneous
opening of the circuit breakers at both ends and hence the faulty line is removed.
270 Electrical Power Systems: Analysis, Security and Deregulation
1.05 ¥ 1
Prefault power, Pmax 1 = = 2.333 p.u.
0.45
1.05 ¥ 1
During fault power, Pmax 2 = = 0.875 p.u.
1.2
1.05 ¥ 1
Postfault power, Pmax 3 = = 1.75 p.u.
0.6
We know that
Dt 2
dn = Dd n -1 + Pa ( n -1)
M
Dt 2 (0.05) 2
= =9
M 2.7778 ¥ 10-4
dn = dn–1 + 9Pa(n–1)
Pi = 1 p.u.
Pa = 1 – Pe
In general, Pe = Pmax sind
For the prefault condition
Pmax 1 sin d = Pm
2.333 sin d = 1
1
sin d = = 0.4286
2.333
d = 25.3808
For the during fault condition, Pe = Pmax2 sin d
For the postfault condition, Pe = Pmax3 sin d
Power System Stability 271
The calculation of d versus t is done in the following table.
t (s) Pmax sind Pe = Pa 9Pa d d degree
Pmax sin d degree
0– 2.333 0.4286 1 — — — 25.3808
0+ 0.875 0.4286 0.375 0.625 — — 25.3808
0av — — — 0.3125 2.8125 — —
2.8125
0.05 0.875 0.4724 0.4134 0.5866 5.2795 — 28.1933
8.0920
0.1 0.875 0.5918 0.5178 0.4822 4.3395 — 36.2853
12.4315
0.15 0.875 0.7515 0.6575 0.3425 3.0823 — 48.7168
15.5138
0.2 0.875 0.9006 0.7880 0.2120 1.9082 — 64.2306
17.4220
0.25 0.875 0.9894 0.8657 0.1343 1.2084 — 81.6526
18.6304
0.3 0.875 0.9839 0.8609 0.1391 1.2515 — 100.2830
19.8819
0.35 0.875 0.8646 0.7565 0.2435 2.1914 — 120.1649
22.0733
0.4– 0.875 0.6124 0.5358 — — — 142.2382
0.4+ 1.75 0.6124 1.0717 — — — 142.2382
0.4av — — 0.8038 0.1962 1.7658 —
23.8391
0.45 1.75 0.2406 0.4211 0.5789 5.2104 166.0773
29.0495
0.5 1.75 –0.2610 –0.4567 1.4567 13.1100 195.1268
42.1595
0.55 1.75 –0.8414 –1.4724 2.4724 22.2518 237.28
64.4113
0.6 1.75 –0.8508 –1.489 2.489 22.4006 301.6976
86.8119
0.65 1.75 0.4773 0.8353 0.1647 1.4825 388.5095
88.2944
0.7 1.75 0.8926 1.5620 –0.5620 –5.0577 476.8039
83.2367
0.75 1.75 –0.3427 –0.5997 1.5997 14.3973 560.0406
97.634
0.8 1.75 –0.8856 –1.5498 2.5498 22.9482 657.67
120.5822
(Contd...)
272 Electrical Power Systems: Analysis, Security and Deregulation
(Contd...)
if (t>tc)
Pa=Pi-Pm(3)*sin(delta);sd=sin(delta);delta=(delta*180)/pi;
disp([t Pm(3) sd Pm(3)*sd Pa (y*Pa) Del_delta delta])
end
Del_delta=Del_delta+(y*Pa);
time(i)=t;
del(i)=delta;
delta=delta+Del_delta;
i=i+1;
end
critical_clearing_angle=del_cri
critical_clearing_time=cri_t
plot(time,del);
title(‘SWING CURVE’);
xlabel(‘t, sec’);
ylabel(‘delta, elec. deg’);
Results:
Enter the power limits for prefault, during fault, and postfault conditions:
[2.333 .875 1.75]
Enter the frequency: 50
Enter the fault clearing time: 0.4
Enter the value of inertia constant: 2.5
Enter the initial displacement angle: 25.3808
Enter the initial steady state power: 1.
the generator is 0.3 p.u. The transmission circuits have negligible resistances
and each has a reactance of 0.2 p.u. on a 25 MVA base. The voltage behind
the transient reactance of the generator is 1.03 p.u. and the voltage of the
metropolitan system is 1.0 p.u. A three-phase short circuit occurs at the middle
of the transmission line circuit and is cleared in 0.4 s by the simultaneous
opening of the circuit breaker at both ends of the line. Calculate and plot the
swing curve of the generator for 1 s.
Solution:
1.03 ¥ 1
Prefault power, Pmax 1 = = 2.58 p.u.
0.4
1.03 ¥ 1
During fault power, Pmax 2 = = 0.936 p.u.
1.1
1.03 ¥ 1
Postfault power, Pmax 3 = = 2.06 p.u.
0.5
SH 1 ¥ 2.76
M = = = 2.56 ¥ 10 -4
180 f 180 ¥ 60
20
Output, Pe = = 0.8 p.u.
25
In general Pe = Pmaxsind
For the prefault condition,
Pe = Pmax 1 sin d
0.8 = 2.58 sin d
0.8
sin d =
2.58
0.8
d = sin -1 = 18.2∞
2.58
For the during fault condition, Pe = Pmax2 sin d
For the postfault condition, Pe = Pmax3 sin d
Immediately after the fault, d remains unchanged momentarily but output
changes.
Pe = Pmax 2 sin d = 0.936 sin 18.2∞
= 0.292 p.u.
Dt 2
Dd n = Dd n -1 + Pa ( n -1)
M
(0.05) 2
= Dd n -1 + Pa ( n -1)
2.564 ¥ 10-4
= Dd n -1 + 9.76 Pa ( n -1)
278 Electrical Power Systems: Analysis, Security and Deregulation
If we plot the swing curve, the value of d increases and then decreases,
this makes the system stable.
Results:
Enter the power limits for prefault, during fault, and postfault conditions:
[2.58 .936 2.06]
Enter the frequency: 60
Enter the fault clearing time: 0.4
Enter the value of inertia constant: 2.76
Enter the initial displacement angle: 18.2
Enter the initial steady state power: 0.8.
280 Electrical Power Systems: Analysis, Security and Deregulation
In the swing curve, the value of d increases and then decreases, this makes
the system stable.
Power System Stability 281
5.8.3 Euler’s Method
The Euler’s method is the simplest and the least accurate of all numerical
methods. It is presented here because of its simplicity. By studying this method,
we will be able to grasp the basic ideas involved in the numerical solutions of
one-dimensional equation (ODE) and can easily understand the comparatively
complex method such as the Runge–Kutta procedure.
Let us consider the first order differential equation
dx
= f ( x, t ) (5.41)
dt
Figure 5.29 illustrate the principles of applying the Euler’s method at initial
condition x = x0 at t = t0.
dx
where is the slope of the curve at (t0, x0), which can be determined
dt x = x0
from Eq. (5.42).Thus, the value of x at t = t1 = t0 + t is given by
dx
x1 = x0 + Dx = x0 + ◊ Dt (5.44)
dt x = x0
282 Electrical Power Systems: Analysis, Security and Deregulation
The Euler’s method is equivalent to using the first two terms of the Taylor
series expansion for x around the point (t0, x0)
Dt 2 Dt 3
x1 = x0 + Dt ( x0 ) + (
x0 ) + (
x0 ) + (5.45)
2! 3!
After using the Euler’s technique for determining x = x1 corresponding to
t = t1, we can take another short time step t and determine x2 corresponding
to t2 = t1 + t as follows:
dx
x2 = x1 + ◊ Dt (5.46)
dt x = x1
The subsequent values of x can be similarly determined. Hence, the computational
algorithm is
dx
xi +1 = xi + ◊ Dt (5.47)
dt x = xi
Ê dx dx ˆ
Á dt + ˜
x = x0 dt x = x1
p
x1 = x0 + Á ˜ ◊ Dt
c
(5.50)
Ë 2 ¯
Similarly,
Ê dx dx ˆ
Á dt + ˜
x = x1 dt p
x = x2
x2 = x1 + Á ˜ ◊ Dt
c
(5.51)
Ë 2 ¯
Ê dx dx ˆ
Á dt + ˜
x = xi dt p
x = xi +1
= xi + Á ˜ ◊ Dt
c
xi + (5.52)
Ë 2 ¯
This process can be repeated until the successive steps converge with
the desired accuracy.
Solution:
(iii) During the fault condition (fault at the middle of the transmission
line 2)
Using delta to star conversion, the circuit becomes
d Dw pf p ¥ 50
= ( Pm - Pe (d 0 )) = (0.8 - 0.742 sin 0.427) = 15.479
dt d0 H 5
End of the first step, t = 0.05
Predicted values are
p dd
d 0.05 = d 0 + ¥ Dt = 0.427 + 0 ¥ 0.05 = 0.427 rad
dt Dw 0
p d Dw
Dw 0.05 = Dw 0 + ¥ Dt = 0 + 15.479 ¥ 0.05 = 0.774 rad/s
dt d0
d Dw pf p p ¥ 50
= ( Pm - Pe (d 0.05 )) = (0.8 - 0.742 sin 0.427) = 15.479 rad/s
dt p
d 0.05 H 5
Power System Stability 287
The corrected values are
È dd ˘ dd
Í dt +
˙ dt
˙ ¥ Dt = 0.427 + È 0 + 0.774 ˘ ¥ 0.05 = 0.446 rad
p
d 0.05 = d 0 + Í
c Dw 0 Dw 0.05
ÍÎ 2 ˙˚ ÍÎ 2 ˙˚
È d Dw d Dw ˘
Í dt +
c Í d0 dt d 0.05 ˙˙
p
È15.479 + 15.479 ˘
Dw 0.05 = Dw 0 + ¥ Dt = 0 + Í ˙˚ ¥ 0.05
ÍÎ 2 ˙
˚ Î 2
= 0.774 rad/s
Iteration 2: Beginning of the first step at t = 0.05,
dd c
= Dw 0.05 = 0.774
dt c
Dw 0.05
d Dw pf c p ¥ 50
= ( Pm - Pe (d 0.05 )) = (0.8 - 0.742 sin 0.446) = 15.077
dt c
d 0.05 H 5
End of the second step, t = 0.1
Predicted values are
p c dd
d 0.1 = d 0.05 + ¥ Dt = 0.446 + 0.774 ¥ 0.05 = 0.485 rad
dt c
Dw 0.05
p c d Dw
Dw 0.1 = Dw 0.05 + ¥ Dt = 0.774 + 15.077 ¥ 0.05 = 1.5278 rad/s
dt c
d 0.05
Derivatives at the end of t = 0.1,
dd p
= Dw 0.1 = 1.5278 rad/s
dt p
Dw 0.1
d Dw pf p p ¥ 50
= ( Pm - Pe (d 0.1 )) = (0.8 - 0.742 sin 0.485) = 14.265 rad/s
dt p
d 0.1 H 5
The corrected values are
È dd dd ˘
Í dt +
dt p ˙
c
È 0.774 + 1.5278 ˘
= d 0.05 + Í
c c Dw 0.05 Dw 0.1 ˙
d 0.1 ¥ Dt = 0.446 + Í ˙˚ ¥ 0.05
ÍÎ 2 ˙˚ Î 2
= 0.5035 rad
È d Dw d Dw ˘
Í dt c + dt p ˙
= Dw 0.05 + Í
c c d 0.05 d 0.1 ˙
Dw 0.1 ¥ Dt
ÍÎ 2 ˙˚
È15.077 + 14.265 ˘
= 0.774 + Í
Î 2 ˙˚ ¥ 0.05 = 1.5076 rad/s
288 Electrical Power Systems: Analysis, Security and Deregulation
d Dw pf c p ¥ 50
= ( Pm - Pe (d 0.1 )) = (0.8 - 1.378 sin 0.5035) = 4.245
dt c
d 0.1 H 5
p c d Dw
Dw 0.15 = Dw 0.1 + ¥ Dt = 1.5076 + 4.245 ¥ 0.05 = 1.7197 rad/s
dt c
d 0.1
d Dw pf p p ¥ 50
= ( Pm - Pe (d 0.15 )) = (0.8 - 1.378 sin 0.579) = 1.444 rad/s
dt p
d 0.15 H 5
The corrected values are
È dd dd ˘
Í dt c + dt p ˙
c
d 0.15
c
= d 0.1 + Í Dw 0.1 Dw 0.15 ˙
¥ Dt
ÍÎ 2 ˙˚
È1.5076 + 1.7197 ˘
= 0.5035 + Í ˙˚ ¥ 0.05 = 0.584 rad
Î 2
È d Dw d Dw ˘
Í dt c + dt p ˙
= Dw 0.1 + Í
c c d 0.1 d 0.15 ˙
Dw 0.15 ¥ Dt
ÍÎ 2 ˙˚
È 4.245 + 1.444 ˘
= 1.5076 + Í
Î 2 ˙˚ ¥ 0.05 = 1.6498 rad/s
(b) Runge–Kutta method
Pe = 0.742 sin d
Pe0 = Pm0 = 0.8 = 1.93 sind0
Power System Stability 289
-1 Ê 0.8 ˆ
d 0 = sin Á = 24.48∞ = 0.427 rad
Ë 1.93 ˜¯
w 0 = 2p f = 2p ¥ 50 = 314.159
Dtt = 0.05 s
Fourth order method
Iteration 1: at t = 0,
Ist estimate: k1 = w0 t = 0 0.05 = 0
pf p ¥ 50
l1 = H ( Pm - Pe (d 0 )) Dt = 5 (0.8 - 0.742 sin 0.427) ¥ 0.05 = 0.774
Ê l ˆ Ê 0.774 ˆ
IInd estimate: k2 = Á Dw 0 + 1 ˜ ¥ Dt = Á 0 + ˜ ¥ 0.05 = 0.0194
Ë 2¯ Ë 2 ¯
pf Ê Ê k1 ˆ ˆ p ¥ 50 Ê Ê 0ˆˆ
l2 = ÁË Pm - Pe ÁË d 0 + ˜¯ ˜¯ Dt = ÁË 0.8 - 0.742 sin ÁË 0.427 + ˜ ˜ ¥ 0.05
H 2 5 2¯¯
= 0.774
Iteration 2: at t = 0.05s,
Ist estimate: k1 = w0.05 t = 0.767 0.05 = 0.0384
pf p ¥ 50
l = ( Pm - Pe (d 0.05 )) Dt = (0.8 - 0.742 sin 0.446) ¥ 0.05 = 0.754
1 H 5
Ê l ˆ Ê 0.754 ˆ
IInd estimate: k2 = Á Dw 0.05 + 1 ˜ ¥ Dt = Á 0.767 + ˜ ¥ 0.05 = 0.0572
Ë 2¯ Ë 2 ¯
pf Ê Ê k ˆˆ
l2 = Á Pm - Pe Ád 0.05 + 1 ˜ ˜ Dt
H Ë Ë 2 ¯¯
p ¥ 50 Ê Ê 0.0384 ˆ ˆ
= Á 0.8 - 0.742 sin Á 0.446 + ˜ ˜ ¥ 0.05
5 Ë Ë 2 ¯¯
= 0.734
Si* Pi - jQi
Ii = = i = 1, 2, ..., m
Vi* Vi* (5.56)
where m is the number of generators
Vi is the terminal voltage of the ith generator
Pi and Qi are the generators of real and reactive powers.
All unknown values are determined from the initial power flow solution.
The generator armature resistances are usually neglected and the voltages
behind the transient reactance are then obtained as
Ei¢ = Vi + jX d¢ I i (5.57)
Power System Stability 293
Next, all loads are converted to equivalent admittances by using the relation
Si* Pi - jQi
yi 0 = 2
= (5.58)
|Vi | |Vi |2
To include voltages behind the transient reactance, m buses are added to the n
bus power system network. The equivalent network, with all loads converted
to admittances is shown in Figure 5.30.
Nodes n + 1, n + 2, ..., n + m are the internal machine buses, i.e. the
buses behind the transient reactances. The node voltage equation, with node
0 as reference for this network, is
To simplify the analysis, all nodes other than the generator internal nodes
are eliminated using Kron’s reduction formula. To eliminate the load buses,
the bus admittance matrix in Eq. (5.59) is partitioned such that the n buses to
be removed are represented in the upper n rows. Since no current enters or
leaves the load buses, currents in the n rows is zero. The generator current is
denoted by the vector Im and the generator and load voltages are represented
by the vectors Em and Vn, respectively. Then, Eq. (5.59), in terms of sub
matrices becomes
È 0 ˘ È Ynn Ynm ˘ È Vn ˘
ÍI ˙ = Í t ˙ ÍE¢ ˙ (5.61)
Î m ˚ ÎYnm Ymm ˚ Î m˚
where, Yij are the elements of the faulted reduced bus admittance matrix,
Hi is the inertia constant of machine i expressed on the common MVA base SB.
If HGi is the inertia constant of machine i expressed on the machine rated
MVA SGi, then Hi is given by
SGi
Hi = H Gi (5.72)
Scb
Showing the electrical power of the ith generator by Pef and transforming
Eq. (5.71) into state variable mode yield
dd i
= wi (5.73)
dt
d Dw i p f0 f
= ( Pmi - Pei ) (5.74)
dt Hi
In the transient stability analysis problem, we have two state equations for
each generator. When the fault is cleared, which may involve the removal of
the faulty line, the bus admittance matrix is recomputed to reflect the change
in the network. Next the postfault reduced bus admittance matrix is evaluated
pf
and the postfault electrical power of the ith generator shown by Pi readily
pf
determined. Using the postfault power Pi , the simulation is continued to
determine the system stability, until the plots reveal a definite trend as to
stability or instability. Usually the slack generator is selected as the reference
machines are plotted. Usually, the solution is carried out for two swings to show
296 Electrical Power Systems: Analysis, Security and Deregulation
that the second swing is not greater than the first one. If the angle differences
do not increase, the system is stable. If any of the angle differences increase
indefinitely, the system is unstable. The flow chart of transient stability analysis
for a multimachine power is given in Figure 5.31.
Figure 5.31 Flow chart of transient stability analysis for a multimachine power.
Review Questions
Part-A
1. What is the power system stability?
2. How is the power system stability classified?
3. What is the rotor angle stability?
4. What is the steady state stability?
5. What is the steady state stability limit?
6. What is the transient stability?
7. What is the transient stability limit?
8. What is the dynamic stability?
9. What is the voltage stability?
10. State the causes of voltage instability.
11. Write the power angle equation and draw the power angle curve.
12. Write the expression for the maximum power transfer.
13. Write the swing equation for a SMIB (single machine connected to an
infinite bus bar) system.
Power System Stability 299
14. Define the swing curve.
15. In a three machine system having ratings G1, G2 and G3 and inertia
constants M1, M2 and M3 what are the inertia constants M and H of
the equivalent system.
16. State the assumptions made in stability studies.
17. State equal area criterion.
18. Define the critical clearing angle.
19. List the methods of improving the transient stability limit of a power
system.
20. What are the numerical integration methods of power system stability?
Part-B
1. A 400 MVA synchronous machine has H1 = 4.6 MJ/MVA and a
1200 MVA machine H2 = 3.0 MJ/MVA. Two machines operate in
parallel in a power plant. Find out Heq relative to a 100 MVA base.
2. A 100 MVA, two pole, 50 Hz generator has moment of inertia
40 103 kg · m2. What is the energy stored in the rotor at the rated
speed? What is the corresponding angular momentum? Determine the
inertia constant H.
3. The sending end and the receiving end voltages of a three-phase
transmission line at a 200 MW load are equal at 230 kV. The per phase
line impedance is j14 . Calculate the maximum steady state power
that can be transmitted over the line.
4. A single line diagram of a system is shown in the figure below.
All the values are in per unit on a common base. The power delivered
into bus 2 is 1.0 p.u. at 0.80 power factor lagging. Obtain the power
angle equation and the swing equation for the system. Neglect all losses.
phase fault occurs and the effective terminal voltage of the generator
becomes 25% of its value before the fault. When the fault is cleared, the
generator is delivering 70% of the original maximum value. Determine
the critical clearing angle.
7. Determine the critical clearing angle of the power system as shown in
the figure below for a three-phase fault at the point F. The generator
is supplying 1.0 p.u. MW power under the prefault condition.