Transient stability analysis

ArnabNandi4 6,943 views 73 slides Jun 13, 2017
Slide 1
Slide 1 of 73
Slide 1
1
Slide 2
2
Slide 3
3
Slide 4
4
Slide 5
5
Slide 6
6
Slide 7
7
Slide 8
8
Slide 9
9
Slide 10
10
Slide 11
11
Slide 12
12
Slide 13
13
Slide 14
14
Slide 15
15
Slide 16
16
Slide 17
17
Slide 18
18
Slide 19
19
Slide 20
20
Slide 21
21
Slide 22
22
Slide 23
23
Slide 24
24
Slide 25
25
Slide 26
26
Slide 27
27
Slide 28
28
Slide 29
29
Slide 30
30
Slide 31
31
Slide 32
32
Slide 33
33
Slide 34
34
Slide 35
35
Slide 36
36
Slide 37
37
Slide 38
38
Slide 39
39
Slide 40
40
Slide 41
41
Slide 42
42
Slide 43
43
Slide 44
44
Slide 45
45
Slide 46
46
Slide 47
47
Slide 48
48
Slide 49
49
Slide 50
50
Slide 51
51
Slide 52
52
Slide 53
53
Slide 54
54
Slide 55
55
Slide 56
56
Slide 57
57
Slide 58
58
Slide 59
59
Slide 60
60
Slide 61
61
Slide 62
62
Slide 63
63
Slide 64
64
Slide 65
65
Slide 66
66
Slide 67
67
Slide 68
68
Slide 69
69
Slide 70
70
Slide 71
71
Slide 72
72
Slide 73
73

About This Presentation

The power at the consumer end is often subjected to changes due to the variation of load or due to disturbances induced within the length of transmission line. For this reason the term power system stability is of utmost importance. The aim of this project is to develop an algorithm for determining ...


Slide Content

A
Project Report
On

TRANSIENT STABILITY ANALYSIS
Project Submitted To

Maulana Abul Kalam Azad University of Technology





In Partial Fulfillment of the Academic Requirements for the Degree of
Bachelor of Technology in Electrical Engineering


By








Under the Guidance of
Prof. Amitava Sil
Principal, Academy of Technology

Department of Electrical Engineering






2017
Arnab Nandi
Debayan Dhar
Indrajit Bose
Indrajit Das
Kaushik Sadhukhan
Anwesha Nag
Roll No: 16901613025
Roll No: 16901613042
Roll No: 16901613047
Roll No: 16901613048
Roll No: 16901613049
Roll No: 16901614127
ACADEMY OF TECHNOLOGY
Adisaptagram, Hooghly
Hooghly: 712121, West Bengal

This is to certify that the project work entitled “TRANSIENT STABILITY ANALYSIS ”
submitted to Maulana Abul Kalam Azad University of Technology in the partial
fulfillment of the requirement for the award of Bachelor of Technology in Electrical
Engineering, is original work carried out under my guidance by:
Arnab Nandi with University Roll No: 16901613025
Debayan Dhar with University Roll No: 16901613042
Indrajit Bose with University Roll No: 16901613047
Indrajit Das with University Roll No: 16901613048
Kaushik Sadhukhan with University Roll No: 16901613049
Anwesha Nag with University Roll No: 16901614127
The matter embodied in this project is genuine work done by the students and has not been
submitted whether to this University or to any other University/Institute for the fulfillment
of the requirement of any course of study.










ACADEMY OF TECHNOLOGY
DEPARTMENT OF ELECTRICAL ENGINEERING




Signature of H.O.D
Prof. Sandip Saha Chowdhury
Date:
Signature of Project Guide
Prof. Amitava Sil
Date:
Certificate of Originality

ACKNOWLEDGEMENT


We express our great deal of thanks to our college for providing us an
opportunity to improve our skills by assigning us a project. Our special thanks
to Prof. Amitava Sil for guiding us for such a project.

We also express our gratitude towards Dr. Soma Biswas for guiding us through
various topics related with this project.

We are greatly thankful to Prof. Sandip Saha Chowdhury,
HOD, Department of Electrical Engineering, AOT and all other faculty
members who have guided us thoroughly and have given important details about
this project and helping us in making it a success.


ARNAB NANDI
Roll No: 16901613025

DEBAYAN DHAR
Roll No: 16901613042

INDRAJIT BOSE
Roll No: 16901613047

INDRAJIT DAS
Roll No: 16901613048
KAUSHIK SADHUKHAN
Roll No: 16901613049

ANWESHA NAG
Roll No: 16901614127

CONTENTS

Genesis of the Project ................................................................................................................. 1

Introduction to Power System Stability .................................................................................. 3
Definition ........................................................................................................................ 4
Classification of Power System Stability ........................................................................ 5

Overview of Synchronous Generator ...................................................................................... 7
Synchronous Machine under Short Circuit Conditions .................................................. 9

Analysis of Swing Equation ................................................................................................... 12
Derivation of Swing Equation ...................................................................................... 13
Linearization of Swing Equation .................................................................................. 17
Numerical Solution of Swing Equation ........................................................................ 19

Damping and Damper Windings ........................................................................................... 22
Purposes of Damper Windings ..................................................................................... 23
Effect of Damping on Stability ..................................................................................... 24

System Under Study ............................................................................................................... 25

Mathematical Formulation .................................................................................................... 27
Mathematical Approach ................................................................................................ 28
Calculation of Damping Constant ................................................................................. 29
Critical Clearing Angle ................................................................................................. 30
Critical Clearing Time .................................................................................................. 32
Effect of Switching Operation on Stability ................................................................... 33

Developing the MATLAB Programs ..................................................................................... 38
Algorithm to Plot Power-Angle Curve and calculate Critical Clearing Angle ............. 39
Algorithm to Plot Delta v/s Time when Damping is Not Considered .......................... 41
Algorithm to Plot Delta v/s Time when Damping is Considered ................................. 43
Program to Plot Power-Angle Curve and calculate Critical Clearing Angle ............... 45
Program to Plot Delta v/s Time when Damping is Not Considered ............................. 47
Algorithm to Plot Delta v/s Time when Damping is Considered ................................. 48

Results and Discussions .......................................................................................................... 49
Variation of Damping Factor ........................................................................................ 50
Critical Clearing Angles ............................................................................................... 51
Fault Clearing Times with Damping Neglected ........................................................... 55
Fault Clearing Times with Damping Considered ........................................................ 60
Effect of Inertia Constant .............................................................................................. 64

Conclusion ............................................................................................................................... 66
Aspects of Future Work ................................................................................................ 68

References ................................................................................................................................. 69

Page | 1

Genesis of the Project




The importance of power system stability is increasingly becoming one of the most limiting
factors for system performance. By the stability of a power system, we mean the ability of
the system to remain in operating equilibrium, or synchronism, while disturbances occur on
the system. There are three types of stability, namely, steady-state, dynamic and transient
stability.

In the study of electric power systems, depending upon the magnitude and nature of the
disturbance, stability studies are often classified into the three following categories:

 Steady-state stability – refers to the stability of a power system subject to small and
gradual changes in load.

 Dynamic stability – refers to the stability of a power system subject to a relatively small
and sudden disturbance, the system can be described by linear differential equations, and
the system can be stabilized by a linear and continuous supplementary stability control.

 Transient stability – refers to the stability of a power system subject to a sudden and
severe disturbance beyond the capability of the linear and continuous supplementary
stability control, and the system may lose its stability at the first swing unless a more
effective countermeasure is taken. For transient stability analysis and control design,
since the disturbances are large, the aforementioned linearization of the differential
equations cannot be used.

The problem of interest is one where a power system operating under a steady load condition
is perturbed, causing the readjustment of the voltage angles of the synchronous machines. If
such an occurrence creates an unbalance between the system generation and load, it results in
the establishment of a new steady – state operating condition, with the subsequent
adjustment of the voltage angles. The perturbation could be a major disturbance such as the
loss of a generator, a fault or the loss of a line, or a combination of such events. It could also
be a small load or random load changes occurring under normal operating conditions.
Adjustment to the new operating condition is called the transient period. The system
behavior during this time is called the dynamic system performance, which is of concern in
defining system stability. The main criterion for stability is that the synchronous machines
maintain synchronism at the end of the transient period.

So we can say that if the oscillatory response of a power system during the transient period
following a disturbance is damped and the system settles in a finite time to a new steady
operating condition, it can be inferred that the system is stable, else the system is considered
unstable.

This primitive definition of stability requires that the system oscillations be damped. This
condition is sometimes called asymptotic stability and means that the system contains
inherent forces that tend to reduce oscillations. This is a desirable feature in many systems
and is considered necessary for power systems. The definition also excludes continuous.

Page | 2


The problem for this project is to investigate the stability of a synchronous machine
connected to finite bus-bar through transformer and a double circuit line. A synchronous
machine is initially operated on steady load with a constant load angle ??????. On sudden change
(unloading) as it happens in case of system short circuits, the speed rises so that the rotating
mass maintains torque balance by taking some mechanical energy. To have the speed
synchronous, the supplied mechanical energy must be withdrawn from the mass. With
increase in speed ?????? increases. Thus due to the unbalance created between mechanical input
and electrical output a new operating state is necessary. The transient following system
perturbation is oscillatory in nature but if the system is stable these oscillations will be
damped towards a new quiescent operating condition.

From the aforesaid it is clear, that under normal operating conditions, the relative position of
the rotor axis and the resultant magnetic field axis is fixed. The angle between the two is
known as the power angle or torque angle ??????. During any disturbance, rotor will decelerate or
accelerate with respect to the synchronously rotating air gap mmf, a relative motion begins.
The equation describing the relative motion is known as the swing equation.

The transient stability analysis requires the solution of a system of coupled non-linear
differential equations. In general, no analytical solution of these equations exists. However,
techniques are available to obtain approximate solution of such differential equations by
numerical methods and one must therefore resort to numerical computation techniques
commonly known as digital simulation.

This project presents a cross-platform application which provides a better comprehension of
power system stability using numerical solution of swing equation. A distinctive feature of
this application is interactive plotting of power angle curves with on-line parameter
adjustment showing accelerating / decelerating areas, operating and clearing angles, stability
margins, etc. It computes the critical clearing angle (CCA) and other stability parameters. It
provides the user with a choice to select the computational method and its step size for
calculating the solution to the swing equation.

In general a power system has large number of generators. The system can be reduced by
considering an equivalent single machine connected to infinite bus. In this project we
adapted a linearized single-machine infinite bus model for design and study of stability is
done with the solution of swing equation using numerical methods.

MATLAB
®
programming and Graphical User Interface (GUI) have been used to study the
stability of a single machine infinite bus (SMIB) system.

Page | 3





Chapter 1

Introduction to
Power System Stability

Page | 4

Definition


At present the demand for electricity is rising phenomenally especially in developing
country like India. This persistent demand is leading to operation of the power system at its
limit. On top of this the need for reliable, stable and quality power is also on the rise due to
electric power sensitive industries like information technology, communication, electronics
etc. In this scenario, meeting the electric power demand is not the only criteria but also it is
the responsibility of the power system engineers to provide a stable and quality power to
the consumers.

Power System Stability, its classification, and problems associated with it have been
addressed by many CIGRE and IEEE publications. The CIGRE study committee and IEEE
power systems dynamic performance committee defines power system stability as:

"Power system stability is the ability of an electrical power system, for given operating
conditions, to regain its state of operating equilibrium after being subjected to a physical
disturbance, with the system variables bounded, so that the entire system remains intact
and the service remains uninterrupted”.

The first electric power system was a dc system built by Edison in 1882. The subsequent
power systems that were constructed in the late 19th century were all dc systems. However
despite the initial popularity of dc systems by the turn of the 20th century ac systems
started to outnumber them. The ac systems were thought to be superior as ac machines
were cheaper than their dc counterparts and more importantly ac voltages are easily
transformable from one level to other using transformers. The early stability problems of ac
systems were experienced in 1920 when insufficient damping caused spontaneous
oscillations or hunting. These problems were solved using generator damper winding and
the use of turbine-type prime movers.

The stability of a system refers to the ability of a system to return back to its steady state
when subjected to a disturbance.

The disturbances mentioned in the definition could be faults, load changes, generator
outages, line outages, voltage collapse or some combination of these. Power system
stability can be broadly classified into rotor angle, voltage and frequency stability. Each of
these three stabilities can be further classified into large disturbance or small disturbance,
short term or long term. The classification is depicted in Fig 1.1.

Page | 5

Classification of Power System Stability

 Rotor angle stability: It is the ability of the system to remain in synchronism when
subjected to a disturbance. The rotor angle of a generator depends on the balance
between the electromagnetic torque due to the generator electrical power output and
mechanical torque due to the input mechanical power through a prime mover.
Remaining in synchronism means that all the generators electromagnetic torque is
exactly balanced by the mechanical torque. If in some generator the balance
between electromagnetic and mechanical torque is disturbed, due to disturbances in
the system, then this will lead to oscillations in the rotor angle. Rotor angle stability
is further classified into small disturbance angle stability and large disturbance
angle stability.
 Small-disturbance or small-signal angle stability: It is the ability of the system to
remain in synchronism when subjected to small disturbances. Small disturbances
can be small load changes like switching on or off of small loads, line tripping,
small generators tripping etc.

 Large-disturbance or Transient angle stability: It is the ability of the system to
remain in synchronism when subjected to large disturbances. Large disturbances
can be faults, switching on or off of large loads, large generators tripping etc.
The time domain of interest in case of large-disturbance as well as small-
disturbance angle stability is anywhere between 0.1- 10 s.

Figure 1.1: Classification of Power System Stability

Page | 6

 Voltage stability: It is the ability of the system to maintain steady state voltages at
all the system buses when subjected to a disturbance. If the disturbance is large then
it is called as large-disturbance voltage stability and if the disturbance is small it is
called as small-disturbance voltage stability.

Unlike angle stability, voltage stability can also be a long term phenomenon. In
case voltage fluctuations occur due to fast acting devices like induction motors,
power electronic drive, HVDC etc then the time frame for understanding the
stability is in the range of 10-20 s and hence can be treated as short term
phenomenon. On the other hand if voltage variations are due to slow change in
load, over loading of lines, generators hitting reactive power limits, tap changing
transformers etc then time frame for voltage stability can stretch from 1 minute to
several minutes.
The main difference between voltage stability and angle stability is that voltage
stability depends on the balance of reactive power demand and generation in the
system where as the angle stability mainly depends on the balance between real
power generation and demand.

 Frequency stability: It refers to the ability of a power system to maintain steady
frequency following a severe disturbance between generation and load. It depends
on the ability to restore equilibrium between system generation and load, with
minimum loss of load.

Frequency instability may lead to sustained frequency swings leading to tripping of
generating units or loads. During frequency excursions, the characteristic times of
the processes and devices that are activated will range from fraction of seconds like
under frequency control to several minutes, corresponding to the response of
devices such as prime mover and hence frequency stability may be a short-term
phenomenon or a long-term phenomenon. Though, stability is classified into rotor
angle, voltage and frequency stability they need not be independent isolated events.
A voltage collapse at a bus can lead to large excursions in rotor angle and
frequency. Similarly, large frequency deviations can lead to large changes in
voltage magnitude.

Page | 7






Chapter 2

Overview of Synchronous
Generator

Page | 8


In this chapter, we give a brief overview on the synchronous generator. The
commercial birth of the alternator (synchronous generator) can be dated back to August 24,
1891. On that day, the first large-scale demonstration of transmission of ac power was
carried out. The transmission extended from Lauffen, Germany, to Frankfurt, about 110
miles away.

Synchronous machines are principally used as alternating current (AC) generators.
According to the arrangement of the field and armature windings, synchronous machines
may be classified as rotating-armature type where the armature winding is on the rotor and
the field system is on the stator or rotating-field type where the armature winding is on the
stator and the field system is on the rotor.

Rotor winding (field winding) (produce rotor magnetic field) is excited by DC supply by
means of slip rings and brushes. The rotor is turned by prime mover producing a rotating
magnetic field. The rotating magnetic field produces three phase sets of voltages within the
stator coil (Armature winding) separated by an angle 120
0
in space. Synchronous means
that the electrical frequency produced is locked with the mechanical rate of rotation of the
generator. Two types of rotor constructions are used in asynchronous alternator: (i) Salient
pole type for low and medium speed alternators. Alternators featuring this type of rotor are
large in diameters and short in axial length and (ii) Cylindrical type for high speed
alternators, especially in turbo alternators. This type of rotor consists of a smooth and solid
steel cylinder having slots along its outer periphery. Field windings are placed in these
slots.

In the cylindrical rotor synchronous machine, the air gap is uniform. The pole structure of
the rotor of a salient pole machine makes the air gap highly non-uniform. The axis along
the axis of the rotor is called the direct or the d axis. The axis perpendicular to d axis is
known as the quadrature or q axis. The direct axis flux path involves two small air gaps and
is the path of the minimum reluctance. The synchronous reactance associated with the
direct axis is therefore and is known as direct axis synchronous reactance ??????
�. The
minimum synchronous reactance is called quadrature axis synchronous reactance ??????
�.




Figure 2.1: Cross-sectional view of salient-pole synchronous machine

Page | 9


Synchronous Machine under Short Circuit Conditions




The concept of Subtransient, Transient and Steady State arises in case of fault in an
Alternator. Let us assume a sudden short circuit in three phase of alternator.













 When the alternator is short-circuited, the currents in all the three-phases rise
rapidly to a high value of about 10 to 18 times of full load current, during the first
quarter cycle. The flux crossing the air gap is large during a first couple of cycles.
The reactance during these first two or three cycle is least and the short circuit
current is high. This reactance is called subtransient reactance and is denoted by
X”. The first few cycles come under subtransient state.

 After a first few cycles, the decrement in the r.m.s. value of short circuit current is
less rapid than the decrements during the first few cycles. This state is called the
Transient State and the reactance in this state is called transient reactance X’.
The circuit breaker contacts separate in the transient state.

 Finally the transient dies out and the current reaches a steady sinusoidal state called
the Steady State. The reactance in this state is called steady state reactance ??????
�.
Since the short circuit current of the alternator lags behind the voltage by 90 degree,
the reactance involved are direct axis reactance.






Figure 2.2: Symmetrical Short Circuit current in Alternator

Page | 10



To understand the behavior of an alternator under transient conditions, the armature and
field resistance is assumed to be negligibly small. Thus, constant flux linkage theorem can
be applied. As per this theorem, in purely inductive circuit, the total flux linkage cannot be
changed instantaneously at the time of any disturbance. Now, if all the three phases of
unloaded alternator with normal excitation are suddenly short circuited there will be short-
circuit current flows in the armature. As the resistance is assumed to be zero, this current
will lag behind the voltage by 90 degree and the mmf produced by this current will be
along the d-axis.

First conclusion is that this current will be affected by d-axis parameters ??????
�, ??????
�

and ??????
�
′′

only.

Further, there will be demagnetizing effect of this current, but as the flux linkage with field
cannot change the effect of demagnetizing armature m.m.f. must be counterbalanced by a
proportional increase in the field current. This additional induced component of field
current gives rise to greater excitation under transient state and results in more short
circuits as compared to the steady state short circuit current.

If field poles are provided with damper bars, then at the instant of three phase short circuit,
the demagnetizing armature m.m.f. induces currents in damper bars, which, in turn,
produces field in the same direction as the main field and hence at this instant, the
excitation further increases and gives rise to further increase in short circuit armature
current. This is for a very short duration, normally 3 to 4 cycles and this period is known as
sub-transient period. Since the field voltage is constant, there is no additional voltage to
sustain these increased excitations during sub transient or transient period.

Consequently the effect of increased field current decreases with a time constant
determined by the field and armature parameters and accordingly the short circuit armature
current also decays with the same time constant.

The reactances offered by the machine during sub transient period are known as sub
transient reactances. Along the direct axis, it is direct axis sub transient reactance, ??????
�
′′
and
along the quadrature axis, it is quadrature axis sub-transient reactance ??????
�
′′
.

Page | 11


Figure 2.3: Flux paths for various reactances of a salient-pole synchronous machine

Page | 12







Chapter 3

Analysis of Swing Equation

Page | 13

Derivation of Swing Equation



The differential equation that relates angular momentum, its acceleration and the rotor
angle is called swing equation
[3]
. Solution of swing equation will show how the rotor
angle changes with respect to time following a disturbance.

Kinetic energy of rotor ??????�=(1/2)�??????
�
2
, where �= moment of inertia in ��/�
2
and
??????
�= angular speed in mechanicalrad/sec. Angular momentum ??????=�??????
�; so ??????�=
(1/2)????????????
�.The value of ?????? changes as the angular velocity changes.However, assuming
??????
� is relatively constant during thetransient analysis time, ?????? may be considered constant.
We know when torque is applied to a body, the body experiences angular acceleration and
the torque is equal to the product of angular acceleration and moment of inertia: �=�??????,
where � is the net torque and ?????? is the angular acceleration. The angular acceleration, a,
may be expressed in terms of ??????
�, a mechanical angle which is measured from a non-
rotating reference frame. So,

??????=�
2
??????
�/��
2
.

So,

�=�(�
2
??????
�/��
2
).

The net torque, which produces acceleration,' is the algebraic difference of mechanical
torque applied and electrical torque output =�
�−�
�=�
??????. Thus, the equation of motion
of the machine is governed by �(�
2
??????
�/��
2
)=�
�−�
�=�
??????; �
� is the mechanical
torque supplied by the prime mover in N-m, �
� is the electrical torque output of the
alternator in N-m and ??????
� is the angular position of the rotor in mechanical rad.
Multiplying both sides with ??????
�, we get

�??????
�
�
2
??????
�
��
2
=�
�??????
�−�
�??????
�=�
????????????
�,�� ??????
�
2
??????
�
��
2
=??????
�−??????
�=??????
??????

where ??????
� , ??????
� and ??????
?????? respectively are the mechanical, electrical and accelerating power in
MW. Neglecting the losses, the difference between the mechanical and electrical torque
gives the net accelerating torque �
?????? . In the steady state, the electrical torque is equal to
the mechanical torque, and hence the accelerating power will be zero. During this period
the rotor will move at synchronous speed ??????
� in mechanical rad/s.

As the rotor is continuously rotating at synchronous speed in steady state ??????
� will also be
continuously varying with respect to time. To make the angle θ
mconstant in steady state
we can measure this angle with respect to a synchronously rotating reference instead of a
stationary reference. Hence, we can write the angular position ??????
m with respect to the
synchronously rotating frame, we define ??????
�=??????
���+??????
�.

�??????
�
��
=??????
��+
�??????
�
��
���
�
2
??????
�
��
2
=
�
2
??????
�
��
2

Page | 14

So, the equation, which describes the behavior of the rotor dynamics, is the swing
equation, which is as follows:

??????
�
2
??????
�
��
2
=??????
�−??????
�=??????
??????

In order to convert them in to electrical radians and electrical radians per second
respectively we have to take the number of poles (??????) of the synchronous machine rotor
into consideration. Hence, the electrical angle and electrical speed can be represented as
??????=(??????/2)??????
� elect. rad and ??????=(??????/2)??????
� elect. rad per second.Rather than giving the
moment of inertia of a synchronous machine for dynamic studies, it is more convenient to
use per unitized quantity called inertia constant � defined as the stored kinetic energy per
unit KVA of machine rating, the swing equation can be represented as in pu system

2�
??????
�
�
2
??????
��
2
=??????
�−??????
�=??????
??????; where �=
1
2
????????????
�

2�
??????
�
�
2
??????
��
2
=??????
�−??????
�=??????
??????;or
�
2
??????
��
2
=
??????�
�
�
[??????
�−??????
�]=
??????�
�
�
[??????
�−??????
�????????????sin??????]

In can be observed from above equation that, if ??????
�=??????
�????????????sin??????, then there will be no
speed change and there will be no angle change. But, if ??????
�≠??????
�????????????sin??????due to
disturbance in the system then either the speed increase or decrease with respect to time.
Let us take the case of ??????
�>??????
�????????????sin??????, there is more input mechanical power than the
electrical power output. In this case, as the energy has to be conserved difference
between the input and output powers will lead toincrease in the kinetic energy of the
rotor and speed increases.

Similarly, if ??????
&#3627408474;<??????
&#3627408474;????????????sin?????? then, the input power is less than the required electrical
power output. Again the balance power, to meet the load requirement, is drawn from the
kinetic energy stored in the rotor due to which the rotor speed decreases.

Since, the mechanical power input P
mand the maximum power output of the generator
max P
max are known for a given system topology and load, we can find the rotor angle δ
from

??????=&#3627408480;??????&#3627408475;
−1
[
??????
&#3627408474;
??????
&#3627408474;????????????
] &#3627408476;&#3627408479; ??????−&#3627408480;??????&#3627408475;
−1
[
??????
&#3627408474;
??????
&#3627408474;????????????
]


It can be observed from the above that rotor angle δhas two solutions at which ??????
&#3627408474;=
??????
&#3627408474;????????????&#3627408480;??????&#3627408475;??????. Now, let us label


??????
0=&#3627408480;??????&#3627408475;
−1
[
??????
&#3627408474;
??????
&#3627408474;????????????
] &#3627408462;&#3627408475;&#3627408465; ??????
&#3627408474;????????????=??????−&#3627408480;??????&#3627408475;
−1
[
??????
&#3627408474;
??????
&#3627408474;????????????
]=??????−??????
0

Page | 15

There is a very important implication of the two solutions (??????
0,??????
&#3627408474;????????????) of swing equation
on the stability of the system that can be written as:

&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
&#3627408474;−??????
&#3627408474;????????????sin??????


This can be clearly understood from what is called as swing curve or ??????−??????curve.




The swing curve shows the plot of electrical power output ??????
&#3627408466;with respect to the rotor
angle ??????. As can be observed the curve is a sine curve with ?????? varying from 0 to ??????. The
maximum power output of the generator ??????
&#3627408474;???????????? occurs at an angle ??????=π/2. On the same
curve the input mechanical power ??????
&#3627408474; can also be represented.

Since P
m is assumed to be constant and does not vary with respect toδ, a straight line is
drawn which cuts the P−δcurve at points A and B. It can be seen that at points A and B
the mechanical power input P
mis equal to P
e. The rotor angles at point A and B are δ
0and
δ
max, the solutions of equation P
m=P
maxsinδ.

Let us try to understand which of the solutions (??????
0, ??????
&#3627408474;????????????) leads to a stable operation. Let
us take first point A (corresponding rotor angle is δ
0). If we perturb the rotorangle δ
0 by a
small positive angle ∆δ, so that the new operating point is at C, thenelectrical power
output will also increase to P
maxsin(δ
0+∆δ).

Since, in steady state P
m=P
maxsinδ
0, after perturbation P
m<P
maxsin(δ
0+∆δ), for a
positive value of ∆δ. Now the output electrical power will become more than the input
mechanical power and hence the rotor starts decelerating due to which the angle δwill be
pulledback to the point A. But since the rotor has certain inertia it cannot stop at the point
Aand decelerate further due to which the angle δmoves to the point say D.

At point D, P
m>P
maxthat is input mechanical power will become more than the electrical
power output and hence the rotor starts accelerating.

Figure 3.1: Swing Curve or P-δ Curve

Page | 16

Again the rotor angle δ starts decreasing and will reach point A but due to inertia it cannot
stop there and it will again move to point C. This phenomenon repeats itself indefinitely if
there is no damping to the rotor oscillations.

Now take the case of second operating point δ
max. Again if we perturb δ
max by a positive
angle ∆δ then the electrical power output will be P
maxsin (δ
max+∆δ).

In this case however, unlike the earlier case P
m>P
maxsin(δ
max+∆δ) and the rotor
starts accelerating due to which the angle will increase further and this will lead to further
decrease in the electrical power output.

Hence, for a small positive perturbation in the rotor angle at the operating point B leads to
continuous increase in the speed of the rotor there by leading to unstable operation of the
generator.

This discussion leads us to an important conclusion that out of the two operating points A
and B, with rotor angles δ
0 and δ
max, operating point A is stable and operating point B is
unstable for small disturbances. Hence, point A is called stable equilibrium point and
operating point B is called as unstable equilibrium point. Similarly, δ
0 is a stable steady
state rotor angle and δ
max is an unstable rotor angle.

Page | 17

Linearization of Swing Equation


From ??????−??????curve is can be observed that if we perturb the rotor angle by a small angle
∆?????? then the rotor angle starts oscillates. If the rotor angle perturbation ∆?????? is very small
then we can assume that the system behaves in a linear fashion around the operating
point. So,


&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
(??????
0+∆??????
0)
&#3627408465;&#3627408481;
2
=??????
&#3627408474;−??????
&#3627408474;????????????sin(??????
0+∆??????
0)


Expanding &#3627408480;??????&#3627408475;(??????
0+∆??????), we get


&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
??????
0
&#3627408465;&#3627408481;
2
+
&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
∆??????
&#3627408465;&#3627408481;
2
=??????
&#3627408474;−??????
&#3627408474;????????????(&#3627408480;??????&#3627408475;??????
0&#3627408464;&#3627408476;&#3627408480;∆??????+&#3627408464;&#3627408476;&#3627408480;??????
0&#3627408480;??????&#3627408475;∆??????)


Since, ∆??????is a very small perturbation angle, we can approximately write &#3627408464;&#3627408476;&#3627408480;∆??????≅1 and
&#3627408480;??????&#3627408475;∆??????≅∆??????. With these approximations substituted in above, we can get


&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
??????
0
&#3627408465;&#3627408481;
2
+
&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
∆??????
&#3627408465;&#3627408481;
2
=??????
&#3627408474;−??????
&#3627408474;????????????&#3627408480;??????&#3627408475;??????
0−??????
&#3627408474;????????????&#3627408464;&#3627408476;&#3627408480;??????
0∆??????


Since at initial operating point,

&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
??????
0
&#3627408465;&#3627408481;
2
=??????
&#3627408474;−??????
&#3627408474;????????????&#3627408480;??????&#3627408475;??????
0=0

So,

&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
∆??????
&#3627408465;&#3627408481;
2
+??????
&#3627408474;????????????&#3627408464;&#3627408476;&#3627408480;??????
0∆??????=0;

or

&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
∆??????
&#3627408465;&#3627408481;
2
+??????
?????? ∆??????=0


Almost all synchronous generators are equipped with damper windings (also called
amortisseur windings).The reason for providing damper windings is to damp oscillations
started by aperiodic shocks, such as short circuits or switching or load variation.

Page | 18

When damper action is neglected, the accelerating power is the difference between the
mechanical power input and the electric power, i.e. ??????
??????=??????
&#3627408474;−??????
&#3627408466;.


When damper winding is considered, another term is included. Thus, ??????
??????=??????
&#3627408474;−??????
&#3627408439;−??????
&#3627408466;,
where ??????
&#3627408439; is the damping torque, where


P
&#3627408439;=&#3627408439;
&#3627408465;∆??????
&#3627408465;&#3627408481;



where, D is the damping coefficient. After including damping torque, the swing equation
can be written as


&#3627408443;
??????&#3627408467;
&#3627408480;
&#3627408465;
2
∆??????
&#3627408465;&#3627408481;
2
+&#3627408439;
&#3627408465;∆??????
&#3627408465;&#3627408481;
+??????
?????? ∆??????=0


Applying Laplace transformation we get


&#3627408480;
2
+&#3627408439;
??????&#3627408467;
&#3627408480;
&#3627408443;
&#3627408480;+
??????&#3627408467;
&#3627408480;
&#3627408443;
??????
?????? =0

Comparing above with the characteristic equation of a standard second order system, we
get

&#3627408480;
2
+2????????????
&#3627408475;&#3627408480;+??????
&#3627408475;
2
=0;&#3627408484;ℎ&#3627408466;&#3627408479;&#3627408466; ??????
&#3627408475;=√
??????&#3627408467;
&#3627408480;
&#3627408443;
??????
?????? &#3627408462;&#3627408475;&#3627408465; ??????=
&#3627408439;
2

??????&#3627408467;
&#3627408480;
&#3627408443;??????
??????


Where, ??????
&#3627408475; is the natural frequency and ?????? is damping ratio of the oscillations of the
system.

Page | 19

Numerical Solution of Swing Equation



From swing equation, we get

??????
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????=??????
&#3627408474;−??????
&#3627408466;=??????
&#3627408474;−??????
&#3627408474;????????????&#3627408480;??????&#3627408475;??????

A direct solution of this equation is very difficult to achieve. For small disturbances, as is
covered in steady state stability studies, it is possible to linearize the equation. This allows
us to obtain a condition, satisfaction of which ensures the stability of the system.

Such a simplistic approach is not applicable to the solution of swing equation for large
disturbances as in the case of transient stability. The method of solution is more involved.
It requires use of such tools such as numerical techniques
[5]
to approximate the non-linear
system to a discrete form with arbitrarily close approximation. Thus, the general approach
we use consists of taking into account the disturbances in the system and solving the
equation in the presence of said disturbances.

Euler’s method: Let a non-linear differential equation of the form &#3627408465;&#3627408485;/&#3627408465;&#3627408481;=f(&#3627408485;) exist,
where f(&#3627408485;) is a nonlinear function of &#3627408485;, we can find the solution of said equation through
Euler’s method
[6]
. The idea is to integrate the equation between the time instants (&#3627408481;
0,&#3627408481;
&#3627408467;)
with an initial value of &#3627408485;=&#3627408485;
0
at &#3627408481;=&#3627408481;
0 with a small time step &#3627408481;=??????&#3627408481; that can be
expressed as

&#3627408485;
1
=&#3627408485;
0
+[
&#3627408465;&#3627408485;
&#3627408465;&#3627408481;
|
??????=??????
0
]??????&#3627408481;

Here, &#3627408485;
1
is the variable value at time instant &#3627408481;=&#3627408481;
0+??????&#3627408481;. This process has to be repeated
iteratively until the final time t
f is reached. The equation is nothing but Taylor series
expansion of the variable &#3627408485; around the operating point (&#3627408481;
0, &#3627408485;
0) with higher order terms
discarded. Since, the higher terms are discarded the above equation may introduce error.

Hence, to reduce the error modified Euler’s method can be used. Modified Euler’s method
has two steps: predictor step and corrector step. These steps in generalized form are as
follows:

&#3627408485;
&#3627408477;
&#3627408472;+1
=&#3627408485;
&#3627408472;
+[
&#3627408465;&#3627408485;
&#3627408465;&#3627408481;
|
??????=??????
??????
]??????&#3627408481;

and

&#3627408485;
&#3627408464;
&#3627408472;+1
=&#3627408485;
&#3627408472;
+
1
2
[
&#3627408465;&#3627408485;
&#3627408465;&#3627408481;
|
??????=??????
??????
??????+1
+
&#3627408465;&#3627408485;
&#3627408465;&#3627408481;
|
??????=??????
??????
]??????&#3627408481;


where &#3627408485;
&#3627408477;
1
and &#3627408485;
&#3627408464;
1
is the predicted and the corrected value of the variable &#3627408485; at the time
instant &#3627408481;=&#3627408481;
0+??????&#3627408481;. The time step ??????&#3627408481; has to be carefully chosen, the smaller the values the
better the accuracy.

Page | 20

Ideally ??????&#3627408481; should approach zero but when the number of steps required to integrate above
equation between the limits(&#3627408481;
0,&#3627408481;
&#3627408467;)will become infinite. Hence there is a trade-off
between the accuracy and the number of steps required.

The Euler’s method can be applied to get a solution of swing equation between the time
limits (&#3627408481;
0,&#3627408481;
&#3627408467;).




Figure 3.2: Graphical Representation of Modified Euler's Method

For the solution of the swing equation, let us assume the initial conditions as ??????
0=0,
??????
0=??????
??????&#3627408475;??????&#3627408481;????????????&#3627408473; &#3627408462;&#3627408475;&#3627408465; ??????
??????=??????
&#3627408474;−??????
&#3627408474;????????????&#3627408480;??????&#3627408475;??????
0. For the n
th
iteration using Euler’s method, the
steps are ??????
0=??????
&#3627408475;−1+ Δt (
1
??????
)??????
????????????
&#3627408475;−1 and ??????
&#3627408475;=??????
&#3627408475;−1+??????
&#3627408475;Δt.

For solving using Modified Euler’s method, we use

??????
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????;
&#3627408465;??????
&#3627408465;&#3627408481;
=?????? &#3627408462;&#3627408475;&#3627408465;
&#3627408465;??????
&#3627408465;&#3627408481;
=
??????
??????
??????


For this set of equations, at any point in the iterative solution, the following sequence of
operations take place ??????
&#3627408477;=??????+ ??????&#3627408481; [(1/??????)??????
??????(??????)] and ??????
&#3627408477;=??????+??????
&#3627408477;??????&#3627408481;. Thus, we get the
predicted values of the functions after a time interval Δt .



t
x=f(t)
(x
0
, t0)
t0 t
1

Δt
??????&#3627408481; f(x
0
, t0)

tf

Exact
Solution
Approximate
Solution

Page | 21

Now, we shall get the values of the slope at the predicted values of the function &#3627408465;??????/&#3627408465;&#3627408481;=
??????
&#3627408477;. We can omit the step, since we already know the value of ω
p and
&#3627408465;??????/&#3627408465;&#3627408481;=[(1/??????)??????
??????(??????)].

Thus we have the slopes at the end of the time interval. Now, in one go we shall evaluate
the corrected values of δ and ω by averaging their slopes –

??????
&#3627408464;=??????+(1/2)(??????+??????
&#3627408477;)??????&#3627408481;;

??????
&#3627408464;=??????+ ??????&#3627408481; [(1/2??????)[??????
??????(??????)+??????
??????(??????
&#3627408477;)]]

These values are updated values for the next iteration.

Page | 22






Chapter 4

Damping and Damper Windings

Page | 23

Almost all synchronous motors, synchronous condensers, and synchronous converters and
many salient-pole synchronous generators are equipped with damper windings (also called
amortisseur windings). In a few installations damper windings have been designed with
regard to their effect upon system stability, but in most installations they have been
provided for other purposes.

Damper windings are not used on turbo generators, but the solid steel rotor cores of such
machines provide paths for eddy currents and thus produce the same effects as dampers.
Even in machines with laminated salient poles and no dampers, eddy currents in the pole
faces have a very small damping effect.


 Purposes of damper windings


The chief reasons for providing damper windings on salient-pole synchronous machines
are the following:

 To provide starting torque for synchronous motors, condensers, and converters.

 To suppress hunting. This was the first application of dampers and the one for
which they were named. Generators driven by reciprocating engines tend to hunt
because of the pulsating torque of the engines. Motors driving loads of pulsating
torque, such as compressors, likewise tend to hunt. Hunting from these causes is
"forced" hunting. There is also a spontaneous hunting that occurs when
synchronous machines are connected together by circuits having a high ratio of
resistance to reactance. This condition is likely to occur if the frequency is low, if
the conductor size is small, or if the inductive reactance is partially cancelled by use
of series capacitors. Both forced and spontaneous types of hunting are greatly
decreased in amplitude by low-resistance dampers.

 To damp oscillations started by aperiodic shocks, such as short circuits or
switching. A possible cause of oscillations between generators in the same station is
inequality of speeds of response of their several exciters. Oscillations from any of
these causes are damped more rapidly if damper windings are added.

 To prevent distortion of voltage wave shape caused by unbalanced loading; in other
words, to suppress harmonics. Dampers have been installed in some generators
which did not originally have them in order to decrease the harmonic voltages to
safe values.

 To balance the terminal voltages of the several phases during unbalanced loading,
that is, to decrease the negative-sequence voltage. Damper windings decrease the
negative-sequence reactance. The lower the negative-sequence impedance, the
lower is the negative-sequence voltage.

 To provide a braking torque on a generator during an unsymmetrical fault, and thus
to reduce the accelerating torque during the fault.

Page | 24

 To provide additional torque for synchronizing generators, especially those which
are synchronized automatically. They also help to pull the generator back into step
after synchronism has been lost because of a fault.


 Effect of Damping on Stability:


Positive-sequence damping results from the torque caused by interaction of the damper
currents with the positive-sequence (forward rotating) magnetic field in the air gap. Except
during starting or pulling out of step, the slip of the rotor with respect to the positive
sequence field is low and oscillatory. Positive-sequence damping causes the oscillations of
the machine rotors, after an aperiodic shock that does not cause loss of synchronism, to
decrease in amplitude and ultimately to die out entirely. Low-resistance amortisseurs
greatly increase the amount of positive-sequence damping. Positive-sequence damping is
present both during and after a fault.

However, it is much more effective after clearing of a fault than during the fault because
the positive-sequence voltage and flux are greater after clearing. By absorbing energy from
the oscillation, positive-sequence damping may prevent a machine which has survived the
first swing from going out of step on the second or subsequent swings because of reduction
of field flux linkage.

The expression of damping coefficient
[2]
is given as:


&#3627408439;=??????
&#3627408481;
2
&#3627408480;??????[
(??????
&#3627408465;

−??????
&#3627408465;
′′
)&#3627408455;
&#3627408465;0
′′
(??????
&#3627408466; +??????
&#3627408465;

)
2
&#3627408480;??????&#3627408475;
2
?????? +
(??????
&#3627408478;

−??????
&#3627408478;
′′
)&#3627408455;
&#3627408478;0
′′
(??????
&#3627408466; +??????
&#3627408478;

)
2
&#3627408464;&#3627408476;&#3627408480;
2
?????? ]


Where,

 ??????
&#3627408481; = Voltage of Infinite Bus in per unit

 &#3627408480; = Slip in per unit =
1
360f


 ??????=2??????&#3627408467;

 ??????
&#3627408466; = External resistance in per unit in series with armature

 ?????? = Power Angle

 &#3627408455;
&#3627408465;0
′′
, &#3627408455;
&#3627408478;0
′′
= Open-circuit sub-transient time constants of d and q axis respectively

 ??????
&#3627408465;

, ??????
&#3627408478;

= Transient reactance of d and q axis respectively

 ??????
&#3627408465;
′′
, ??????
&#3627408478;
′′
= Sub-transient reactance of d and q axis respectively

Page | 25












Chapter 5

System Under Study

Page | 26

System under Study


The system under study consists of a turbo alternator connected to an infinite bus through
transformers and double circuit line as shown below:





If a three phase fault occurs at the mid-point of one of the lines, the machine will not be
able to supply the full amount of power with the full prime mover input being present and
with reduced output, the rotor will accelerates and the load angle increases. If it increases
indefinitely, the machine loses synchronism and the stability will be lost. Now, if the fault
can be cleared within a certain time, the machine can restore synchronism and the value of
rotor angle at that instant is called critical clearing angle, which is the edge of instability.

The synchronous generator considered is having the following parameters:
120 MW, 13.8 KV, 0.8 p.f. (lag), star connected

??????
&#3627408465;= 1.460 p.u. ??????
&#3627408465;

= 0.2960 p.u. ??????
&#3627408465;
′′
= 0.2000 p.u.
??????
&#3627408478;= 1.4000 p.u. ??????
&#3627408478;
′′
= 0.592 p.u ??????
&#3627408478;
′′
= 0.2000 p.u.
&#3627408455;
&#3627408465;0

= 4.4 s &#3627408455;
&#3627408465;0
′′
= 0.492 s &#3627408455;
&#3627408478;0
′′
= 0.07 s

Inertia constant ??????= 0.01878978. Saliency is negligible.

The ratings of both the transformers are identical. Each 150 KVA, 11.8 / 220 KV, 3 –
phase, 50 Hz with a vector group Yd 11. The leakage reactances of transformers 1 and 2
are 0.12 and 0.11 p.u. respectively both referred to alternator base. Each line is capable of
transmitting about 75 MW and its positive sequence reactance is 0.2 p.u. with reference to
alternator base.

Figure 5.1: Single Line Diagram of SMIB System

Page | 27







Chapter 6

Mathematical Formulation

Page | 28


Mathematical Approach


The swing equation considering damping is given as

M
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????−
&#3627408440;
&#3627408467; ??????
&#3627408481; &#3627408480;??????&#3627408475;??????
??????
−&#3627408439;
&#3627408465;??????
&#3627408465;&#3627408481;



M
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
+
&#3627408440;
&#3627408467; ??????
&#3627408481; &#3627408480;??????&#3627408475;??????
??????
+&#3627408439;
&#3627408465;??????
&#3627408465;&#3627408481;
=??????
??????


&#3627408462;
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
+&#3627408463;
&#3627408465;??????
&#3627408465;&#3627408481;
+&#3627408464;(??????)=d


Let &#3627408485;=?????? and y =
&#3627408465;??????
&#3627408465;&#3627408481;
. Then the above equation reduces to: &#3627408462;
&#3627408465;&#3627408486;
&#3627408465;&#3627408481;
+&#3627408463;&#3627408486;+&#3627408464;(&#3627408485;)=&#3627408465;

&#3627408465;&#3627408486;
&#3627408465;&#3627408481;
=
&#3627408465;−&#3627408463;&#3627408486;−&#3627408464;(&#3627408485;)
&#3627408462;


For an ideal un-damped system, &#3627408463;=0. To find the values of ??????(=&#3627408485;), we need to solve the
two differential equations i.e.

&#3627408486;=
&#3627408465;??????
&#3627408465;&#3627408481;


and

&#3627408465;&#3627408486;
&#3627408465;&#3627408481;
=
&#3627408465;−&#3627408463;&#3627408486;−&#3627408464;(&#3627408485;)
&#3627408462;


simultaneously, using numerical methods.

The value of X will vary depending on the location of fault and the response of the curve
will be plotted for different values of fault clearing time.

In pre-fault state the machine is delivering rated current at 0.8 pf to the infinite bus (at
rated voltage). Therefore, the voltage behind transient reactance
[1]
is given as:

(&#3627408440;)
2
=(??????
&#3627408481;)
2
+(2??????
&#3627408481;&#3627408444;
????????????
&#3627408465;

sin∅)+(&#3627408444;
????????????
&#3627408465;

)
2


&#3627408440;=1.464 p.u

Page | 29

Calculation of Damping Constant


For numerical solution of the swing equation, we take

&#3627408465;??????
&#3627408465;&#3627408481;
=
Δ??????
Δ&#3627408481;
and assume Δ&#3627408481;=0.01 sec

Thus for the given system, the damping power

&#3627408439;
&#3627408465;??????
&#3627408465;&#3627408481;
results to &#3627408439;
Δ??????
Δ&#3627408481;


Damping power

??????
&#3627408439;=
2??????
&#3627408481;
2
??????&#3627408467;∆??????
360&#3627408467;∆&#3627408481;
[
(??????
&#3627408465;

−??????
&#3627408465;
′′
)&#3627408455;
&#3627408465;0
′′
(??????
&#3627408466;+ ??????
&#3627408465;

)
2
&#3627408480;??????&#3627408475;
2
?????? +
(??????
&#3627408478;

−??????
&#3627408478;
′′
)&#3627408455;
&#3627408478;0
′′
(??????
&#3627408466;+ ??????
&#3627408478;

)
2
&#3627408464;&#3627408476;&#3627408480;
2
?????? ]
















??????
&#3627408466;=0.12+(0.2||0.2)+0.11=0.33 &#3627408477;&#3627408482;


Assume the smallest time interval for numerical solution ∆&#3627408481;=0.01 &#3627408480;&#3627408466;&#3627408464;.



??????
&#3627408439;=
2??????
&#3627408481;
2
πf∆δ
360&#3627408467;∗0.01
[
(0.296−0.2 )∗0.492
(0.33+0.296)
2
&#3627408480;??????&#3627408475;
2
?????? +
(0.592−0.2 )∗0.07
(0.33+0.592)
2
&#3627408464;&#3627408476;&#3627408480;
2
?????? ]


= ??????
&#3627408481;
2
(0.178 &#3627408480;??????&#3627408475;
2
δ + 0.0172&#3627408464;&#3627408476;&#3627408480;
2
δ) ∆δ


Damping coefficient &#3627408439;=??????
&#3627408481;
2
(0.178 &#3627408480;??????&#3627408475;
2
δ + 0.0172&#3627408464;&#3627408476;&#3627408480;
2
δ) = ??????
&#3627408481;
2
(0.097 + 0.0804cos2δ)

Figure 6.1: Single Line Diagram with Per-Unit Reactances

Page | 30

Critical Clearing Angle


The critical clearing angle is defined as the maximum change in the load angle curve
before clearing the fault without loss of synchronism. In other words, when the fault occurs
in the system the load angle curve begin to increase, and the system becomes unstable. The
angle at which the fault becomes clear and the system becomes stable is called critical
clearing angle.

Let ??????
1, ??????
2 and ??????
3 be respectively the reactance of the system before fault, during fault, and
after the fault is cleared. Then the electrical power that can be transferred is:

P1 =
&#3627408440;
&#3627408467;
??????
&#3627408481;

??????
1
&#3627408480;??????&#3627408475;δ = P
m1&#3627408480;??????&#3627408475;δ ---------- Before fault

P2 =
&#3627408440;
&#3627408467;
??????
&#3627408481;
??????
2
&#3627408480;??????&#3627408475;δ = P
m2&#3627408480;??????&#3627408475;δ ---------- During occurrence of fault

P3 =
&#3627408440;
&#3627408467;
??????
&#3627408481;
??????
3
&#3627408480;??????&#3627408475;δ = P
m3&#3627408480;??????&#3627408475;δ ---------- After fault is cleared


Here the curves for P1, P2, P3 are represented by A, B
and C respectively.

Using the principle of Equal Area Criteria i.e. net
accelerating area equals to zero, we get:

∫ ??????
??????
????????????????????????
??????0
&#3627408465;?????? =0


∫??????
??????
????????????
??????0
&#3627408465;?????? + ∫ ??????
??????
????????????????????????
????????????
&#3627408465;?????? =0


∫(??????
&#3627408480; −??????
2 )
????????????
??????0
&#3627408465;?????? + ∫(??????
&#3627408480; −??????
3 )
????????????????????????
????????????
&#3627408465;?????? =0

Solving this, we get
??????
&#3627408464;=&#3627408464;&#3627408476;&#3627408480;
−1
[
??????
&#3627408480;(??????
&#3627408474;????????????−??????
0)+??????
&#3627408474;3 &#3627408464;&#3627408476;&#3627408480;??????
&#3627408474;????????????−??????
&#3627408474;2 &#3627408464;&#3627408476;&#3627408480;??????
0
??????
&#3627408474;3 −??????
&#3627408474;2
]


Figure 6.2: P-δ curve showing Accelerating and De-
accelerating Areas

Page | 31

Initial power angle is expressed as:
??????
0=&#3627408480;??????&#3627408475;
−1
(
??????&#3627408480;
??????
&#3627408474;1
)
??????
180


Maximum swing angle is expressed as:
??????
&#3627408474;????????????=180
0
−&#3627408480;??????&#3627408475;
−1
(
??????
&#3627408480;
??????
&#3627408474;3
)

Page | 32

Critical Clearing Time

The Critical Clearing Time is the maximum time during which a disturbance can be
applied without the system losing its stability. The aim of its calculation is to determine the
characteristics of protections required by the power system. Consider the undamped swing
equation:

M
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????−??????
&#3627408466;

During the period of occurrence of fault, ??????
&#3627408466;=0. Hence the equation reduces to

M
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????

2&#3627408443;
??????
&#3627408480;
&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=??????
??????

&#3627408465;
2
??????
&#3627408465;&#3627408481;
2
=
??????
&#3627408480;
2&#3627408443;
??????
??????

&#3627408465;??????
&#3627408465;&#3627408481;
=∫
??????
&#3627408480;
2&#3627408443;
??????
??????&#3627408465;&#3627408481;
&#3627408481;
&#3627408476;
=
??????
&#3627408480;
2&#3627408443;
??????
??????&#3627408481;

δ=∫
??????
&#3627408480;
2&#3627408443;
??????
??????&#3627408481;&#3627408465;&#3627408481;
&#3627408481;
&#3627408476;
=
??????
&#3627408480;
4&#3627408443;
??????
??????&#3627408481;
2
+ ??????
0

Now δ = ??????
&#3627408438; and t = &#3627408481;
&#3627408438;
δ
C−δ
0=
??????
&#3627408480;
4&#3627408443;
??????
??????&#3627408481;
2


&#3627408481;
&#3627408464;= √
2&#3627408443;(??????
&#3627408464;−??????
0)
πf??????
&#3627408480;


If the time of clearance of the fault is more than the critical clearing time, the system does
not regain its stability. If damping is absent and fault is cleared within the critical clearing
time, the system will exhibit sustained oscillation. The clearing time of is derived based on
the assumption that the electrical power ??????
&#3627408466; becomes zero during the fault
[3]
. This need not
be the case always. In that even we have to resort to finding the clearing time using the
numerical integration of the swing equation.

Page | 33

Effect of Switching Operation on Stability


Consider a generator supplying power to an infinite bus over two transmission lines
operating in parallel.









The power angle curve for the two lines operating in
parallel is represented by A. Curve B represents the
power-angle diagram when a three phase fault occurs
on any one of the lines. Curve C represents the post-
fault power-angle diagram when the faulty line is
removed after clearance of the fault.
The value of reactance of the system during
occurrence of fault and after clearance depends on the
distance of fault location from the alternator. In the
load-angle characteristics, the value of ??????
&#3627408464; and hence
the stability limit will be determined by the position of
the fault since the electrical power is dependent on
reactance.

The load angle v/s time characteristics will
however depend on both the reactance and the fault clearing time. If &#3627408481;
&#3627408467; represents the fault
clearing time, then output electrical power is ??????
2 for t <&#3627408481;
&#3627408467;and is ??????
3 for t ≥ &#3627408481;
&#3627408467;. Again the
magnitude of ??????
2 and ??????
3 is dependent on the distance of location of the fault from the
alternator end.

??????
2=
&#3627408440;
&#3627408467;
??????
&#3627408481;
??????
2
sin?????? &#3627408462;&#3627408475;&#3627408465; ??????
3=
&#3627408440;
&#3627408467;
??????
&#3627408481;
??????
3
sin??????

This information will be used to determine whether a system will regain back its stability
for different fault locations and different values of fault clearing time.





Figure 6.3: Single Line diagram of system under study
Figure 6.4: P-δ curve for Pre, During and Post
Fault

Page | 34

Fault occurs on the Bus-Bar






 Reactance before fault: X
1 = (0.296 + 0.12) + (0.2||0.2) + 0.11 = 0.626 p.u

 Reactance during fault: ??????
2 = Infinity

During occurrence of fault, the bus-bar is protected by the relays which are again
restored after clearance of fault. Hence,

 Reactance after clearance of fault: X
3 = (0.296 + 0.12) + (0.2||0.2) + 0.11 = 0.626
p.u


Fault occurs on any line close to the Bus-Bar







 Reactance before fault: X
1= (0.296 + 0.12) + (0.2||0.2) + 0.11 = 0.626 p.u

As the fault occurs very close to the bus-bar, the protective relays of the bus-bar trips
during fault, preventing any power flow from source to load.
 Reactance during fault: ??????
2 = Infinity

 Reactance after clearance of fault: ??????
3 = (0.296 + 0.12) + 0.2 +0.11 = 0.726 p.u


Figure 6.5: Single-Line diagram for fault on Bus-Bar
Figure 6.6: Single Line diagram for fault on a line very close to Bus-Bar

Page | 35


Fault occurs at 1/4
th
Distance from Alternator









 Reactance before fault: X
1 = (0.296 + 0.12) + (0.2||0.2) + 0.11= 0.626 p.u














 Reactance during fault: X
2 = 5.8065 p.u

 Reactance after clearance of fault:??????
3 = (0.296+0.12) + 0.2 + 0.11 = 0.726 p.u


Figure 6.7: Single Line diagram for fault on T.L at 1/4th distance
Figure 6.8: Step wise calculation of the reactance during fault

Page | 36


Fault occurs on the middle of a line












 Reactance before fault: X
1 = (0.296 + 0.12) + (0.2||0.2) + 0.11 = 0.626 p.u

 The step-by-step calculation of Reactance after fault, X
2 using ∆ - Y and Y-∆ is
shown below: X
2 = 3.6084 p.u























 Reactance after clearance of fault: X
3 = (0.296+0.12) + 0.2 + 0.11 = 0.726 p.u


Figure 8. 9: Single Line diagram for fault at the middle of T.L
Figure 8.10: Step wise calculation of the reactance during fault

Page | 37


Fault occurs at 3/4
th
Distance from Alternator

 Reactance before fault: X
1= (0.296 + 0.12) + (0.2||0.2) + 0.11 = 0.626 p.u
















 Reacance during fault: X
2= 3.587 p.u

 Reactance after clearance of fault: X
3 = (0.296+0.12) + 0.2 +0.11 = 0.726 p.u


Figure 8.11: Single Line diagram for fault on T.L at 3/4th distance
Figure 8.12: Step wise calculation of the reactance during fault

Page | 38






Chapter 7

Developing the MATLAB
Programs

Page | 39

Algorithms

 Algorithm to plot Power-Angle Curve and calculate Critical Clearing Angle


























START
INPUT
Pm, E, V
X1, X2, X3, H, f



Pe1max = E*V/X1
Pe2max = E*V/X2
Pe3max = E*V/X3

delta = 0
delta< pi

Pe1 =Pe1max*sin(delta)
Pe2 = Pe2max*sin(delta)
Pe3 = Pe3max*sin(delta)


d0=&#3627408532;??????&#3627408527;
−&#3627409359;
(Pm/Pe1max)
dmax= pi-&#3627408532;??????&#3627408527;
−&#3627409359;
(Pm/Pe3max)

cosdc= (Pm*(dmax-d0)+Pe3max*cos(dmax)-Pe2max*cos(d0))/(Pe3max-Pe2max)


|cosdc| >1
No Critical
Clearing Angle
found
1
delta = delta + 0.01
YES
NO
YES
2
NO

Page | 40































dc = &#3627408516;&#3627408528;&#3627408532;
−&#3627409359;
(&#3627408516;&#3627408528;&#3627408532;&#3627408517;&#3627408516;)
1
dc>dmax
No Critical
Clearing Angle
found
X2 = infinity
d0r=d0*pi/180
dcr=dc*pi/180;
tc = sqrt(2*H*(dcr-d0r)/(pi*f*Pm))

STOP
YES NO
plot (delta, Pe1,Pe2,Pe3)
YES
NO
2

Page | 41

 Algorithm to obtain Delta v/s Time Curve when Damping is Not Considered



























START
INPUT
Pm, E, V
X1, X2, X3, H, f
tc = Fault Clearing time
tf = Final Time


Pe1max = E*V/X1
Pe2max = E*V/X2
Pe3max = E*V/X3
∆t = Dt = 0.01
np = tf /Dt
Pemax=Pe2max
ck=pi*f/H


Take 3 Arrays:
t ( ) to store time
x1( ) to store values of δ
x2 ( ) to store values of
∆??????
∆&#3627408533;
= ω


1
k = 1
YES
d0=&#3627408532;??????&#3627408527;
−&#3627409359;
(Pm/Pe1max)

t (1) = 0
x1 (1)= d0
x2 (1)= 0

k < np
2
3
NO

Page | 42











ωω

x2(k+1)=x2(k)+Dx2b*Dt















t (k) >= tc
plot (delta v/s t )
NO
Pemax=Pe3max


YES
1
t (k+1) = t (k)+ Dt


Solve ω: Dx1b= x2(k)
Solve
∆??????
∆&#3627408533;
: Dx2b=ck*(Pm-Pemax*sin(x1(k)))


From Modified Euler’s Method
x1(k+1)=x1(k)+Dx1b*Dt
x2(k+1)=x2(k)+Dx2b*Dt


Dx1e=x2(k+1)
Dx2e=ck*(Pm-Pemax*sin(x1(k+1)))


Calculate Average Slope
Dx1=(Dx1b+Dx1e)/2
Dx2=(Dx2b+Dx2e)/2


Compute δ: x1(k+1)=x1(k)+Dx1*Dt
Compute
∆??????
∆&#3627408533;
: x2(k+1)=x2(k)+Dx2*Dt


k = k +1


2 3
delta = 180*x1/pi


STOP

Page | 43

 Algorithm to obtain Delta v/s Time Curve when Damping is Considered



























START
INPUT
Pm, E, V
X1, X2, X3, H, f
tc = Fault Clearing time
tf = Final Time


Pe1max = E*V/X1
Pe2max = E*V/X2
Pe3max = E*V/X3
∆t = Dt = 0.01
np = tf /Dt
Pemax=Pe2max
ck=pi*f/H


Take 4 Arrays:
t ( ) to store time
x1( ) to store values of δ
x2 ( ) to store values of
∆??????
∆&#3627408533;
= ω
b ( ) to store Damping Power

1
k = 1
YES
d0=&#3627408532;??????&#3627408527;
−&#3627409359;
(Pm/Pe1max)

t (1) = 0
x1 (1)= d0
x2 (1)= 0

k < np
2
3
NO

Page | 44












ωω

x2(k+1)=x2(k)+Dx2b*Dt














t (k) >= tc
plot (delta v/s t)
NO
Pemax=Pe3max


YES
1
t (k+1) = t (k)+ Dt


Solve ω: Dx1b= x2(k)
Solve b (k)=??????
&#3627409360;
*(0.0976+(0.0804*cos(2*x1(k))))*Dx1b)
Solve
∆??????
∆&#3627408533;
: Dx2b=ck*(Pm - Pemax*sin(x1(k)) - b(k))


From Modified Euler’s Method
x1(k+1)=x1(k)+Dx1b*Dt
x2(k+1)=x2(k)+Dx2b*Dt


Dx1e=x2(k+1)
b(k+1) = ??????
&#3627409360;
*(0.0976+(0.0804*cos(2*x1(k+1))))*Dx1e)
Dx2e=ck*(Pm - Pemax*sin(x1(k+1)) – b(k+1))


Calculate Average Slope
Dx1=(Dx1b+Dx1e)/2
Dx2=(Dx2b+Dx2e)/2


Compute δ: x1(k+1)=x1(k)+Dx1*Dt
Compute
∆??????
∆&#3627408533;
: x2(k+1)=x2(k)+Dx2*Dt


k = k +1


2 3
delta = 180*x1/pi


STOP

Page | 45

MATLAB Programs


 Program to plot Power-Angle Curve and calculate Critical Clearing Angle

functionloadangle(Pm, E, V, X1, X2, X3)
% This program obtains the power angle curves for a one -
machine system before fault, during fault and after the fault
clearance.
% The equal area criterion is applied to find the critical
clearing angle

Pm = input('Generator output power in p.u. Pm = ' );
E = input('Generator e.m.f. in p.u. E = ' );
V = input('Infinite bus-bar voltage in p.u. V = ' );
X1 = input('Reactance before Fault in p.u. X1 = ' );
X2 = input('Reactance during Fault in p.u. X2 = ' );
X3 = input('Reactance aftere Fault in p.u. X3 = ' );
H=input('Inertia Constant H = ' );
f=input('Frequency f = ');
Pe1max = E*V/X1; Pe2max=E*V/X2; Pe3max=E*V/X3;
delta = 0:.01:pi;
Pe1 = Pe1max*sin(delta); Pe2 = Pe2max*sin(delta); Pe3 =
Pe3max*sin(delta);
d0 =asin(Pm/Pe1max); dmax = pi -asin(Pm/Pe3max);
cosdc = (Pm*(dmax-d0)+Pe3max*cos(dmax) -
Pe2max*cos(d0))/(Pe3max -Pe2max);
if abs(cosdc) > 1
fprintf('No critical clearing angle could be found. \n')
fprintf('system can remain stable during this
disturbance.\n\n')
return
else, end
dc=acos(cosdc);
if dc >dmax
fprintf('No critical clearing angle could be found. \n')
fprintf('System can remain stable during this
disturbance.\n\n')
return
else, end
Pmx=[0 pi-d0]*180/pi; Pmy=[Pm Pm];
x0=[d0 d0]*180/pi; y0=[0 Pm]; xc=[dc dc]*180/pi; yc=[0
Pe3max*sin(dc)];
xm=[dmaxdmax]*180/pi; ym=[0 Pe3max*sin(dmax)];
d0=d0*180/pi; dmax=dmax*180/pi; dc=dc*180/pi;
x=(d0:.1:dc);
y=Pe2max*sin(x*pi/180);
y1=Pe2max*sin(d0*pi/180);
y2=Pe2max*sin(dc*pi/180);
x=[d0 x dc];
y=[Pm y Pm];
xx=dc:.1:dmax;

Page | 46

h=Pe3max*sin(xx*pi/180);
xx=[dc xx dmax];
hh=[Pm h Pm];
delta=delta*180/pi;
if X2 == inf

if H ~= 0
d0r=d0*pi/180; dcr=dc*pi/180;
tc = sqrt(2*H*(dcr-d0r)/(pi*f*Pm));
else, end
else,
end

fprintf('\nInitial power angle = %7.3f \n', d0)
fprintf('Maximum angle swing = %7.3f \n', dmax)
fprintf('Critical clearing angle = %7.3f \n\n', dc)
if X2==inf& H~=0
fprintf('Critical clearing time = %7.3f sec. \n\n', tc)
else, end
h = figure; figure(h);
fill(x,y,'m')
hold;
fill(xx,hh,'c')
plot(delta, Pe1,'-', delta, Pe2,'r-', delta, Pe3,'g-', Pmx,
Pmy,'b-', x0,y0, xc,yc, xm,ym), grid
Title('Application of equal area criterion to a critically
cleared system')
xlabel('Power angle, degree' ), ylabel(' Power, per unit')
text(5, 1.07*Pm, 'Pm')
text(50, 1.05*Pe1max,[ 'Critical clearing angle =
',num2str(dc)])
axis([0 180 0 1.1*Pe1max])
holdoff;

Page | 47

 Program to obtainDelta v/s Time Curve when Damping is Not Considered

%Plots Delta v/s Time for Undamped System.
%Modified Euler's Method has been used

function [delta] = undamped(Pm, E, V, X1, X2, X3, H, f, tc,
tf)

%tc = fault clearing time
%tf = final time
%Dt = smallest time interval
Dt = 0.01
Pe1max = E*V/X1; Pe2max=E*V/X2; Pe3max=E*V/X3;
d0 =asin(Pm/Pe1max);
t(1) = 0;
x1(1)= d0;
x2(1)=0;
np=tf /Dt;
Pemax=Pe2max;
ck=pi*f/H;
for k = 1:np
if t(k) >= tc
Pemax=Pe3max;
else, end
t(k+1)=t(k)+Dt;
%Modified Euler's Method begins
Dx1b=x2(k);
Dx2b=ck*(Pm-Pemax*sin(x1(k)));
x1(k+1)=x1(k)+Dx1b*Dt;
x2(k+1)=x2(k)+Dx2b*Dt;
Dx1e=x2(k+1);
Dx2e=ck*(Pm-Pemax*sin(x1(k+1)));
Dx1=(Dx1b+Dx1e)/2;
Dx2=(Dx2b+Dx2e)/2;
x1(k+1)=x1(k)+Dx1*Dt;
x2(k+1)=x2(k)+Dx2*Dt;
end
delta=180*x1/pi;
%delta is an array

Page | 48

 Program to obtainDelta v/s Time Curve when Damping is Considered

%Plots Delta v/s Time for Damped System.
%Modified Euler's Method has been used

function [delta] = damping(Pm, E, V, X1, X2, X3, H, f, tc,
tf)

%tc = fault clearing time
%tf = final time
%Dt = smallest time interval

Dt = 0.01
Pe1max = E*V/X1; Pe2max=E*V/X2; Pe3max=E*V/X3;
d0 =asin(Pm/Pe1max);
t(1) = 0;
x1(1)= d0;
x2(1)=0;
np=tf /Dt;
Pemax=Pe2max;
ck=pi*f/H;
for k = 1:np
if t(k) >= tc
Pemax=Pe3max;
else, end
t(k+1)=t(k)+Dt;

%b() is an array which contains value of damping power
%Modified Euler's Method begins
Dx1b=x2(k);
b(k)=(power(V,2)*(0.0976+(0.0804*cos(2*x1(k))))*Dx1b);
Dx2b=ck*(Pm-Pemax*sin(x1(k))-b(k));
x1(k+1)=x1(k)+Dx1b*Dt;
x2(k+1)=x2(k)+Dx2b*Dt;
Dx1e=x2(k+1);
b(k+1)= (power(V,2)*(0.0976+(0.0804*cos(2*x1(k+1))))*Dx1e)
Dx2e=ck*(Pm-Pemax*sin(x1(k+1))-b(k+1));
Dx1=(Dx1b+Dx1e)/2;
Dx2=(Dx2b+Dx2e)/2;
x1(k+1)=x1(k)+Dx1*Dt;
x2(k+1)=x2(k)+Dx2*Dt;
end
delta=180*x1/pi;

Page | 49






Chapter 8

Results and Discussions

Page | 50
0
20
40
60
80
100
120
140
160
180
200
1
2
3
4
5
6
7
8
9
10
11
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
Delta(degree)
Variation of Damping Factor with Load Andle and Series Reactance
Series Reactance (p.u)
Damping Factor
Variation of Damping Factor

Damping factor

&#3627408439;=
2??????
&#3627408481;
2
πf
360&#3627408467;
[
(??????
&#3627408465;

−??????
&#3627408465;
′′
)&#3627408455;
&#3627408465;0
′′
(??????
&#3627408466;+ ??????
&#3627408465;

)
2
&#3627408480;??????&#3627408475;
2
?????? +
(??????
&#3627408478;

−??????
&#3627408478;
′′
)&#3627408455;
&#3627408478;0
′′
(??????
&#3627408466;+ ??????
&#3627408478;

)
2
&#3627408464;&#3627408476;&#3627408480;
2
?????? ]

For the given system, we assume that the alternator transient and sub transient reactances,
and the time constants remain unchanged. The load angle δ and the series reactance ??????
&#3627408466; is
considered to be variable. Then:

&#3627408439;=
2??????
&#3627408481;
2
πf
360&#3627408467;
[
(0.296−0.2 )∗0.492
(??????
&#3627408466;+0.296)
2
&#3627408480;??????&#3627408475;
2
?????? +
(0.592−0.2 )∗0.07
(??????
&#3627408466;+0.592)
2
&#3627408464;&#3627408476;&#3627408480;
2
?????? ]

Thus &#3627408439;=f(??????,??????
&#3627408466;)










Figure 8.1 shows the variation of damping factor when both load angle and series reactance
are variable. For a particular network, ??????
&#3627408466; is constant and δ is variable. If it is found that the
value damping constant is not sufficient enough to bring the system back to steady state,
then the value of ??????
&#3627408466; can be increased by inserting additional reactance in series with the
alternator.
Figure 8.1: Variation of Damping factor with ??????
&#3627408466; and δ

Page | 51


Critical Clearing Angles

For the given system, Load-Angle curves for before fault, during fault and after clearance
of fault at different locations are shown. The curves have been obtained using the
MATLAB function loadangle( )by entering suitable values of system reactances for
before, during and after fault situations





















Figure 8.2: P v/s δ curve for fault occurring on Bus-Bar
Figure 8.3: P v/s δ curve for fault occurring on a line close to Bus-Bar

Page | 52
0 20 40 60 80 100 120 140 160 180
0
0.5
1
1.5
2
2.5
Application of equal area criterion to a critically cleared system
Power angle, degree
Power, per unit
Pm
Critical clearing angle = 95.8669


A1
A2
P1
P2
P3
P1max = 2.339 p.u
P2max = 0.2521 p.u
P3max = 2.017 p.u














Figure 8.4: P v/s δ curve for fault occurring at 1/4
th
distance from Alternator
Figure 8.5: P v/s δ curve for fault occurring at the center of a transmission line

Page | 53




Observations:

 Figure 8.2 shows the Power-Angle curve when a 3-phase fault occurs on the bus-
bar. During fault, the bus-bar is isolated and the alternator gets disconnected from
the system. As a result, the alternator runs on no-load.
Since no-load operation of alternator for longer time will lead to instability, the bus-
bar is reconnected using &#3627408454;&#3627408441;
6 circuit breaker and microprocessor based relays. Thus
electrical power before fault (??????&#3627408474;
1 ) is equal to electrical power after clearance of
fault (??????&#3627408474;
3 ).

Since electrical power during occurrence of fault is zero, the critical clearing time
using the conventional formula is obtained as 0.252sec.

 Figure 8.3 shows the Power-Angle curve when a 3-phase fault occurs on a
transmission line very close to the bus-bar. Since the fault is very close to the bus-
bar, the bus-bar is isolated during the fault. However after clearance of fault, the
faulty line is removed from the system. Hence ??????&#3627408474;
1 ≠??????&#3627408474;
3

 Figure 8.4 shows the Power-Angle curve when a fault occurs at 1/4
th
distance from
the alternator on a transmission line. Since the reactance during the fault is not zero
like the above two cases, electrical power during fault is also sinusoidal.


Figure 8.6: P v/s δ curve for fault occurring at 3/4th distance from Alternator

Page | 54

 Figure 8.5 shows the Power-Angle curve when a fault occurs at the center of any
one transmission line and Figure 8.6 shows the Power-Angle curve when a fault
occurs at 3/4
th
distance from the alternator on a transmission line. The values of
??????&#3627408474;
1 and??????&#3627408474;
3 for faults in the transmission line is same since after fault, the faulty
line is removed. However values of ??????&#3627408474;
2 are different since line reactance during
fault (??????
2 ) depends on the distance of occurrence of fault.


Fault Location
Maximum
Power During
Fault
(??????&#3627408526;
&#3627409360; ) in p.u
Maximum
Power After
Fault
(??????&#3627408526;
&#3627409361;) in p.u
Initial
Load
Angle ??????
&#3627409358;
(degree)
Maximum
Load Angle
??????
&#3627408526;????????????
(degree)
??????
&#3627408526;????????????− ??????
&#3627409358;
Critical
Clearing
Angle
??????
??????(degree)
On Bus-Bar 0 2.339 20.0034 159.996 139.99 95.9616
On a line very
close to Bus-Bar
0 2.016 20.0034 156.6266 136.6232 88.3925
1/4
th
Distance
From Alternator
0.2521 2.017 20.0034 156.6266 136.6232 95.8669
Middle of
Transmission
Line
0.4047 2.017 20.0034 156.6266 136.6232 101.6282
3/4
th
Distance
From Alternator
0.4081 2.017 20.0034 156.6266 136.6232 101.7287

Table 1: Stability Parameters for Fault at Different Locations

From the above data, it can be observed that when ??????&#3627408474;
2 = 0, ??????
&#3627408464;∝(??????
&#3627408474;????????????−??????
0)


And for the last three cases, when ??????&#3627408474;
1 ≠ 0 and the value of ??????
&#3627408474;????????????−??????
0 is same for all the
cases, then ??????
&#3627408464;∝
1
??????&#3627408474;
3
− ??????&#3627408474;
2

Page | 55

Fault Clearing Time with damping neglected


When damping is neglected, the system behavior is undamped and the variation of δ with
time is sinusoidal. For analysis we plot the δ v/s time curve for only one cycle or a half
cycle.
The variation of δ with time is plotted for different fault clearing times using MATLAB
function undamped( ). The main objective behind this is to determine the critical clearing
time for faults at different locations. When the fault clearing time is more than critical
clearing time, the system becomes unstable and δ reaches infinity.
The conventional formula for calculating critical clearing time is based on the assumption
that electrical power during fault is zero. But this condition is not always true. For non-zero
electrical power transfer during fault, critical clearing time can be found out using complex
integration techniques.
Here we have determined the critical clearing time for 3-phase fault at different locations
using bisection method. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.1sec
tc = 0.2sec
tc = 0.25sec
tc = 0.26sec
Figure 8.7: δ v/s Time curve for fault on Bus-Bar

Page | 56
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.1 sec
tc = 0.2 sec
tc = 0.27 sec
tc = 0.271 sec
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.1sec
tc = 0.2sec
tc = 0.23sec
tc = 0.24sec
Figure 8.8: δ v/s Time curve for fault on a transmission line close to Bus-Bar
Figure 8.9: δ v/s Time curve for fault at 1/4th distance from Alternator

Page | 57
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.25sec
tc = 0.3sec
tc = 0.31sec 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.25sec
tc = 0.3sec
tc = 0.31sec













Figure 8.10: δ v/s Time curve for fault at center of any one Transmission Line
Figure 8.11: δ v/s Time curve for fault at 3/4th distance from Alternator

Page | 58

Observations:

 Figure 8.7 shows the variation of δ with time for a 3-phase fault on the bus-bar, for
different fault clearing times. Since electrical power transmitted during the fault is
zero, the critical clearing time as calculated from the formula is 0.252 seconds.

 Figure 8.8 shows the variation of δ with time for a 3-phase fault on a transmission
line very close the bus-bar, for different fault clearing times. Since electrical power
transmitted during the fault is zero, the critical clearing time as calculated from the
formula is 0.239 seconds.

 Figure 8.9 shows the variation of δ with time for a 3-phase fault at 1/4
th
distance
from alternator, for different fault clearing times. Since the system becomes
unstable for fault clearing times more than 0.27 seconds, it is the critical clearing
time of the system.

 Figure 8.10 and 8.11 shows the variation of δ with time for a 3-phase fault at the
center of transmission line, and at 3/4
th
distance from alternator respectively. The
critical clearing time for both the cases has been found to be 0.3 seconds.



Table 2: Computerized result of Critical Clearing Time for fault at different locations

From the obtained data, we can conclude that critical clearing time increases as δ
C− δ
0
increases. For the given system, inertia constant H = 3 MJ/MVA and frequency f = 50Hz.

Fault Location
Initial Load
Angle ??????
&#3627409358;
(degree)
Critical
Clearing
Angle ??????
??????
(degree)
??????
??????−??????
&#3627409358;
Critical
Clearing
Time&#3627408533;
?????? (sec)
On Bus-Bar 20.0034 95.9616 75.9582 0.25
On a line very close to
Bus-Bar
20.0034 88.3925 68.3891 0.23
1/4
th
Distance From
Alternator
20.0034 95.8669 75.8635 0.27
Middle of
Transmission Line
20.0034 101.6282 81.6248 0.3
3/4
th
Distance From
Alternator
20.0034 101.7287 81.7253 0.3

Page | 59

The above graphs have been obtained by numerical solution of the swing equation using
modified Euler’s method where time interval ∆t = 0.01 sec. The ambiguity in the critical
clearing time for faults at center of transmission line and at 3/4
th
distance from alternator
can be resolved by reducing ∆t, which in turn will reduce the truncation error associated
with modified Euler’s method.

 For fault on Bus-Bar, &#3627408481;
&#3627408464; &#3627408462;&#3627408464;&#3627408481;&#3627408482;&#3627408462;&#3627408473; = 0.252sec
&#3627408481;
&#3627408464; &#3627408464;&#3627408462;&#3627408473;&#3627408464;&#3627408482;&#3627408473;&#3627408462;&#3627408481;&#3627408466;&#3627408465; = 0.25sec
% error = 0.79%

 For fault very close to Bus-Bar, &#3627408481;
&#3627408464; &#3627408462;&#3627408464;&#3627408481;&#3627408482;&#3627408462;&#3627408473; = 0.239sec
&#3627408481;
&#3627408464; &#3627408464;&#3627408462;&#3627408473;&#3627408464;&#3627408482;&#3627408473;&#3627408462;&#3627408481;&#3627408466;&#3627408465; = 0.23sec
% error = 3.76%
When ∆&#3627408481; = 0.001 sec, the data obtained is:
 For fault on Bus-Bar, &#3627408481;
&#3627408464; &#3627408462;&#3627408464;&#3627408481;&#3627408482;&#3627408462;&#3627408473; = 0.252sec
&#3627408481;
&#3627408464; &#3627408464;&#3627408462;&#3627408473;&#3627408464;&#3627408482;&#3627408473;&#3627408462;&#3627408481;&#3627408466;&#3627408465; = 0.215sec
% error = 0.39%

 For fault very close to Bus-Bar, &#3627408481;
&#3627408464; &#3627408462;&#3627408464;&#3627408481;&#3627408482;&#3627408462;&#3627408473; = 0.239sec
&#3627408481;
&#3627408464; &#3627408464;&#3627408462;&#3627408473;&#3627408464;&#3627408482;&#3627408473;&#3627408462;&#3627408481;&#3627408466;&#3627408465; = 0.238sec
% error = 0.418%

The error arising by numerical solution of differential equation using modified Euler’s
method can be improved by reducing the value of ∆t, but it cannot be eliminated. Since the
error is very small, it will not cause a big difficulty is setting the circuit breaker operating
time in a practical system.

Page | 60
0 0.5 1 1.5 2 2.5
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2 sec
tc = 0.3 sec
tc = 0.36 sec
tc = 0.361 sec 0 0.5 1 1.5 2 2.5
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.3 sec
tc = 0.34 sec
tc = 0.341 sec
Fault Clearing Time with damping considered


In the previous case, the effect of internal damping was neglected. Now the swing equation
considering damping is solved using modified Euler’s method. The variation of δ with time
is now plotted using the MATLAB function damping( ).

We now observe the critical clearing time for 3-phase fault at different locations and
compare the values with that obtained when damping was not considered.
Figure 8.13: δ v/s Time curve for fault on Bus-Bar

Figure 8.13: δ v/s Time curve for fault close to Bus-Bar

Page | 61
0 0.5 1 1.5 2 2.5
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.3sec
tc = 0.41sec
tc = 0.42sec 0 0.5 1 1.5 2 2.5
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.3sec
tc = 0.48sec
tc = 0.49sec

Figure 8.14: δ v/s Time curve for fault at 1/4th distance from Alternator
Figure 8.15: δ v/s Time curve for fault at center of a Transmission Line

Page | 62
0 0.5 1 1.5 2 2.5
-50
0
50
100
150
200
250
Variation of Load-angle with Time for different values of Fault Clearing Time (tc)
Time (sec)
Delta (degree)


tc = 0.2sec
tc = 0.3sec
tc = 0.48sec
tc = 0.49sec

Observations:

 Figure 8.12 shows the variation of δ with time for a 3-phase fault on the bus-bar,
for different fault clearing times. The critical clearing time is found to be 0.36
seconds.

 Figure 8.13 shows the variation of δ with time for a 3-phase fault on a transmission
line very close the bus-bar, for different fault clearing times. The critical clearing
time is found to be 0.34 seconds.

 Figure 8.14 shows the variation of δ with time for a 3-phase fault at 1/4
th
distance
from alternator, for different fault clearing times. Since the system becomes
unstable for fault clearing times more than 0.41 seconds, it is the critical clearing
time of the system.


Figure 8.16: δ v/s Time curve for fault at 3/4th distance from Alternator

Page | 63

 Figure 8.15 and 8.16 shows the variation of δ with time for a 3-phase fault at the
center of transmission line, and at 3/4
th
distance from alternator respectively. The
critical clearing time for both the cases has been found to be 0.49 seconds.


Fault Location
Critical
Clearing Time
of Undamped
System
&#3627408533;
?????? (sec)
Critical
Clearing Time
of Damped
System
&#3627408533;
??????

(sec)
&#3627408533;
??????

− &#3627408533;
?????? (sec)
On Bus-Bar 0.25 0.36 0.11
On a line very close to Bus-
Bar
0.23 0.34 0.11
1/4
th
Distance From
Alternator
0.27 0.41 0.14
Middle of Transmission Line 0.3 0.48 0.18
3/4
th
Distance From
Alternator
0.3 0.48 0.18

Table 3: Effect of Damping on Critical Clearing Time

From the above data, it is clear that the value of critical clearing time increases as damping
is included in the system. As the system is stable when fault clearing time is less than the
critical clearing time, damping gives the system more time to clear the fault and regain its
stability. In short, we can say that the stability limit of the system is increased.

Page | 64

Effect of Inertia Constant


To study the effect of change of inertia constant H, we plot the δ v/s time curve for a fault
on the bus-bar with a fault clearing time of 0.2sec. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
-100
-50
0
50
100
150
200
Variation of Load-angle with Time for different values of Inertia Constant (H)
Time (sec)
Delta (degree)


H = 2 MJ/MVA
H = 3 MJ/MVA
H = 4 MJ/MVA
H = 5 MJ/MVA 0 0.5 1 1.5 2 2.5
-10
0
10
20
30
40
50
60
70
80
Variation of Load-angle with Time for different values of Inertia Constant (H)
Time (sec)
Delta (degree)


H = 2 MJ/MVA
H = 3 MJ/MVA
H = 4 MJ/MVA
H= 5 MJ/MVA
Figure 8.17: Effect of change of H on &#3627408481;
&#3627408464;when damping is neglected
Figure 8.18: Effect of change of H on &#3627408481;
&#3627408464; when damping is considered

Page | 65

When damping is neglected, it is observed that increase of H of alternator reduces the angle
through which it swings in a given time interval. The number of cycles in a given time is
also reduces, as well as the peak value of δ is also reduced.
Thus by increasing the inertia constant, the stability of the system can be improved.
For the same fault location when damping is considered, the effect of change of H on δ v/s
time characteristics is obtained.
From the response, it is observed that with increase in value of H, the maximum overshoot
of δ decreases. If we observe for the time period 0 to 1sec, the number of cycles reduces as
value of H increases. However, if we observe for the entire time period till the system
attains stability, the number of cycles remains unaltered. The time at which δ attains
overshoots and undershoots increases with increase in H.
The value of δ for first undershoot and second overshoot however remains unaltered with
increase in value of H. It is also observed that the critical clearing time increases with
increase in H.


Table 4: Effect of change of Inertia Constant on Critical Clearing Time


Fault Location

&#3627408533;
?????? (sec) for
H=2
MJ/MVA
&#3627408533;
?????? (sec) for
H=3
MJ/MVA
&#3627408533;
?????? (sec) for
H=4
MJ/MVA
&#3627408533;
?????? (sec) for
H=5
MJ/MVA
On Bus-Bar 0.31 0.36 0.39 0.43
On a line very close to
Bus-Bar
0.3 0.34 0.38 0.41
1/4
th
Distance From
Alternator
0.37 0.41 0.45 0.49
Middle of Transmission
Line
0.43 0.49 0.52 0.56
3/4
th
Distance From
Alternator
0.43 0.49 0.52 0.56

Page | 66







Chapter 9

Conclusion

Page | 67

In the following project, we have determined the effect of different parameters on the
transient stability of a single machine infinite bus system. But in practice, a power
system consists of a number of alternators. These alternators are interconnected to each
other and can be represented by a single equivalent alternator.

The study subsystem is modeled in details whereas approximate modeling is carried out
for the rest of the subsystem. The qualitative conclusions regarding system stability
drawn an equivalent one-machine infinite bus system can be easily extended to a multi-
machine system.

It can be seen that transient stability is greatly affected by the type and location of a
fault so that a power system analyst must at the very outset of a stability study decide
on these two factors.

For the case of one-machine system connected to infinite bus it can be seen that an
increase in the inertia constant H of the machine reduces the angle through which it
swings in a given time interval offering a method of improving stability. But this
cannot be employed in practice because of economic reasons and for the reason of
slowing down of the response of the speed-governor loop apart from an excessive rotor
weight.

Our aim should be to improvise methods to increase transient stability. A stage has
been reached in technology whereby the methods of improving stability have been
pushed to their limits. With the trend to reduce machine inertias there is a constant need
to determine availability, feasibility and applicability of new methods for maintaining
and improving stability.
We have explored a more detailed model for transient stability analysis taking into
account the effect of damping which is clearly visible from the dynamic response of the
system. We have included a damping factor in the original swing equation which
accounts for the damping taking place at various points within the system. We have
calculated the critical clearing time of the system with and without damping.
We have also discussed how the damping factor can be changed in order to get a
desired result.

Although throughout the project we have worked with a given system, the MATLAB
programs that we have developed can be applied to determine the stability of any single
machine system, where the per unit reactances and the alternator parameters can be
entered by the user.

The response obtained using the programs can be further improved by using any
alternate method of numerical solution of differential equation, which has lesser error.
Here we have used modified Euler’s method due to our convenience and for the ease of
developing the programs.

Page | 68

Aspects of Future Work



To date the computational complexity of transient stability problems have kept them
from being run in real-time to support decision making at the time of a disturbance. If a
transient stability program could run in real time or faster than real time. Then power
system control room operators could be provided with a detailed view of the scope of
cascading failure. This view of unfolding situation could assist an operator in
understanding the magnitude of the problem and its ramifications so that proactive
measure could be taken to limit the extent of the incident. Faster transient stability
simulation implementations may significantly improve power system reliability which
in turn will directly or indirectly affect.

 electrical utility company profits
 environmental impact
 customer satisfaction

In addition to real time analysis, there are other areas where transient stability analysis
could become an integral part of daily power system operations.

 system restoration analysis
 economic / environmental dispatch
 expansion planning

Page | 69

References


[1] Edward Wilson Kimbark, “Power System Stability, Volume I: Elements of
Stability Calculations”. New York: I.E.E.E Press, 1995. Print.

[2] Edward Wilson Kimbark, “Power System Stability, Volume III: Synchronous
Machines”. New York: I.E.E.E Press, 1995. Print.

[3] J.B. Gupta, “A Course in Power Systems”. Ludhiana, S.K. Kataria & Sons,
1968. Print.

[4] D.P. Kothari, I.J Nagrath, “Power System Engineering”. New Delhi, McGraw
Hill Education (India) Private Limited, 2007. Print.

[5] A. Chakrabarti, Sunita Halder, “Power System Analysis, Operation and
Control”. New Delhi, Prentice-Hall of India Private Limited, 2006. Print.

[6] S.A. Mollah, S. Chakrabarty, “Computing Systems”. Kolkata, J.B Books and
Learning, 2012. Print.

[7] Rene Rossi, Mohammad A.S. Masoum, “A Practical Approach Based on Equal
Area Criteria for Stability Analysis of Weak Networks with Wind Generation
Penetration”, Journal of Energy and Power Engineering, vol. 6, pp. 629-637,
April 2011.

[8] S. Muknahallipatna, Stanislaw Legowski, Sadrul Ula, Jason Kopas, “Power
System Transient Stability Analysis Software Tool for an Undergraduate
Curriculum”, John Wiley & Sons, Inc. Computer Applications in Engineering
Education, vol. 9, pp. 37-48, March 2001.

[9] D. P. Koester, S. Ranka, and G. C. Fox, “Power Systems Transient Stability -
A Grand Computing Challenge”. The Northeast Parallel Architectures Center
(NPAC) Technical Report - SCCS 549, August 1992.

[10] Richard Johnson, “MATLAB
®
Programming Style Guidelines”, October 2002.

[11] The MathWorks, Inc., “MATLAB
®
Creating Graphical User Interfaces”,
November 2000.