ABSTRACT
Title of dissertation: PROBLEMS IN SELECTIVE MODAL
ANALYSIS AND CONTROL
Mohamed S. Saad, Doctor of Philosophy, 2007
Dissertation directed by: Professor Eyad H. Abed
Department of Electrical and Computer Engineering
In this dissertation, we develop monitoring and control systems for improv-
ing the performance of systems that are required to operate at the edge of their
stability envelopes. The concept of modal participation factors, which is an essen-
tial construct in the theory of Selective Modal Analysis, is used extensively in this
work. The basic deflnition of modal participation factors that was originally given
for unforced linear time-invariant systems is revisited, and related notions of out-
put participation factors and input-to-output participation factors are introduced,
studied and applied to models of electrical power systems.
A signal-based approach for real-time detection of impending instability in
nonlinear systems is considered. The main idea pursued involves using a small ad-
ditive white Gaussian noise as a probe signal and monitoring the spectral density
of one or more measured states or outputs for certain signatures of impending in-
stability. Input-to-state and input-to-output participation factors are used as tools
to aid in selection of locations for probe inputs and states or outputs to be mon-
itored, respectively. Since these participation factors are model-based, the work
presented combines signal-based and model-based ideas toward achieving a robust
methodology for instability monitoring. Case studies from power systems are used
to illustrate the developed monitoring system, one of which involves the WSCC
3-generator, 9-bus network.
Feedback algorithms are developed for assigning modal participation factors
in general linear time-invariant systems using eigenvector assignment-based tech-
niques. The goal is to reduce the interaction between a selected group of states (the
high-value group) and an undesirable mode (for example a critical mode, i.e., one
corresponding to an eigenvalue or pair of eigenvalues approaching the imaginary axis
in the complex plane). In particular, we address two cases, one in which the mode
of interest corresponds to a real eigenvalue approaching zero, and the other in which
the mode involves a complex conjugate pair of eigenvalues that may be approaching
the imaginary axis. A novel procedure for computing the desired closed-loop right
eigenvector(s) associated with the critical mode (based on given constraints on the
desired closed-loop participation factors) is presented. An example from power sys-
tems is presented to demonstrate the efiectiveness of the controller. The example
used is the WSCC 3-generator, 9-bus network.
PROBLEMS IN SELECTIVE MODAL ANALYSIS AND
CONTROL
by
Mohamed S. Saad
Dissertation submitted to the Faculty of the Graduate School of the
University of Maryland, College Park in partial fulflllment
of the requirements for the degree of
Doctor of Philosophy
2007
Advisory Committee:
Professor Eyad H. Abed, Chair/Advisor
Professor William S. Levine
Professor Gilmer L. Blankenship
Professor Amr Baz
Assistant Professor Richard J. La
c Copyright by
MOHAMED S. SAAD
2007
DEDICATION
To my parents, brothers
To my wife Maha and my kids Karim and Malak
ii
ACKNOWLEDGMENTS
First and foremost I?d like to thank my advisor, Professor Eyad H. Abed for
his endless support and advise throughout my graduate studies at the University of
Maryland. He has been a great source of inspiration throughout the years.
I would also like to thank Professors William Levine, Gilmer Blankenship,
Amr Baz and Richard La for agreeing to serve on my thesis committee and for
sparing their invaluable time reviewing the manuscript.
I would like to thank Professor Andre Tits for the helpful discussions on
my research. I would also like to thank Professor Robert Israel of the University
of British Columbia for the help with my research questions. Thanks are due to
Professor David Elliott for his advice on research. Thanks are also due to Dr.
Monther Hasouneh for the great time we spent on discussing research and other
general life topics.
I owe my deepest thanks to my family - my mother and father who have
always stood by me, gave me endless patience and love and guided me through my
career. Words cannot express the gratitude I owe them. Deep thanks are also due
to my dear brothers Ahmed and Maged.
Special and deep thanks to my wife Maha. She has always been surrounding
iii
me with care and love. I will never forget her nice encouraging words throughout
the time of my PhD. Also, I would like to thank my lovely twins, Karim and Malak
for the joy they have been fllling my life with.
I am grateful to all my colleagues at the University of Maryland; Chao Wu,
Amirali Sharifl, Maben Rabi, Nassir BenAmmar and Guojian Lin. I will never forget
our nice discussions on every aspect of life.
I am also grateful to my close friends during this whole PhD journey; Mo-
hamedAbdallah, MohamedZahran, Abdel-hameedBadawy, AhmedElsaid, Mostafa
Arafa, Tarek Massoud and Mohamed Fahmi. I hope that our friendship will remain
forever.
It is impossible to remember all, and I apologize to those I?ve inadvertently
left out.
iv
Contents
List of Tables ix
List of Figures xi
1 Introduction 1
1.1 Thesis Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2 Mathematical Background 12
2.1 Modal Participation Factors . . . . . . . . . . . . . . . . . . . . . . . 12
2.1.1 Conventional participation factors . . . . . . . . . . . . . . . . 13
2.1.2 Generalized participation factors . . . . . . . . . . . . . . . . 16
2.2 Bifurcation Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.2.1 Stationary bifurcation . . . . . . . . . . . . . . . . . . . . . . 19
2.2.2 Hopf bifurcation . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.3 Eigenstructure Assignment for
Control System Design . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.3.1 Eigenstructure assignment using state feedback . . . . . . . . 26
v
3 Output Participation Factors 32
3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.2 Output Participation Factors . . . . . . . . . . . . . . . . . . . . . . 35
3.3 Input-to-Output Modal Participation Factors . . . . . . . . . . . . . 45
4 Detecting Impending Instability in Nonlinear Systems 47
4.1 Precursor-Based Monitoring . . . . . . . . . . . . . . . . . . . . . . . 48
4.1.1 Developed monitoring system . . . . . . . . . . . . . . . . . . 49
4.1.2 Instability proximity index . . . . . . . . . . . . . . . . . . . . 53
4.2 Power System Examples . . . . . . . . . . . . . . . . . . . . . . . . . 55
4.2.1 Single generator with dynamic load . . . . . . . . . . . . . . . 55
4.2.2 Single generator connected to an inflnite bus . . . . . . . . . . 57
4.2.3 Three-generator nine-bus power system . . . . . . . . . . . . . 63
5 Participation Factor Assignment Using Eigenvector Assignment-
Based Feedback Control 74
5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
5.2 Participation Factor Assignment: Case of
Real Eigenvalue Approaching Zero . . . . . . . . . . . . . . . . . . . 77
5.2.1 Computing desired closed-loop right eigenvector ri0d . . . . . . 78
5.2.2 Feedback control design: exactly assignable case . . . . . . . . 80
5.2.3 Feedback control design: not exactly assignable case . . . . . . 84
vi
5.2.4 Feedback control design: \sacriflcial" case . . . . . . . . . . . 85
5.3 Participation Factor Assignment: Case of
Eigenvalue Pair Approaching Imaginary Axis . . . . . . . . . . . . . . 90
5.3.1 Computing desired closed-loop right eigenvectors
fri0d; ri+10d = ri0d?g . . . . . . . . . . . . . . . . . . . . . . . . 91
5.3.2 Feedback control design: exactly assignable case . . . . . . . . 95
5.3.3 Feedback control design: not exactly assignable case . . . . . . 98
6 Application to Electric Power Networks 101
6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102
6.2 Three-Generator Nine-Bus Power System . . . . . . . . . . . . . . . . 102
6.2.1 Feedback control design . . . . . . . . . . . . . . . . . . . . . 107
6.3 Remarks on Use of State Feedback in Power Networks . . . . . . . . . 112
7 Conclusions and Future Directions 114
7.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114
7.2 Future Directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116
A Proof of Propositions 5.1 and 5.2 118
A.1 Proof of Proposition 5.1 . . . . . . . . . . . . . . . . . . . . . . . . . 118
A.2 Proof of Proposition 5.2 . . . . . . . . . . . . . . . . . . . . . . . . . 120
B Table of Symbols 122
vii
Bibliography 123
viii
List of Tables
4.1 Input-to-state participation factors and spectral peaks at !c. . . . . . 60
4.2 Parameter values for the single generator connected to an inflnite bus
model. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.3 Input-to-state participation factors and spectral peaks at !c ? 7:8
for the single generator connected to an inflnite bus system. . . . . . 63
4.4 Three-generator nine-bus power system: machine data. . . . . . . . . 68
4.5 Three-generator nine-bus power system: exciter data. . . . . . . . . . 68
4.6 Input-to-state participation factors for the 3-machine nine-bus system
(partial listing). The load at bus 5 is 4.4 pu. . . . . . . . . . . . . . . 68
4.7 Input-to-output participation factors for the 3-machine nine-bus sys-
tem (partial listing). The load at bus 5 is 4.4 pu. . . . . . . . . . . . 69
5.1 Closed-loop ratios fz1; z2g for coe?cient vectors fk1; k2; kcadg. The
desired closed-loop ratios are fzd1 = 1:2; zd2 = 2:24g. . . . . . . . . . . 100
ix
6.1 The states with the highest norms of open-loop participation factors
in the critical mode associated with ?fi;i+1g = ?0:04602 ? j2:1151.
The load at bus 5 is 4.4 pu. . . . . . . . . . . . . . . . . . . . . . . . 105
6.2 Closed-loop ratios fz13; z23g for coe?cient vectors fk1; k2; kcadg.
The desired closed-loop ratios are fz13 = 2:92; z23 = 2:03g. . . . . . . 108
6.3 The states with the highest norms of closed-loop participation factors
in the critical mode associated with ?fi;i+1g = ?0:04602 ? j2:1151.
The load at bus 5 is 4.4 pu. . . . . . . . . . . . . . . . . . . . . . . . 109
x
List of Figures
2.1 Transcritical bifurcation. . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.2 Pitchfork bifurcation. (a) supercritical, (b) subcritical . . . . . . . . . 22
2.3 Andronov-Hopf bifurcation. (a) supercritical, (b) subcritical . . . . . 24
3.1 Illustrating diagonal transformation T in two dimensions, (a) Elliptic
original uncertainty set ?, (b) Rectangular original uncertainty set ?. 43
4.1 Precursor-based instability monitor with external probe signal. . . . . 49
4.2 Power spectral densities of the states of the model given in (4.5)-(4.8).
The bifurcation parameter was set to Q1 = 2:9. Band-limited AWGN
of zero mean and (0:001)2 power was added to Pm. . . . . . . . . . . 58
4.3 Variation of the peak value of the power spectral density of ! as a
function of the bifurcation parameter Q1. Band-limited AWGN of
zero mean and (0:001)2 power was added to Pm. . . . . . . . . . . . . 59
xi
4.4 Power spectral densities of the states of the single generator connected
to an inflnite bus system. The bifurcation parameter was set to KA =
185. Band-limited AWGN of zero mean and (0:000032)2 power was
added to Pm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
4.5 Power spectral densities of the states of the single generator connected
to an inflnite bus system. The bifurcation parameter was set to KA =
185. Band-limited AWGN of zero mean and (0:000032)2 power was
added to Vref. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.6 Variation of the peak value of the power spectral density of VR as a
function of the bifurcation parameter KA. Band-limited AWGN of
zero mean and (0:000032)2 power was added to Pm. . . . . . . . . . . 66
4.7 Plot of the inverse of the power spectral density index of ! versus the
bifurcation parameter KA. The threshold value is equal to 5. . . . . . 67
4.8 WSCC 3-machine, 9-bus system. . . . . . . . . . . . . . . . . . . . . . 69
4.9 Power spectral densities of the states Efd1, Efd2 and Efd3. The load
on bus 5 was used as a bifurcation parameter. The load value is 4.0
pu. Band-limited AWGN of zero mean and 0:05 power was added to
Pm3, the mechanical power of generator number 3. . . . . . . . . . . . 70
xii
4.10 Power spectral densities of the states Efd1, Efd2 and Efd3. The load
on bus 5 was used as a bifurcation parameter. The load value is 4.4
pu. Band-limited AWGN of zero mean and 0:05 power was added to
Pm3, the mechanical power of generator number 3. . . . . . . . . . . . 71
4.11 Power spectral density of Efd1 for three values of PL5 (the load on bus
5): PL5 = 4:0 pu (dash-dotted line), PL5 = 4:25 pu (dashed line) and
PL5 = 4:4 pu (solid line). Band-limited AWGN of zero mean and 0:05
power was added to Pm3, the mechanical power of generator number 3. 72
4.12 Power spectral densities of the state Efd1 and outputs Q14 and Q27.
The load on bus 5 was used as a bifurcation parameter. The load
value is 4.4 pu. Band-limited AWGN of zero mean and 0:05 power
was added to Pm3, the mechanical power of generator number 3. . . . 73
6.1 Exciter for the three-generator, nine-bus system. . . . . . . . . . . . . 104
6.2 Open-loop response for a perturbation of 2% in the initial condition of
the q-axis transient voltage of generator 2 (e00q2) of the WSCC 3-gen.,
9-bus power system. (a) e0q3, (b) e0q2 and (c) e0q1. . . . . . . . . . . . . 107
6.3 Plot of q-axis transient voltage of generator 1 (e0q1) versus time for
a perturbation of 2% in the initial condition of the q-axis transient
voltage of generator 2 (e00q2) before and after control. . . . . . . . . . . 110
xiii
6.4 Plot of q-axis transient voltage of generator 2 (e0q2) versus time for
a perturbation of 2% in the initial condition of the q-axis transient
voltage of generator 2 (e00q2) before and after control. . . . . . . . . . . 111
6.5 Plot of q-axis transient voltage of generator 3 (e0q3) versus time for
a perturbation of 2% in the initial condition of the q-axis transient
voltage of generator 2 (e00q2) before and after control. . . . . . . . . . . 112
xiv
Chapter 1
Introduction
This dissertation focuses on development of monitoring and feedback control sys-
tems for managing the relative impact of certain critical modes on the states of a
dynamical system. The work begins with a review of modal participation factors,
which are quantities measuring the relative impact of a system mode in system
states. Participation factors provide a convenient tool for determining which states
should be monitored in order to detect that a critical mode is approaching the imag-
inary axis, which would lead to instability of the system. They can also be used
to determine which inputs should be used in state feedback control to mitigate the
efiect of a critical mode on system states that may correspond to equipment that
must be protected. We flnd that the conventional notion of participation factors is
insu?cient for monitoring and control when we have access not to system states,
but to system output variables. This motivates us to introduce several new con-
cepts extending the traditional notion of participation factors. These include output
1
participation factors and input-to-output participation factors. We use these new
concepts to develop a signal-based approach for real-time detection of impending
instability in nonlinear systems. The main idea pursued involves using a small ad-
ditive white Gaussian noise as a probe signal and monitoring the spectral density of
one or more measured states or outputs for certain signatures of impending insta-
bility. Input-to-state and input-to-output participation factors are used as tools to
aid in selection of locations for probe inputs and states or outputs to be monitored.
We also apply eigenstructure assignment techniques to develop procedures for as-
signment of participation factors. We provide examples from the fleld of electric
power systems demonstrating the concepts introduced in this thesis, both for modal
monitoring and control.
The conventional deflnition of participation factors was introduced in the
framework of Selective Modal Analysis (SMA) in [2, 34]. SMA is a methodol-
ogy for the modal analysis of linear time-invariant systems that emphasizes scale-
independent notions that facilitate modal analysis, control design, and order reduc-
tion. The number of modes involved in a phenomenon of interest is usually only
a small fraction of those appearing in the dynamics of a system, usually conflned
to a single mode. Participation factors are non-dimensional scalars that measure
the interaction between modes and the state variables of a linear system. Since
their introduction, participation factors have been used for stability analysis, order
reduction, actuator placement and controller design in a variety of flelds, especially
2
electric power systems [9, 10, 11]. Participation factors have recently also been used
in coherency and clustering studies in electrical power systems [31, 32].
In [15], a new approach to deflning modal participation factors for unforced
systems was presented. This approach involves taking an average or a probabilistic
expectation of a quantitative measure of relative modal participation. The average
is performed with respect to all possible values of the initial state vector, assumed
to be unknown but to lie in a known set or satisfy a known probability density.
The deflnitions found in [15] were shown to reduce to the original deflnition of
participation factors of [2, 34] if the uncertainty in the initial condition satisfles a
symmetry condition.
There is a need to extend the participation factors deflnitions referred to
above in several directions. Two extensions that are pursued in this thesis are output
participation factors and input-to-output participation factors. The need for these
extensions becomes apparent when we consider systems for which the states may not
be readily available for measurement or as feedback signals, while other variables
(namely certain system output variables) are readily available. In such cases, it
is useful to have a criterion for determining which output signal carries the most
information about a critical mode of interest. In some applications, for example,
clustering in power system studies, participation of modes in variables that are not
state variables (e.g., bus voltages) are needed.
The concept of participation of modes in output variables has not previously
3
been considered in the literature. Such a concept can be useful in many applica-
tions such as clustering studies and in stability monitoring and detection of closeness
to instability. In clustering in power systems, for example, development of output
participation factors will facilitate studies and understanding of intracluster oscil-
lations (i.e., oscillations resulting from modes within the cluster) and intercluster
oscillations (i.e., oscillations resulting from modes common to more than one clus-
ter [31, 32]. It will be also useful for local and coordinated control design. For
example, in control design to dampen intercluster oscillations, coordinated control
(distributed control) must be used. While local controls guarantee stable operation
within a cluster, poor choices of local controls might lead to instabilities caused by
the coupled operations of the clusters.
Recently, Canizares and co-workers have proposed some extensions to the
participation factors notions available in the literature, also targeting the incor-
poration of output variables in the deflnitions. In [35], a concept of \extended
eigenvector" was introduced and used to identify the dominant output variable as-
sociated with a critical mode. In [36] an observability index was used to determine
the preferred output signal for observing a critical mode. They have used their new
deflnitions in a modal-based analysis technique for monitoring and control of electric
power systems [35, 36].
The concept of \mode observability factor," introduced by [26], has been
used to determine the best output for either monitoring a mode, or possibly as a
4
good candidate feedback signal to control a mode of interest. In particular, the
mode observability factor was used for determining the most suitable buses in a
power system for placing static VAR compensators in order to dampen the critical
modes of oscillation.
A drawback to using the deflnitions referred to above is that it is scale-
dependent, i.e., for a change in the units of system outputs, the factors values will
change. Another drawback to using these deflnitions is that they do not reduce
to the conventional state participation factors when the output is simply a state
variable.
The thesis then focuses on using the newly introduced concept of input-to-
output participation factors and the concept of input-to-state participation factors
given in [40] to develop a signal-based approach for real-time detection of impend-
ing instability in nonlinear systems with examples from power systems. Kim and
Abed [25] developed monitoring systems for detecting impending instability in non-
linear systems that builds on Wiesenfeld?s research on \noisy precursors of bifurca-
tions". Wiesenfeld?s research was originally introduced to characterize and employ
the noise ampliflcation properties of nonlinear systems near bifurcations of various
types [37, 38]. Noisy precursors are features of the power spectral density (PSD)
of a measured output of a system nearing instability that is excited by additive
white Gaussian noise (AWGN). In [25], the noisy precursors concept was extended
from systems operating at limit cycles to systems operating near equilibria, and
5
closed-loop monitoring systems were developed to facilitate use of noisy precur-
sors in revealing impending loss of stability for such systems. It was shown in [25]
that systems driven by white Gaussian noise and operating near an equilibrium
point exhibit sharply growing peaks near certain frequencies as the system nears
a bifurcation. In particular, it was shown that for stationary bifurcation where an
eigenvalue passes through the origin (as in the case of pitchfork or transcritical bi-
furcation), the peak in the PSD occurs at zero frequency. Analogously, for the case
of Hopf bifurcation (complex conjugate pair of eigenvalues crossing the imaginary
axis transversely), the peak in the PSD occurs near !c, the critical frequency of the
Hopf bifurcation.
As noted by Hauer [19], the recurring problems of system oscillations and
voltage collapse in power systems are due in part to system behavior not well cap-
tured by the models used in planning and operation studies. In the face of compo-
nent failures, system models quickly become mismatched to the physical network,
and are only accurate if they are updated using a powerful and accurate failure
detection system. Therefore, it is important to employ nonparametric techniques
for instability monitoring.
The use of probe signals is shown to help reveal an impending loss of stabil-
ity. This is because probe signals propagate in the power system and give certain
signatures near an instability that can be used as a warning signals for possible
impending voltage collapse. Such warning signals are needed to alert system opera-
6
tors of a situation that may require preventive control, and to provide the operators
with valuable additional time to take necessary preventive (rather than corrective)
measures.
The focus of the thesis then shifts to addressing possible actions to be taken
to counter the possible efiects of the impending instability on the system. The need
to carry such actions can be attributed to the fact that in many systems composed
of smaller subsystems, it may be desirable to assign the magnitudes of participation
factors corresponding to various subsystems to desired values. This would ensure
that the energy corresponding to a critical mode minimally afiects subsystems of
high value, directing the energy of this mode to other, possibly purposefully intro-
duced, subsystems of lesser value to the system operator. This idea is based on the
fact that the higher the magnitude of participation factor, the higher the tendency
of a subsystem to be afiected in a serious way by an impending (or true) instability.
Power systems are a good example for this discussion. In power systems operation,
nuclear units, for example, are regarded as high-value subsystems. Reducing partic-
ipation factors of these units in some critical mechanical mode would reduce the risk
of premature tripping, as well as alleviate the efiect on these units of a troublesome
mode being excited.
The goal of the control efiort will not be to alter the position of the eigen-
value(s) (associated with the critical mode), but rather to reduce the participation
factors of a high-value group of states in that mode. This means that the closed-loop
7
eigenvalues will be kept the same as the open-loop ones. Similarly, all open-loop
right eigenvectors will be kept unchanged except for those eigenvector(s) associated
with the critical mode. This is possible using the freedom ofiered by state feedback
in multivariable systems beyond closed-loop eigenvalue assignment [17]. Moreover,
and as noted in the seminal paper of Moore [17], \One advantage of eigenstructure
assignment is that you can still assign eigenvectors associated with uncontrollable
eigenvalues." This is important in participation factor assignment where it is desired
that open-loop eigenvalues are kept unchanged.
A state feedback control system is given for assigning modal participation
factors using eigenvector assignment-based techniques. The idea of using eigenvector
assignment-based techniques to assign participation factors has not previously been
addressed. Some work has addressed the problem of assigning right eigenvector(s)
associated with a mode of interest (usually the critical mode). In [4], the author
alerted the power systems community to the possible use of eigenvector assignment
to make a target area participate highly in an unstable mode. This could be done by
an unscrupulous provider of power to hurt the operation of a competing provider.
In [13], a control scheme was proposed to assign eigenvalues and eigenvectors to
improve system damping and dynamic behavior, respectively. Since right and left
eigenstructure problems are strongly interdependent, assigning participation factors
is not an easy task. Indeed, if the required objective is not exactly attainable,
participation factor assignment may require solving an optimization problem for
8
achieving an objective as close as possible to the desired one.
Inspired by vibration absorbers in mechanical and civil engineering, in this
work we desire to design controllers that reduce participation of certain, \high-
value" system states in a particular mode at the expense of other \low-value" states.
Those low-value states can be viewed as taking a role akin to that of vibration
absorbers in mechanical systems. In power systems for instance, an isolated remote
generator station or  ywheel could be regarded as a vibration absorber that protects
heavilyloadedhigh-capacitygeneratingstations. Wewillrefertotheaforementioned
problem as the participation factor assignment problem.
Reducing the participation factors of a high-value group of states in the
critical mode implies that some other groups of states will be burdened with higher
magnitudes of their participation factors in the critical mode. Deciding which group
of states will be seen as \low-value" is based on system conditions and operator
experience.
Control systems that assign participation factors should work in unison with
monitoring systems that detect closeness to instability. The control becomes active
when the monitoring system provides an alert that instability is near. The objective
of participation factor assignment in this setting is as described above: to lessen the
impact of a critical mode on high-value states by a trade with impact on low-value
states.
9
1.1 Thesis Outline
The dissertation proceeds as follows. In Chapter 2, theoretical background mate-
rial used throughout the thesis is collected. The topics include bifurcation theory,
modal participation factors and eigenstructure assignment-based feedback control
algorithms.
In Chapter 3, a study is given on extending the uncertain initial condition
approach to deflning state participation factors in [15] to facilitate flnding a general
measure of output participation. The new proposed deflnition has the desirable
property that it reduces to the original state participation deflnition when the out-
puts are simply the system states. Also, this new deflnition is unit independent,
i.e., the output participation factors do not depend on the units of system states or
outputs. Another notion of modal participation is presented in this chapter. This
new deflnition is one of input-output participation factors. It involves probe input
signals as well as modes and outputs, and is applicable to a large set of systems.
In Chapter 4, a signal-based approach for real-time detection of impending
instability in nonlinear systems is considered. The developed monitoring system is
based on using a small additive white Gaussian noise as a probe signal and monitor-
ing the spectral density of one or more measured states or outputs for certain sig-
natures of impending instability. Input-to-output participation factors (introduced
in Chapter 3), and input-to-state participation factors [40] (discussed in Chapter
2) are used as tools to aid in selection of locations for probe inputs and outputs or
10
states to be monitored. Examples from power systems are given to demonstrate the
performance of the monitoring system.
In Chapter 5, state feedback control is presented for assigning modal partic-
ipation factors for systems nearing instability using eigenvector assignment-based
techniques. The focus is on assigning participation factors for state variables within
a predeflned group of state variables, where the group can be any subvector of the
state vector, or could consist of all system states. A procedure for computing the
desired closed-loop right eigenvector(s) (based on a given desired closed-loop par-
ticipation factors) is given. Two cases are addressed, flrst the case of systems with
a mode associated with a single real eigenvalue approaching zero, and second the
case of systems with a mode associated with a pair of conjugate complex eigenvalues
approaching the imaginary axis.
In Chapter 6, the results obtained in Chapter 5 are applied to a sample power
system operating close to its stability limits. The considered power system is the
WSCC 3-generator, 9-bus network. As expected, participation factor assignment
has a clear efiect on the relative time domain responses of system states.
Finally, we collect concluding remarks in Chapter 7, along with possible direc-
tions for future research. Some of the results reported in this thesis were published
in a conference paper and a book chapter [22, 21].
11
Chapter 2
Mathematical Background
In this chapter, the theoretical and algorithmic background materials employed
throughout this dissertation are reviewed. Section 2.1 focuses on modal partici-
pation factors, a basic construct in the theory of Selective Modal Analysis [2, 34].
In Section 2.2 a summary on bifurcation theory is provided, and various types of
continuous-time local bifurcations are reviewed. In section 2.3, background material
on eigenstructure assignment for control system design is reviewed.
2.1 Modal Participation Factors
In the theory of Selective Modal Analysis [2, 34], participation factors are non-
dimensional scalars that measure the interaction between the modes and the state
variables of a linear time-invariant system. The notion of participation factors pro-
vides a tool for modal participation analysis where an e?cient method is required to
single out a selected portion of the model related to the dynamics of interest. The
12
conventional deflnition of participation factors has been used as a tool for stability
analysis, modal reduction, actuator placement and controller design in a variety of
flelds, especially electric power systems [9, 11, 10]. The number of modes involved
in the phenomenon of interest is usually only a small fraction of those governing the
dynamics of the system, mostly restricted to only one mode. Participation factors
have also recently been used in coherency and clustering studies in electrical power
systems [31, 32].
Another potential use of modal participation factors is in the area of stability
monitoring of nonlinear systems close to instability. Such systems are also known
to be close to a bifurcation in their dynamics. Input-to-state participation factors
(given in [40]) and input-to-output participation factors (introduced in Chapter 3)
will be employed in Chapter 4 to design stability monitoring systems for detecting
impending instability in nonlinear systems.
2.1.1 Conventional participation factors
Consider a general continuous linear time-invariant system
_x(t) = Ax(t) (2.1)
where x(t) 2 Rn, and A is a real (n ? n) matrix. Suppose that A has a set of
n distinct eigenvalues f?1;?2;:::;?ng. Let fr1;r2;:::;rng be right eigenvectors
of the matrix A associated with the eigenvalues f?1;?2;:::;?ng, respectively. Let
fl1;l2;:::;lng denote left (row) eigenvectors of the matrix A associated with the
13
eigenvalues f?1;?2;:::;?ng, respectively. The right and left eigenvectors are taken
to satisfy the normalization
lirj = ?ij (2.2)
where ?ij is the Kronecker delta:
?ij =
8
>><
>>:
1 i = j
0 i 6= j
The participation factor of the i-th mode in the k-th state is deflned to be the
complex number [2, 34, 15]
pki := likrik (2.3)
This formula also gives the participation of the k-th state in the i-th mode. Partic-
ipation factors is a measure of the relative participation of modes in states and the
relative participation of states in modes. Since rik measures the activity of xk in the
i-th mode and lik weighs the contribution of this activity to the mode, the product
pki measures the net participation.
Properties of participation factors
The following properties of participation factors are recalled from [1, 3] and can be
easily proved through appropriate consideration of the Deflnition in (2.3) and the
normalization assumption in (2.2).
1. Participation factors are dimensionless quantities that are independent of the
units in which state variables are measured.
14
2. In view of the eigenvector normalization (2.2), the sum of the participation fac-
tors associated with any mode (Pni=1 pki) or with any state variable (Pnk=1 pki)
is equal to 1.
3. The participation factor pki represents the sensitivity of the eigenvalue ?i to
the the diagonal element akk of the state matrix A
pki = @?i@a
kk
:
The solution to system dynamic equations in (2.1) for an initial condition x0
is given by:
x(t) = eAtx0 (2.4)
Since the state matrix A is assumed to have distinct eigenvalues, then A could be
written in terms of right and left modal matrices as well as a diagonal matrix of its
eigenvalues. Then x(t) in (2.4) is
x(t) =
nX
i=1
(lix0)e?itri (2.5)
Now suppose the initial condition x0 is ek, the unit vector along the k-th coordinate
axis. Then the evolution of the k-th state becomes
xk(t) =
nX
i=1
(likrik)e?it (2.6)
=
nX
i=1
pkie?it: (2.7)
15
Equation (2.7) indicates that pki can be viewed as the relative participation of the
i-th mode in the k-th state at time t = 0 1.
2.1.2 Generalized participation factors
The concept of participation factors of modes in states and vice versa has been
extended to linear time-invariant systems with inputs [40]
_x(t) = Ax(t)+Bu(t) (2.8)
y(t) = Cx(t): (2.9)
In [40], the case where the input is applied to one component, say the q-th compo-
nent, of the right side of (2.8) and only one state, say the k-th state, is measured is
considered. That is, in equations (2.8)-(2.9), B and C take the form
B = eq = [0 ::: 0 1|{z}
q?th
0:::0]T;
C = ekT = [0 ::: 0 1|{z}
k?th
0:::0]:
1In [15], a new approach to deflning modal participation factors was presented. The new
approach involved taking an average or a probabilistic expectation of a quantitative measure of
relative modal participation over an uncertain initial state vector. The new deflnitions were shown
to reduce to the original deflnition of participation factors of [2, 34] if the initial state obeys a
symmetry condition.
16
With this choice of C and B, in steady state the output in (2.9) (in the frequency
domain) is given by
y(s) = xk(s) = C(sI ?A)?1Bu
=
nX
i=1
CriliB
s??i u(s)
=
nX
i=1
rikliq
s??iu(s)
The generalized participation factor, participation factor of mode i in state k when
excitation is applied to state q, is deflned as follows:
piqk := r
i
kl
i
q
kHkk (2.10)
Here kHkk := jri~kri~qj is the normalization factor where
(~k; ~q) = arg maxp;q (jrikliqj)
attains maximum participation of norm 1. The above normalization makes the
participation factors to be bounded by zero and one in amplitude and independent
of system scaling. In [1], the quantity piqk = rikliq is called a generalized participation.
We will call this quantity the input-to-state participation factor (ISPF) for mode
i, with measurement at state k and input applied to state q. The input-to-state
participation factor (ISPF) will be employed in designing monitoring systems of
impending instabilities in nonlinear systems in Chapter 4.
17
2.2 Bifurcation Theory
The term bifurcation is commonly used to describe a qualitative change in steady
state behavior as the parameters on which the dynamical system depends are var-
ied through certain values [6, 7, 5]. The parameter being varied is referred to as
the bifurcation parameter. In other words, a bifurcation is a change in the num-
ber of candidate operating conditions of a nonlinear dynamical system that occurs
as a parameter is quasi-statically varied. A value of the bifurcation parameter at
which a bifurcation occurs is called a critical value of the bifurcation parameter. A
nonlinear system operating at an equilibrium point undergoes a bifurcation when a
quasistatic change in parameters causes the equilibrium to lose stability. Typically,
new equilibria and/or limit cycles appear from the nominal equilibrium at a bifur-
cation. When these so-called bifurcated solutions are stable, the bifurcation is said
to be supercritical or soft. When they are unstable, the bifurcation is said to be
subcritical or hard. The implication of these concepts for operation of real systems
is easy to describe (see [5] for a detailed exposition). For a subcritical bifurcation,
the system trajectory escapes from a local neighborhood of the nominal equilibrium
and can diverge to a distant (unacceptable) steady state. Therefore, detecting an
impending subcritical bifurcation is highly desirable. For a supercritical bifurcation,
on the other hand, the motion is likely to settle on the nearby bifurcated steady
state, which is stable and does not represent a major departure from the nominal
operation.
18
To flx ideas, consider a general one-parameter family of ordinary difierential
equation systems
_x = F?(x;?) (2.11)
where x 2Rn is the system state, ? 2Rdenotes the bifurcation parameter, and F?
is smooth in x and ?. For any value of ?, the equilibrium points of (2.11) are given
by the solutions for x of the algebraic equations F?(x;?) = 0. Local bifurcations are
those that occur in the vicinity of an equilibrium point. Local bifurcations of equi-
librium points consist of stationary bifurcations and the Andronov-Hopf bifurcation
(Hopf bifurcation for short). Stationary bifurcation is any bifurcation of one or
more equilibrium points from a nominal equilibrium point. Saddle-node, transcriti-
cal, and pitchfork bifurcations are some examples of stationary bifurcation. Unlike
the stationary bifurcation, in Hopf bifurcation a branch of periodic orbits bifurcates
from an equilibrium point. Stationary and Hopf bifurcation of equilibrium points
can be also categorized on the basis of the critical eigenvalues at the corresponding
bifurcation points. Next, a brief summary of various types of continuous-time local
bifurcations is given.
2.2.1 Stationary bifurcation
Stationary bifurcation occurs when a single real eigenvalue goes from being negative
to being positive as the bifurcation parameter ? passes through a critical value
?c. Precisely, the origin of (2.11) undergoes a stationary bifurcation at the critical
19
parameter value ? = 0 if assumptions, (L1), (L2) and (L3) hold.
(L1) F? of system (2.11) is su?ciently smooth in x, ?, and F?(0) = 0 for all ? in a
neighorhood of 0 (i.e., 0 is an equilbrium for all values of ? in the neighborhood);
(L2) The Jacobian A(?) := DxF?(0) possesses a simple real eigenvalue ?(?) such
that at the critical parameter value ?c = 0, ?(0) = 0 and
?0(0) := @?@?j?=0 6= 0;
(L3) All eigenvalues of the critical Jacobian A(0) besides 0 have negative real parts.
Under assumptions (L1), (L2) and (L3), two new equilibrium points of (2.11)
emerge from the origin at ? = 0. Bifurcation stability coe?cients are quantities
that determine the direction of bifurcation, and in particular the stability of the
bifurcated solutions. Locally, the new equilibrium points occur for parameter values
given by a smooth function of an auxiliary small parameter ? (? can be positive or
negative):
?(?) = ?1?+?2?2 +O(?3) (2.12)
One of the new equilibrium points occurs for ? > 0 and the other for ? < 0. Also,
the stability of the new equilibrium points is determined by the sign of an eigenvalue
fl(?) of the system linearization at the new equilibrium points. This eigenvalue is
near 0 and is also given by a smooth function of the parameter ?:
fl(?) = fl1?+fl2?2 +O(?3) (2.13)
20
Stability of the bifurcated equilibrium points is determined by the sign of fl(?).
If fl(?) < 0 the corresponding equilibrium point is stable, while if fl(?) > 0 the
equilibrium point is unstable. The coe?cients fli, i = 1;2;::: in the expansion
above are the bifurcation stability coe?cients mentioned earlier, for the case of
stationary bifurcation. The values of these coe?cients determine the local nature
of the bifurcation, as explained next.
Since ? can be positive or negative, it follows that if fl1 6= 0 the bifurcation
is neither subcritical nor supercritical. (This is equivalent to the condition ?1 6= 0.)
The bifurcation is therefore generically transcritical. Figure 2.1 shows the corre-
sponding bifurcation diagram. If fl1 = 0 and fl2 6= 0, a stationary bifurcation is
known as a pitchfork bifurcation. The pitchfork bifurcation is subcritical if fl2 > 0;
it is supercritical if fl2 < 0. Figure 2.2 shows the corresponding bifurcation diagram.
 
x 
?  
Figure 2.1: Transcritical bifurcation.
21
 
x 
?
(b) (a) 
x 
? 
Figure 2.2: Pitchfork bifurcation. (a) supercritical, (b) subcritical
2.2.2 Hopf bifurcation
Suppose that the origin of (2.11) loses stability as the result of a complex conjugate
pair of eigenvalues of A(?) crossing the imaginary axis. All other eigenvalues are
assumed to remain stable, i.e., their real parts are negative for all values of ?.
Under this simple condition on the linearization of a nonlinear system, the nonlinear
system typically undergoes a bifurcation. To be more precise, assume the following
conditions are satisfled (the critical parameter value is taken to be ?c = 0 without
loss of generality).
(H1) F? of system (2.11) is su?ciently smooth in x, ?, and F?(0;?) = 0 for all
? in a neighborhood of 0. The Jacobian A(?) possesses a complex-conjugate pair
of (algebraically) simple eigenvalues ?1;2(?) = fi(?) ? i!(?), such that fi(0) = 0,
22
fi0(0) 6= 0 and !c := !(0) > 0.
(H2) All eigenvalues of the critical Jacobian A(0) besides ?i!c have negative real
parts.
The Hopf Bifurcation theorem asserts that, under conditions (H1) and (H2), a small-
amplitude non constant limit cycle (i.e., periodic solution) of (2.11) emerges from
the origin at ? = 0. Locally, the limit cycles occur for parameter values given by a
smooth and even function of the amplitude ? of the limit cycles:
?(?) = ?2?2 +?4?4 +::: (2.14)
where the ellipsis denotes higher order terms. The stability of the limit cycle result-
ing from a Hopf bifurcation is determined by the sign of a particular characteristic
exponent fl(?). This characteristic exponent is real and vanishes in the limit as the
bifurcation point is approached. It is given by a smooth and even function of the
amplitude ? of the limit cycles:
fl(?) = fl2?2 +fl4?4 +::: (2.15)
The coe?cients ?2 and fl2 in the expansions above are related by the exchange of
stability formula
fl2 = ?2fi0(0)?2:
Generically, these coe?cients do not vanish. Their signs determine the direction of
bifurcation. The coe?cients fl2, fl4;::: in the expansion (2.15) are the bifurcation
23
stability coe?cients for the case of Hopf bifurcation. If fl2 > 0, then locally the
bifurcated limit cycle is unstable and the bifurcation is subcritical. If fl2 < 0,
then locally the bifurcated limit cycle is stable (more precisely, one says that it is
orbitally asymptotically stable). This is the case of supercritical Hopf bifurcation. If
it happens that fl2 vanishes, then stability is determined by the flrst non vanishing
bifurcation stability coe?cient (if one exists). Figure 2.3 shows a typical Hopf-
bifurcation diagram.
 
y
x 
?
(b) (a) 
y
x 
? 
Figure 2.3: Andronov-Hopf bifurcation. (a) supercritical, (b) subcritical
24
2.3 Eigenstructure Assignment for
Control System Design
Eigenstructure assignment problem is considered the general extension to the well
known eigenvalue assignment problem [23, 17]. This general notion is attributed
to the fact that, apart from the case of single-input single-output (SISO) systems,
the speciflcation of closed-loop eigenvalues does not uniquely deflne the feedback
structure of a closed-loop system. The source of nonuniqueness can be identifled as
that coming from the freedom ofiered by state or output feedback in multi-input
multi-output (MIMO) systems, beyond eigenvalue assignment, in selecting the as-
sociated eigenvectors and generalized eigenvectors from allowable spaces. Most of
the work on eigenstructure assignment has been purely eigenvalue assignment (pole
placement) up to late 60?s. Later, Moore found that degrees of freedom are available
over and above eigenvalue assignment using state feedback control for linear time-
invariant multi-input multi-output (MIMO) systems [17]. Since then, numerous
methods and algorithms have been developed to exercise those extra degrees of free-
dom to give the system some good performance characteristics. Good performance
characteristics span the areas of low sensitivity, robust, multi-objective and fault
detection [23]. In this section, background material on eigenstructure assignment
that is needed in the dissertation is recalled.
25
2.3.1 Eigenstructure assignment using state feedback
Consider the closed-loop LTI system
_x(t) = (A+BF)x(t) (2.16)
associated with an open-loop system
_x(t) = Ax(t)+Bu(t) (2.17)
to which a state feedback u(t) = Fx(t) is applied. Here x(t) 2 Rn; u(t) 2 Rm;
A 2R(n?n); F 2R(m?n). For a given complex scalar ?, which may be viewed as a
closed-loop eigenvalue, deflne the Hautus matrix [4]:
S? :=
?
(?I ?A) B
?
; S? 2Cn?(n+m)
Also deflne the matrix K?, with its row partitions along the same boundary as the
partitions of columns in S?. K? is given by:
K? :=
2
66
4
N?
M?
3
77
5; K? 2C(n+m)?m; (2.18)
and is composed of columns that constitute a basis for null space (or kernel) of S?.
The matrix K? can be calculated by performing the singular value decomposition
(SVD) or orthogonal triangular decomposition (OTD) on S? [23]. To compute K?
through SVD, apply singular value decomposition to S? to yield:
S? = U?V H; (2.19)
26
where U is an (n?n) unitary matrix over the fleld of complex numbers, the matrix
? is (n?(n+m)) with nonnegative numbers (singular values of S?) on the diagonal
and zeros ofi the diagonal, and V H denotes the conjugate transpose of V, an ((n+
m)?(n+m)) unitary matrix over the fleld of complex numbers. The matrix K? is
given by the last m columns of the unitary matrix V. Knowing K? means that for
any z 2C(m+n) such that
S?z = 0; (2.20)
there exists a k 2Cm such that
z = K? k: (2.21)
It is clear that the columns of N? are linearly independent if B has linearly
independent columns, and that N?? = N??. It is known that F may be chosen to
yield any self-conjugate set of closed-loop eigenvalues provided that the pair (A;B)
is controllable, i.e.,
rank [B AB ::: An?1B] = n
Moore in [17] described the available freedom for assigning eigenvectors for
an arbitrary self-conjugate set of eigenvalues using state feedback in the case m > 1.
He gave both necessary and su?cient conditions for a state feedback to exist that
achieves a desired set of eigenvalues and eigenvectors.
Proposition 2.1 [17] Consider f ?1;?2;:::;?n g a set of self-conjugate distinct com-
plex numbers, and a set of complex vectors fv1;v2;::::;vng. There exists a matrix F
27
of real numbers such that ?ivi = (A + BF)vi for all i 2 f1;2;:::;ng if and only if
the following three conditions are satisfled:
1. The vectors fv1;v2;:::;vng form a linearly independent set within Cn (over the
fleld of complex scalars),
2. vi = vj? whenever ?i = ??j, and
3. vi 2 span fN?ig,
where N?i is deflned as given in (2.18). Furthermore, if such a feedback gain matrix
F exists and rank B = m, then F is unique.
A proof of this proposition, along with a procedure for computing the feed-
backgain matrix F, is givenin [17]. Next, werecall Moore?s procedure for computing
the feedback gain matrix F [17].
Computing feedback gain matrix F
Suppose that vi, i 2f1;2;:::;ng are chosen to satisfy the three conditions of Propo-
sition 2.1. Since vi 2 span fN?ig for i 2 f1;2;:::;ng, then vi can be expressed
as
vi = N?i ki: (2.22)
28
for some vector ki 2Rm(Cm) (based on whether ?i is real (complex)), which implies
that
?
(?I ?A) B
?
2
66
4
vi
M?iki
3
77
5 = 0;
i.e.,
(?iI ?A)vi +BM?iki = 0:
If F is chosen such that ?M?iki = Fvi, then [?iI ?(A+BF)]vi = 0. Let
wi = ?M?i ki:
If all n eigenvalues are real numbers, then vi; wi are vectors of real numbers. Also
the flrst condition of Proposition 2.1 implies that the matrix [v1 v2 :::vn] is always
nonsingular. For this case,
F =
?
w1 w2 ::: wn
? ?
v1 v2 ::: vn
??1
(2.23)
For the case where there are complex eigenvalues, assume that ?1 = ??2. The second
condition of Proposition 2.1 states that v1 = v2? which implies that w1 = w2?. The
equation which must be solved is then:
F
?
v1R +jv1I v1R ?jv1I V
?
=
?
w1R +jw1I w1R ?jw1I W
?
(2.24)
where the columns of V and W are vi, i = 3;:::;n and wi, i = 3;:::;n, respectively.
Multiplying both sides of (2.24) from the right by the nonsingular matrix
29
2
66
66
4
1
2 ?j
1
2 0
1
2 +j
1
2
0 I
3
77
77
5
yields the equivalent equation
F
?
v1R v1I V
?
=
?
w1R w1I W
?
:
As vi, i 2 f1;2;:::;ng are independent, the columns of [v1R v1I V] are linearly
independent. Then F is given by:
F =
?
w1R w1I W
? ?
v1R v1I V
??1
: (2.25)
The third condition of Proposition 2.1 implies that the eigenvector vi must
be in the subspace spanned by the columns of N?i. This subspace is of dimension
m which is equal to the rank of B, or equivalently the number of independent
control variables. Therefore, the number of control variables available determines
how large (dimension) the subspace in which achievable eigenvectors must reside.
The orientation of the subspace is determined by the open-loop quantities A; B and
the desired closed-loop eigenvalue ?i. In general, however, a desired eigenvector vid
will not reside in the prescribed subspace and hence cannot be achieved. Instead
a \best possible" choice for an achievable eigenvector is made. Next, we recall the
analysis to compute the best possible achievable eigenvector via associated with a
desired eigenvector vid [8].
30
\Best possible via"
Thisbestpossibleachievableeigenvectorvia istheprojectionofvid ontothesubspace
spanned by the columns of N?i. To compute via associated with a real eigenvalue
?i, recall from (2.22) that
via = N?i ki; ki 2Rm:
To flnd the value of ki corresponding to the projection of vid onto the \achievable
subspace," ki is chosen to minimize
J = kvid ?viak2 = kvid ?N?i kik2:
Now
dJ=dki = 2 NT?i(N?i ki ?vid):
Hence dJ=dki = 0 implies
ki = (NT?i N?i)?1 NT?i vid (2.26)
via = N?i (NT?i N?i)?1 NT?i vid: (2.27)
For complex eigenvalues/eigenvectors the presentation above is formally correct with
minor arithmetic adjustments to accommodate complex quantities [17].
31
Chapter 3
Output Participation Factors
In this chapter, the basic deflnition of modal participation factors that was originally
given for LTI systems without input or output signals [2, 34] is revisited, along with
an extension proposed in [15], and then several new extensions are proposed. The
flrst such extension is a concept of output participation factors that measure the
participation of modes in output variables that are not necessarily system states.
The deflnition is based on computing an average relative contribution of a mode
in an output relative to uncertainty in the system initial condition. This follows
the approach proposed in [15] for deflning state participation factors. The new
deflnition proposed here has the desirable property that it reduces to the original
state participation deflnition when the outputs are simply the system states. Also,
this new deflnition is unit independent, i.e., the output participation factors do
not depend on the units of system states or outputs. It is also found that the
proposed output participation factors depend on the set to which the uncertain
32
initial condition belongs.
Another notion of modal participation is presented in this chapter. This
new deflnition is one of input-output participation factors. It involves probe input
signals as well as modes and outputs, and is applicable to a large set of systems. This
chapter proceeds as follows. In Section 3.1, we recall the work of [15] on deflning
participation factors using averaging over an uncertain initial condition. In Section
3.2, development of a notion of output participation factors is presented. In Section
3.3, input-to-output participation factors are introduced and studied.
3.1 Introduction
Consider a continuous linear time-invariant system described as given in (2.1). The
system solution for a general initial condition x0 (given in (2.5)), and the solution
for the particular initial condition x0 = ek (given in (2.7)) are repeated here for
convenience:
_x(t) = Ax(t) (3.1)
x(t) =
nX
i=1
(lix0)e?itri (3.2)
x(t) =
nX
i=1
(likrik)e?it
=
nX
i=1
pkie?it: (3.3)
33
Abed et al. [15] redeflned the participation factors concept of Verghese et al. [2, 34]
by considering the efiect of uncertainty in the initial condition x0. The reconsid-
eration of modal participation in this light was motivated by the following simple
observations for the linear system (3.1) and its solution (3.2): If the initial condi-
tion x0 lies along the i-th eigenvector, then the only mode that participates in the
evolution of any state is the i-th mode. On the other hand, if the initial condition
lies along the k-th coordinate axis, then the evolution of the k-th state involves all
system modes according to (3.3). Clearly, then, the degree to which a mode par-
ticipates in a state depends on the initial condition. The authors proposed a new
deflnition for state participation factors that reduces to the original deflnition of
participation factors of [2, 34] if the initial condition x0 is taken to lie in a connected
set ? containing the origin and obeys a symmetry condition deflned as follows:
Deflnition 3.1. ([15]) The set ? is symmetric with respect to each of the hy-
perplanes xk = 0, k = 1;2;:::;n if and only if for any k 2 f1; :::;ng and z =
(z1;:::;zk;:::;zn) 2Rn, z 2 ? implies that (z1;:::;?zk;:::;zn) 2 ?.
The symmetry assumption in Deflnition 3.1 is reasonable for typical engineering sys-
tem models although some applications might have restrictions on initial conditions.
The deflnition of state participation factors proposed in [15] is:
Deflnition 3.2. ([15]) The participation factor for the mode associated with ?i in
34
state xk with respect to an uncertainty set ? is
pki := avg
x02?
f(l
ix0)ri
k
x0k g (3.4)
whenever this quantity exists.
Here, avgx02? is an operator that computes the average of a function over the
set ?. This quantity measures the average relative contribution at time t = 0 of the
i-th mode to state xk . In the deflnition, the i-th mode is interpreted as the e?it
term in (3.2). Also, the denominator on the right-hand side of (3.4) is simply the
sum of the contributions of all system modes to xk(t) at t = 0:
x0k =
nX
j=1
(ljx0)rjk: (3.5)
In the next section, we consider the natural extension of the deflnition above to a
concept of participation of modes in outputs. The efiect of uncertainty in the initial
condition x0 is incorporated in a way similar to that used in [15] for calculation of
state participation factors.
3.2 Output Participation Factors
Consider a continuous linear time-invariant system with outputs
_x(t) = Ax(t) (3.6)
y(t) = Cx(t) (3.7)
35
where x 2Rn, y 2Rm, A is a real n?n system matrix and C is a real m?n output
matrix. Using the solution for x(t) in (3.2), the output vector is given by
y(t) = C
nX
i=1
(lix0)e?itri (3.8)
The evolution of the k-th output variable yk is given by
yk(t) = Ck
nX
i=1
(lix0)e?itri (3.9)
where Ck is the k-th row of C. It may be tempting to say that the participation
of the i-th mode in the k-th output is equal to the k-th row of C multiplied by the
state participation factors, i.e., the sum of the participation of the i-th mode in all
states weighted by the elements of the k-th row of C:
pyik =
nX
j=1
Ckj pji (3.10)
This temptation is supported by the fact that the expression in (3.10) reduces to the
participation of state k in mode i when Ck = ek (the k-th unit vector). However,
this would lead to a notion that is unit-dependent and hence unacceptable from the
point of view of this dissertation.
Next, we give a high-level deflnition that is the natural extension of the
state participation factors deflnition in [15] to the case of participation of modes
in outputs. It is easily verifled that this deflnition is independent of units. Below,
the concept is reduced to a simple calculation for systems with initial condition
uncertainty known to lie in an n-ellipsoid or n-rectangle with main axes parallel to
the coordinate axes.
36
Deflnition 3.3. The participation factor for the mode associated with ?i in output
yk with respect to an uncertainty set ? is
pyik := avg
x02?
f(l
ix0)(Ckri)
y0k g (3.11)
whenever this quantity exists.
Here, the expression avgx02? denotes the operator that computes the average of a
function over the set ?.
Note that in computing the implied multidimensional volume integral, the
argument function is undeflned for y0k = 0, i.e., undeflned for x0 2 fx : Ckx = 0g
and the Cauchy principal value of the integral is to be used. The quantity on the
right-hand side of Deflnition 3.3 measures the average relative contribution at time
t = 0 of the i-th mode to output yk . The denominator of this quantity is simply
the sum of the contributions of all modes to yk(t) at t = 0:
y0k =
nX
j=1
(ljx0)Ckrj
= Ckx0:
Writing li = liCk + li?Ck, where liCk is the projection of li in the direction of Ck and
li?Ck is the projection of li in the direction orthogonal to Ck, the expression in (3.11)
becomes
pyik = Ckri avg
x02?
f(l
i
Ck +l
i?
Ck)x
0
Ckx0 g
= Ckri avg
x02?
fl
i
Ckx
0
Ckx0 +
li?Ckx0
Ckx0g (3.12)
37
Writing liCk in terms of li and Ck,
liCk = klik?f l
iCkT
klikkCkkg?f
Ck
kCkkg
= fl
iCkT
kCkk2g?C
k (3.13)
Using (3.13), pyik reduces to
pyik = Ckri avg
x02?
fl
iCkT
kCkk2 +
li?Ckx0
Ckx0g
= C
kriliCkT
kCkk2 +C
kri avg
x02?
fl
i?
Ckx
0
Ckx0g
Let
Vol(?) :=
Z
x02?
dx0 (3.14)
denote the volume of the set ?. The expression on the right-hand side of (3.11) is
pyik = C
kriliCkT
kCkk2 +C
kri
Z
x02?
fl
i?
Ckx
0
Ckx0gdx
0=Vol(?) (3.15)
The orthogonal vectors li?Ck and Ck operate on the initial state vector x0 2 ? to
orthogonally transform it to another set ?0. Generally, even if the set ? is symmetric
according to Deflnition 3.1, the new set ?0 is not symmetric with respect to the axes
x0i = 0, i 2 f1;2;:::;ng, and the integral term doesn?t vanish. Note that in the
case where output yk is simply the k-th state (i.e., Ck = ek), the integral term on
the right hand side of (3.15) vanishes. This follows from the observation that the
integral term reduces to
g ?
Z
x02?
fx
0
j
x0kgdx
0
38
for some j 6= k and a constant g = Ckri=Vol(?). Because ? is symmetric according
to Deflnition 3.1,
Z
x02?
fx
0
j
x0kgdx
0 = 0
where the integral is interpreted in the sense of the Cauchy principal value. The
expression in (3.11) is dimensionless, i.e., it does not depend on the units in which
state variables are measured or the unit in which the output variable is measured.
Also, its simpliflcation (3.15) reduces to the conventional state participation factors
when the output is simply a state variable. Thus, it can be used to provide a good
measure of the the mode that has the highest participation in an output when it is
reasonable to assume uncertainty in the system?s initial conditions. Also note that
the deflnition in its current form is dependent on the set ? to which x0 belongs.
Suppose that the initial condition uncertainty set ? is known to be symmetric
according to Deflnition 3.1. This is not enough to ensure that the integral term in
(3.15) vanishes. However, this integral vanishes if the initial condition uncertainty
set is an n-cube or an n-sphere. For instance, the initial state z0 lies in an n-
sphere if kz0k ? c0 where c0 is a positive constant. Therefore, suppose the initial
condition uncertainty set is an n-ellipsoid or n-rectangle with main axes parallel to
the coordinate axes. Next, consider a change in units (x0 = Tx), T is a positive
diagonal matrix, such that in the new coordinates the state initial condition x00 lies
either in an n-cube or an n-sphere. In this case, the integral term vanishes (a proof
39
is given below) and the expression in (3.15) reduces to
pyik = C
0kr0il0iC0kT
kC0kk2 ; (3.16)
where
C0k = CkT?1; r0i = Tri and l0i = liT?1:
Proof:
Let x00 = Qz0, where Q is an (n?n) orthonormal matrix, i.e.,
Q?1QT = I
The orthonormal similarity transformation matrix Q can be written as
Q =
2
66
4
^C0k
Q1
3
77
5
T
where
^C0k = C0k
kC0kk;
and Q1 is chosen such that C0kQT1 = 0. Note that kx00k = kz0k, dx00 = dz0. Using
the orthonormal similarity transformation Q, the integral term in (3.15) becomes
Z
x00?c0
fl
0i?
C0k x
00
C0k x00 gdx
00 =
Z
z0?c0
fl
0i?
C0k Q z
0
C0k Q z0 gdz
0
=
Z
z0?c0
fl
0i?
C0k ( ^C0
k)T z0
1 +l
i?
C0k Q
T
1 z
00
C0k ( ^C0k)T z01
gdz0
=
Z
z0?c0
f l
0i?
C0kQ
T
1 z
00
C0k ( ^C0k)T z01
gdz0
= 0; (3.17)
40
where z0 = [z01 z00]T. The last step in (3.17) follows from the observation that for
any j 6= 1,
Z
z0?c0
fz
0
j
z01gdz
0 = 0;
where the integral above is interpreted in the sense of Cauchy principal value. ?
Eq. (3.16) is taken as the deflnition of output participation factors, and is
applied to the transformed system in which the initial condition uncertainty set is
an n-cube or n-sphere.
The expression in (3.16) is crucial for computing the value of the high-level
deflnition of output participation factors given in Deflnition 3.3. This value (that
depends on the n-ellipsoid or n-rectangle set ? with main axes parallel to the co-
ordinate axes) is equal to the value of the deflnition of output participation factors
in (3.16) when a diagonal transformation (equivalent to a change in units of system
states) is performed to transfer the set ? into an n-cube or an n-sphere.
Next, we summarize the properties of output participation factors (deflned
in (3.16)). These properties can be easily proved through appropriate consideration
of Deflnition 3.3 and Eq. (3.16).
Properties of output participation factors
1. The value of the high-level deflnition of output participation factors given in
Deflnition 3.3, and which is a function of A; C and the set ?, is equal to the
value of the deflnition of output participation factors given in (3.16) evaluated
41
for ^A and ^C, where
^A = TAT?1
^C = CT?1;
where T is the diagonal transformation matrix required to transform the sym-
metric n-ellipsoid or n-rectangle set ? into a an n-cube or an n-sphere. Figure
3.1 depicts such a diagonal transformation T for the case of a two-dimensional
system.
2. The sum of output participation factors associated with any mode is equal to
1, i.e.,
nX
i=1
CkriliCkT
kCkk2 = 1
3. The output participation factors are dimensionless quantities that are inde-
pendent of the units in which state and/or output variables are measured.
Next, a numerical example is presented to show that the values of the output
participation factors computed according to the expression in (3.10) do not agree
with the values of output participation factors (deflned in (3.16)).
42
 
T 
z1 
z2 x2 
x1 
T 
z1 
z2 x2 
x1 
(a) 
(b) 
Figure 3.1: Illustrating diagonal transformation T in two dimensions, (a) Elliptic
original uncertainty set ?, (b) Rectangular original uncertainty set ?.
Numerical example
Consider an LTI system given by (3.6){(3.7) with
A =
0
BB
BB
BB
@
0 1 0
?1:5 ?0:5 1:225
?0:433 0 0
1
CC
CC
CC
A
and C =
0
BB
@
0 1 1
0 2 1
1
CC
A
The initial conditions of this system are uncertain and lie in a set ?. We assume
that ? is an n-sphere of radius c0, i.e., kx0k ? c0 (c0 is a positive constant). Thus,
no transformation is required and the output participation factors can be computed
using A and C according to the deflnition in (3.16). The eigenvalues of A are
43
?1;2 = ?0:0672?j1:2026, ?3 = ?0:3655. We are interested in the participation of
the critical mode (corresponding to the complex conjugate pair of eigenvalues ?1;2)
in the states and in the outputs. Right and left eigenvectors corresponding to ?1
that satisfy the normalization condition (2.2) are, respectively,
r1 =
?
?0:0348?j0:6216 0:7499 0:2224?j0:0249
?T
l1 =
?
?0:0078+j0:8297 0:6087+j0:2255 0:1944?j0:6308
?
:
The participations of mode 1 in states xi; i = 1;2;3 are p11 = 0:5160 ? j0:0240
(jp11j = 0:5166), p21 = 0:4565+j0:1691 (jp21j = 0:4868) and p31 = 0:0275?j0:1451
(jp31j = 0:1477). The participation of mode 1 in the outputs y1 and y2 can be calcu-
lated using (3.16): py11 = 0:3854?j0:2071 (jpy11j = 0:4375) and py21 = 0:4854?j0:0690
(jpy21j = 0:4903). A simple calculation shows that participation of mode 1 in either
output is not the weighted sum of the state participations (as suggested by the ex-
pression given in (3.10)). Next, we present another notion of modal participation
factors. This new deflnition is one of input-to-output participation factors. It in-
volves probe input signals as well as modes and outputs, and is applicable to a large
set of systems.
44
3.3 Input-to-Output Modal Participation Factors
Consider a linear time-invariant system with inputs and outputs:
_x(t) = Ax(t)+Bu(t) (3.18)
y(t) = Cx(t): (3.19)
The sought deflnition of input-to-output participation factors is considered a natural
extension to the deflnition of input-to-state participation factors [40] given in Section
2.1.2. Recall that the deflnition of input-to-state participation factors assumes a
special form of the input and output matrices B and C, respectively. It is assumed
that the input is applied to one component, say the q-th component, of the right side
of (3.18) and only one state, say the k-th state, is measured. For cases of systems
where availability of measurements of states is not usually possible, and/or inputs
are applied to more than one component of the right side of (3.18), extending the
deflnition in (2.10) to address such cases is needed. We therefore, consider the case
of systems with general Bq and Ck, where Bq and Ck are the q-th column of input
matrix B and k-th row of output matrix C, respectively. This means that the input
is applied through the q-th input channel of the right side of (3.18) and the k-th
output of (3.19) is measured.
With this choice of the input location and the measured output, in steady
state the k-th output in (3.19) (in the frequency domain) is given by
45
yk(s) = Ck(sI ?A)?1Bqu
=
nX
i=1
CkriliBq
s??i u(s)
=
nX
i=1
Rikq
s??iu(s)
Let
Rikq := CkriliBq:
Take
pyiqk = jRikqj (3.20)
as the participation factor of mode i in output k when the input is applied through
the q-th input channel. We call this quantity the input-to-output participation factor
(IOPF) for mode i, with measurement at output k and input applied through the
q-th input channel. Note that the IOPF, similar to the ISPF, is dimensionless. Also,
note that, the quantity Rikq = CkriliBq is numerically equal to the entry on the k-th
row and q-th column of the residue matrix associated with mode i. This deflnition
of input-to-output participation factors is useful in designing monitoring systems to
detect impending instability in nonlinear systems as discussed in the next chapter.
46
Chapter 4
Detecting Impending Instability in
Nonlinear Systems
In this chapter, a signal-based approach for real-time detection of impending insta-
bility in nonlinear systems is considered. The main idea pursued here involves using
a small additive white Gaussian noise (AWGN) as a probe signal and monitoring
the spectral density of one or more measured states or outputs for certain signatures
of impending instability.
Input-to-output participation factors introduced in Chapter 3 along with
Input-to-state participation factors given in [40] are used as a tool to aid in selection
of locations for probe inputs and outputs or states to be monitored. Since these
participation factors are model-based, the chapter combines signal-based and model-
basedideastowardachievingarobustmethodologyforinstabilitymonitoring. Power
systems will be the main application of this methodology where the detection of
impending instabilities is crucial.
47
The chapter proceeds as follows. In Section 4.1, the sought signal-based
approach to instability monitoring using input-to-state and/or input-to-output par-
ticipation factors is presented. A stability index based on power spectral density
measurements is also given. In Section 4.2, three case studies from power systems
are given to demonstrate the proposed approach to instability monitoring.
4.1 Precursor-Based Monitoring
In this section, we show that noisy precursors can be used as a warning signal
indicating that a system is operating dangerously close to instability. We also show
that the spectrum of a measured output (state) of the system is proportional to the
square of the input-to-output (input-to-state) participation factors. Thus, IOPFs
(ISPFs) can be used to determine the best location for applying the probe signal and
for choosing which output (state) to measure where the noisy precursor would be
most apparent. A stability index based on power spectral density measurements that
can be used to aid operators in assessing proximity to the approaching instability is
given. Figure 4.1 shows a schematic diagram of our instability monitoring system.
The Figure depicts the idea of having the monitoring system to trigger a preventive
control action to protect the system against the efiects of the approaching instability.
One possible preventive action is discussed in Chapter 5.
48
Measurements
Control
Power System
Probe Signal
Inputs
Instability Detector
Signal Analysis/
Preventive 
Figure 4.1: Precursor-based instability monitor with external probe signal.
4.1.1 Developed monitoring system
Consider a nonlinear dynamic system (\the plant")
_x = f(x;?)+?(t) (4.1)
where x 2 Rn, ? is a bifurcation parameter, and ?(t) 2 Rn is a zero-mean vector
white Gaussian noise process [25]. Let the system possess an equilibrium point x0.
For small perturbations and noise, the dynamical behavior of the system can be
described by the linearized system in the vicinity of the equilibrium point x0. The
linearized system corresponding to (4.1) with a small noise forcing ?(t) is given by
_x = Df(x0;?)x+?(t) (4.2)
49
where Df(x0;?) represents the gradient of f(x;?) with respect to x evaluated at
x = x0. In (4.2), x now denotes x ? x0 (the state vector referred to x0). For the
results of the linearized analysis to have any bearing on the original nonlinear model,
we must assume that the noise is of small amplitude.
The noise ?(t) can occur naturally or can be injected using available controls.
We consider the case where the noise is applied through one input channel and the
power spectral density of one output is calculated. That is, we consider the case
where ?(t) = Bq?(t) where Bq is the q-th column of the input matrix B, ?(t) is a
scalar white Gaussian noise with zero mean and power  2, applied through the q-th
input channel and the output is given by y = Cpx where Cp is the p-th row of the
output matrix C.
In steady state, the p-th output of system (4.2) forced by a small AWGN is
given by
yp(s) =
nX
i=1
CpriliBq
s??i ?(s)
=
nX
i=1
Ripq
s??i?(s)
= H(s) ?(s):
The absolute value of Ripq was deflned in Chapter 3 as the participation factor of
mode i in output p when the input is applied through the q-th input channel. The
power spectral density (PSD) of the output of a linear system with transfer function
50
H(j!) is related to the PSD of the input by [24]
Sy(!) = H(j!)H(?j!)S?(!)
Thus, the power spectrum of the p-th output is given by
Syp(!) =
? nX
i=1
Ripq
j! ??i
!? nX
k=1
Rkpq
?j! ??k
!
 2
=  2
nX
i=1
nX
k=1
Ripq
j! ??i
Rkpq
?j! ??k
Suppose that the system is nearing a Hopf bifurcation. Speciflcally, assume that a
complex conjugate pair of eigenvalues is close to the imaginary axis, and has rela-
tively small negative real part in absolute value compared to other system eigenval-
ues. Denote this pair as ?1;2 = ???j!c, with ? > 0 small and !c > 0:
jRe(?i)j ?; i = 3;:::;n:
Under this assumption, and for values of ! close to !c, Syp(!) can be approximated
as
Syp(!) ?  2
2X
i=1
2X
k=1
Ripq
j! ??i
Rkpq
?j! ??k
=  2
 jR
1pqj2
?2 +(! ?!c)2 +
jR1pqj2
?2 +(! +!c)2
+ 2Re
? 1
?+j(! ?!c)
(R1pq)2
??j(! +!c)
 ?
:
Note that all terms containing ?i, i = 3;:::;n have been neglected and only terms
containing the critical eigenvalues ?1 and ?2 have been retained. After algebraic
manipulation and substituting (R1pq)2 = fi + jfl where fi = jR1pqj2 cos(2 pq) and
51
fl = jR1pqj2 sin(2 pq), with  pq = tan?1
?
ImfR1pqg
RefR1pqg
?
, the PSD of yp can be rewritten
as
Syp(!) =  2jR1pqj2
 1
?2 +(! ?!c)2 +
1
?2 +(! +!c)2
?
+  2(fl?+fi!c)(! ?!c)+?(?fi?!cfl)(?2 +!2
c)(?2 +(! ?!c)2)
?  2(fl?+fi!c)(! +!c)??(?fi?!cfl)(?2 +!2
c)(?2 +(! +!c)2)
=  2jR1pqj2
?
(1+G1(!)) 1?2 +(! ?!
c)2
+ (1?G2(!)) 1?2 +(! +!
c)2
?
(4.3)
where
G1(!) = 1?2 +!2
c
[(?sin(2 pq)+!c cos(2 pq))(! ?!c)
+ ?(?cos(2 pq)?!c sin(2 pq))];
G2(!) = 1?2 +!2
c
[(?sin(2 pq)+!c cos(2 pq))(! +!c)
? ?(?cos(2 pq)?!c sin(2 pq))]:
For ! = !c and su?ciently small ? (? ? !c), the power spectral density of yp is
given by
Syp(!c) =  2jR1pqj2
 1
?2 +O
 1
?
?
+O(1)
?
:
Note that the IOPFs are related to the spectral densities of the outputs of a system
driven by small AWGN as in (4.3). The amplitude of the spectrum is proportional
to the square of the IOPFs. The input-to-output participation factors can be used
52
to determine the best location for applying the probe signal and also the output
that will have the highest spectral peak.
The special case where only system states are monitored and input is applied
to one component of the right side of (4.2) bears similar results [21]. In [21], we show
that the ISPFs are related to the spectral densities of the states of a system driven
by small AWGN. The amplitude of the spectrum is proportional to the square of the
ISPFs. Thus, the input-to-state participation factors can also be used to determine
the best location for applying the probe signal and also the state that will have the
highest spectral peak.
4.1.2 Instability proximity index
In this section, an instability proximity index that helps predict closeness to insta-
bility based on power spectral density measurements is given. Performance indices
have been used to predict proximity to voltage collapse in power systems [44]. Sev-
eral stability indices have been proposed in the literature, some of these indices are
model-based while others rely mainly on online measurements and operators expe-
rience. Sensitivity factors, singular values and eigenvalues are among those indices
used in the literature [44].
In this paper, a sensitivity-based index is proposed. This index is based on
online measurements of power spectral density peaks at certain critical frequencies.
53
The proposed index is given by
PSDI = dfPSDvgd? : (4.4)
Here, PSDI stands for power spectral density index, PSDv is the power spectral
density peak value at the critical frequency of the measured variable v (state or
output) corresponding to some value of the bifurcation parameter ?. The PSDI
value at time t is approximated using ?fPSDvg=?f?g at time t, where ?fPSDvg
and ?f?g are the changes in the values of PSDv and ? between time instances
t and t ? ?t, respectively. The sampling time ?t is determined based on oper-
ator?s experience and system?s conditions. The instability index is calculated for
the variable v (state or output) which has the highest ISPF or IOPF, respectively.
The instability index grows signiflcantly as the system approaches criticality, where
it is theoretically inflnity at criticality. The reciprocal of the index is more infor-
mative than the index itself as the reciprocal value approaches zero as the system
approaches criticality. This index can be used to assist a system operator or an
automatic controller in taking a preventive action. For example, a threshold can be
used such that if the reciprocal of the index drops below that threshold, a preventive
action is triggered. The efiectiveness of the proposed instability index in detecting
closeness to instability in a single generator power system model is demonstrated in
Section 4.2.2.
54
4.2 Power System Examples
In this section, we apply results from the previous section to monitor and detect
impending instability of examples from power systems. These examples show how
ISPF and IOPF can be used to better monitor an impending instability.
Although white Gaussian noise was used to derive the results presented so
far, band-limited white Gaussian noise signals will be used to probe the monitored
power system examples. This is because it is practically impossible to generate
white noise signals. Band-limited white noise is a useful theoretical approximation
when the noise disturbance has a correlation time that is very small relative to the
natural bandwidth of the system [43].
4.2.1 Single generator with dynamic load
Consider the single generator power system model with induction motor load [12]:
_?m = ! (4.5)
M _! = ?dm! +Pm ?EmVYm sin(?m ??) (4.6)
Kqw _? = ?Kqv2V 2 ?KqvV +Q(?m;?;V)?Q0 ?Q1 (4.7)
TKq!Kpv _V = Kp!Kqv2V 2 +(Kp!Kpv ?Kq!Kpv)V
+ Kq!(P(?m;?;V)?P0 ?P1)
? Kp!(Q(?m;?;V)?Q0 ?Q1) (4.8)
55
The state variables are ?m (the generator phase angle, closely related to the me-
chanical angle of the generator rotor), ! (the rotor speed), ? (the load voltage phase
angle) and V (the magnitude of the load voltage). The load includes a constant PQ
load in parallel with an induction motor. The real and reactive powers supplied to
the load by the network are
P(?m;?;V) = ?E00VY 00 sin(?)+EmVYm sin(?m ??);
Q(?m;?;V) = E00VY 00 cos(?)+EmVYm cos(?m ??)
? (Y 00 +Ym)V 2
where
E00 = E0p1+C2Y ?2
0 ?2CY
?1
0 cos 0
Y 00 = Y0
q
1+C2Y ?20 ?2CY ?10 cos 0
 00 =  0 +tan?1
? CY ?1
0 sin 0
1?CY ?10 cos 0
 
The values of the parameters for this model are:
M = 0:01464, C = 3:5, Em = 1:05, Y0 = 3:33,  0 = 0,  m = 0, Kp! = 0:4,
Kpv = 0:3, Kq! = ?0:03, Kqv = ?2:8, Kqv2 = 2:1, T = 8:5, P0 = 0:6, P1 = 0:0,
Q0 = 1:3, E0 = 1:0, Ym = 5:0, Pm = 1:0, dm = 0:05. All values are in per unit
except for angles, which are in degrees.
It has been shown that a supercritical Hopf bifurcation occurs in this power
system model as the reactive load Q1 is increased through the critical value Q?1 =
2:9801438 [12]. Next, we consider the system operating at loads close to the Hopf
56
bifurcation, say at Q1 = 2:9. The corresponding operating point is x0 = [0:2473,
0; 0:0398; 0:9248]. The Jacobian of the system at this operating point is
A =
2
66
66
66
66
66
4
0 1 0 0
?324:5254 ?3:4153 324:5254 ?73:8611
33:3333 0 ?29:2479 72:7220
?3:3656 0 1:5180 ?11:1529
3
77
77
77
77
77
5
The eigenvalues of A are f?0:7923?j6:6318; ?21:1157?j10:9959g.
To monitor the system, a band-limited AWGN probe signal is applied to the
mechanical power Pm. Figure 4.2 depicts the spectral densities for the four states
?m, !, ? and V for Q1 = 2:9 and  = 0:001. As it is clear from this flgure, the state
! has a higher peak than all other states. Figure 4.3 demonstrates the variation of
the spectral density peak near ! = !c ? 6:6 rad=s as a function of the bifurcation
parameter Q1. The values of the input-to-state participation factors of the critical
mode in all states are given in Table 4.1. As predicted by the analysis in Section 4.1,
the ordering of the peaks of the spectral densities of all states at !c can be predicted
from the values of the ISPFs.
4.2.2 Single generator connected to an inflnite bus
Consider a synchronous machine connected to an inflnite bus together with excita-
tion control [14]. It was shown in [14] that this system undergoes a Hopf bifurcation
as the control gain in the excitation system is increased beyond a critical value. The
57
0 5 10 15 20 25 30 35 40
10?7
10?6
10?5
10?4
10?3
10?2
10?1
Freq [rad/s]
Power spectral density
?m
??
V
Figure 4.2: Power spectral densities of the states of the model given in (4.5)-(4.8).
The bifurcation parameter was set to Q1 = 2:9. Band-limited AWGN of zero mean
and (0:001)2 power was added to Pm.
58
2.5 2.6 2.7 2.8 2.9 30
0.02
0.04
0.06
0.08
0.1
0.12
Q1
Value of PSD of 
? at critical freq 
? c
Figure 4.3: Variation of the peak value of the power spectral density of ! as a
function of the bifurcation parameter Q1. Band-limited AWGN of zero mean and
(0:001)2 power was added to Pm.
dynamics of the generator are given by:
_? = ! (4.9)
2H _! = ?D! +!0(Pm ?Pe) (4.10)
?0d0 _E0q = EFD ?E0q ?(Xd ?X0d)id (4.11)
59
States Spectral peak Input-to-state participation
at !c ? 6:6318 factors (ISPFs)
?m 9:528?10?4 p121 = 2:8992
! 40:38?10?4 p122 = 19:364
? 6:616?10?4 p123 = 2:4013
V 0:305?10?4 p124 = 0:4981
Table 4.1: Input-to-state participation factors and spectral peaks at !c.
with the following algebraic equations:
Pe = Eqiq
Eq = E0q +(Xq ?X0d)id
id = x(Eq ?Ecos?)?rEsin?
iq = r(Eq ?Ecos?)+xEsin?
x = Xl+XqR2
l+(Xl+Xq)2
r = RlR2
l+(Xl+Xq)2
The subscripts d and q refer to the direct and quadrature axes, respectively. The
dynamics of the excitation control are given by
?E _EFD = ?KEEFD +VR ?EFDSE(EFD) (4.12)
?F _V3 = ?V3 + KF?
E
(?KEEFD +VR ?EFDSE(EFD)) (4.13)
?A _VR = ?VR +KA(VREF ?Vt ?V3) (4.14)
Here Vt is the terminal voltage and is given by
V 2t = v2d +v2q
60
where
?vd = ?q = ?Xqiq
vq = ?d = E0q ?X0did:
The saturation function SE(EFD) is usually approximated as
SE(EFD) = AEX exp(BEXEFD):
An equilibrium point of this system is denoted by x0 = (?0;!0;E0q0;E0FD;V 03 ;V 0R).
The values of the parameters that appear in this power system model are given in
Table 4.2.
Synchronous machine Exciter Transmission line
H = 2:37 s KE = ?0:05 R0l = 0:02
D = 1 pu KF = 0:02 X0l = 0:40
Xd = 1:7 ?E = 0:50 s Rl = ?R0l
X0d = 0:245 ?F = 0:60 s Xl = ?X0l
Xq = 1:64 ?A = 0:10 s
!0 = 377:0 rad/s AEX = 0:09
?0d0 = 5:9 s BEX = 0:50
Table 4.2: Parameter values for the single generator connected to an inflnite bus
model.
For Pm = 0:937, VREF = 1:130, and ? = 2, it has been shown that a
subcritical Hopf bifurcation occurs at K?A = 193:74 [14]. Next, we consider the
system operating before the Hopf bifurcation, say at KA = 185. The corresponding
operating point is given by x0 = [1:3515, 0, 1:1039, 2:3150, 0, 0:5472]T. The Jacobian
of the system at this operating point is
61
A =
2
66
66
66
66
66
66
4
0 1 0 0 0 0
?62:2 ?0:2 ?79:7 0 0 0
?0:2 0 ?0:4 0:2 0 0
0 0 0 ?1:1 0 2
0 0 0 0 ?1:7 0:1
125:9 0 ?1157:6 0 ?1850 ?10
3
77
77
77
77
77
77
5
:
The eigenvalues of A are f?0:0139 ? j7:7707; ? 4:5832 ? j12:6178;?2:1029 ?
j0:9417g.
Note that for this model, there are two physically feasible locations for ap-
plying the probe signal. The probe signal can be either applied to Vref or to Pm.
The input-to-state participation factors are used to determine the best location for
applying the probe signal. From the values of the ISPFs (see Table 4.3), it is clear
that mode 1 has higher participation in other states when the probe signal is applied
to Pm than when applied to Vref. This can be also seen from the power spectral
densities shown in Figures 4.4 and 4.5. Also, the ISPFs give an indication of which
state to monitor. The higher the participation factor of the critical mode in a state,
the higher the peak of the spectrum for that state. Figure 4.6 shows the variation
of the power spectral peak at the critical frequency as a function of the bifurcation
parameter when noise is added to Pm. Figure 4.7 shows the inverse of the PSDI
based on measurements from the state VR versus the exciter gain KA. Note that,
for example, a threshold value of 5 can be used such that if the reciprocal of the
index drops below that value, a preventive action is triggered.
62
State Spect. peak at !c ISPFs Spect. peak at !c ISPFs
(noise added to Pm) (noise added to Vref)
? 0.0226 p121 = 0:0648 0.0019 p161 = 0:0024
! 1.2880 p122 = 0:4923 0.1084 p162 = 0:0185
E0q 0:75938?10?4 p123 = 0:0038 0:68651?10?5 p163 = 0:0001
EFD 0.2326 p124 = 0:2084 0.0210 p164 = 0:0078
V3 2:4644?10?4 p125 = 0:0068 2:2288?10?5 p165 = 0:0003
VR 3.3923 p126 = 0:8006 0.3064 p166 = 0:0301
Table 4.3: Input-to-state participation factors and spectral peaks at !c ? 7:8 for
the single generator connected to an inflnite bus system.
4.2.3 Three-generator nine-bus power system
Below, we consider the Western System Coordinating Council (WSCC) 3{machine,
9{bus power system model, which is widely used in the literature [33, pp. 170{
177],[16]. The dynamics of this model includes three identical IEEE-Type I exciters
for the three machines. The machines and exciters data are given in Tables 4.4 and
4.5, respectively [33, 16].
In this model, a subcritical Hopf bifurcation occurs as the load on bus 5 is
increased beyond 4.5 pu [33]. Our goal in this case study is to detect this impending
loss of stability by using a band-limited AWGN probe signal and continuously mon-
itoring the power spectral densities of certain states or outputs. This would give
the system operator (or an automatic controller) valuable time to take appropriate
preventive measures (e.g., shedding loads at certain buses or trigger a controller
63
2 4 6 8 10 12 14
10?10
10?8
10?6
10?4
10?2
100
Freq [rad/sec]
Power spectral density
??
Eq?E
FDV
3V
R
Figure 4.4: Power spectral densities of the states of the single generator connected
to an inflnite bus system. The bifurcation parameter was set to KA = 185. Band-
limited AWGN of zero mean and (0:000032)2 power was added to Pm.
that assigns the participation factors of system states in the critical mode as will
be discussed in Section 6.2). The simulations of this model were conducted us-
ing the software package PSAT [16]. For values of the load on bus 5 close to 4.0
pu, the linearization of the system at the operating point has two complex con-
jugate pair of eigenvalues close to the imaginary axis, ?1;2 = ?0:17665 ? j8:184
and ?3;4 = ?0:3134?j1:7197. As the load on bus 5 is increased further, the pair
?3;4 approaches the imaginary axis, while the other pair ?1;2 changes only slightly.
For example, when the load at bus 5 is 4.4 pu, ?1;2 = ?0:18231 ? j8:0978 and
?3;4 = ?0:04602?j2:1151. Increasing the load on bus 5 beyond 4.5 pu causes the
pair ?3;4 to cross the imaginary axis from left to right.
64
2 4 6 8 10 12 14
10?8
10?6
10?4
10?2
100
Freq [rad/sec]
Power spectral density
??
E?qE
FDV
3V
R
Figure 4.5: Power spectral densities of the states of the single generator connected
to an inflnite bus system. The bifurcation parameter was set to KA = 185. Band-
limited AWGN of zero mean and (0:000032)2 power was added to Vref.
From the values of the ISPFs calculated for this system, we found that both
of the critical modes have higher participation when the probe signal is applied to
Pm3, the mechanical power of generator number 3. Also, we found that these modes
have high participation in the fleld voltage of the exciters. Therefore, in the following
simulations, the probe signal is added to Pm3 and the power spectral densities of the
fleld voltages of the three exciters (i.e., Efdi;i = 1;2;3) are monitored. Table 4.6
shows input-to-state participation factors for the 3-machine nine-bus system (partial
listing) when the load at bus 5 is 4.4 pu. Figure 4.9 and Figure 4.10 show the power
spectral densities of Efdi, i = 1;2;3 when the load on bus 5 (PL5) is 4.0 pu and
4.4 pu, respectively. It is clear from Figure 4.9 that when the load on bus 5 is
65
170 175 180 185 1900
0.2
0.4
0.6
0.8
1
1.2
1.4
KA
Value of PSD of V
R at critical freq 
? c
Figure 4.6: Variation of the peak value of the power spectral density of VR as a
function of the bifurcation parameter KA. Band-limited AWGN of zero mean and
(0:000032)2 power was added to Pm.
4.0 pu, the spectrum has two peaks at 0.28 Hz and 1.3 Hz. These two frequencies
correspond to the complex eigenvalues ?3;4 and ?1;2, respectively. Note that the
peak at 1.3 Hz that corresponds to the pair of complex eigenvalues ?1;2 is higher
than the peak at 0.28 Hz. However, when the load at bus 5 is increased to 4.4 pu,
the peak at 0.28 Hz becomes much larger than the one at 1.3 Hz (see Figure 4.10),
which is an indicator that an instability is being approached. Figure 4.11 shows the
power spectral density of Efd1 for three values of PL5: 4.0 pu, 4.25 pu and 4.4 pu.
Next, we consider using IOPFs to decide which is the best output to mea-
sure. Table 4.7 shows partial listings of the IOPFs when the load at bus 5 is 4.4
66
178 180 182 184 186 188 190 192 1940
10
20
30
40
50
60
70
80
Bifurcation parameter Ka
1/PSDI
 
Thershold value=5 
A 
Figure 4.7: Plot of the inverse of the power spectral density index of ! versus the
bifurcation parameter KA. The threshold value is equal to 5.
pu. The values listed correspond to cases where the probe signal is applied to
one of the mechanical powers of the generators fPm1; Pm2;Pm3g and one of the
outputs fV5;V6;V7;Q14;Q39;Q27g is measured. fV5;V6;V7g are voltages at buses
f5;6;7g, respectively and fQ14;Q39;Q27g are reactive powers  owing between buses
f(1;4);(3;9);(2;7)g, respectively. Figure 4.12 shows the power spectral density for
state Efd1 and outputs Q14 and Q27. The load value is 4.4 pu. The noise is added
to Pm3, the mechanical power of generator number 3. It is clear from the flgure that
IOPFs can be used to aid in the selection of signals to monitor the system other
than signals of state variables.
67
Parameters Machine 1 Machine 2 Machine 3
H(secs) 23.64 6.4 3.01
xd(pu) 0.146 0.8958 1.3125
x0d(pu) 0.0608 0.1198 0.1813
xq(pu) 0.0969 0.8645 1.2578
x0q(pu) 0.0969 0.1969 0.25
T0do(sec) 8.96 6.0 5.89
T0qo(sec) 0.31 0.535 0.6
Table 4.4: Three-generator nine-bus power system: machine data.
Parameters KA TA(sec) KE TE(sec) KF TF(sec)
20 0.2 1.0 0.314 0.063 0.35
Table 4.5: Three-generator nine-bus power system: exciter data.
States measured
Input noise Efd1 Efd2 Efd3 !1 !2 !3
added to
Pm1 3.0017 2.6973 2.1357 0.0033 0.0028 0.0031
Pm2 2.6113 2.3465 1.858 0.0029 0.0024 0.0027
Pm3 4.7816 4.2967 3.4022 0.0052 0.0044 0.0049
Vref1 0.0155 0.014 0.0111 0.0000169 0.0000143 0.000016
Vref2 0.0233 0.021 0.0166 0.0000255 0.0000215 0.00002409
Vref3 0.0475 0.0427 0.0338 0.0000519 0.000043 0.000049
Table 4.6: Input-to-state participation factors for the 3-machine nine-bus system
(partial listing). The load at bus 5 is 4.4 pu.
68
 
230 KV 
Gen 1 
230 KV 
Gen 3 
13.8 KV 
1.025 pu 
Gen 2 
18 KV 
1.025 pu 
90 MW 
30 MVAR 
125 MW 
50 MVAR 
Gen 1 
16.5 KV 
1.04 pu 
1 
4 
5 6 
8 9 7 
2 3 
85 MW 
Slack Bus 
100 MW 
35 MVAR 
230 KV 
163 MW 
Figure 4.8: WSCC 3-machine, 9-bus system.
Outputs measured
Input noise V5 V6 V7 Q14 Q39 Q27
added to
Pm1 0.8255 0.6490 0.6922 4.2022 0.9882 1.4460
Pm2 0.7182 0.5646 0.6021 3.6557 0.8597 1.2579
Pm3 1.3152 1.0339 1.1027 6.6948 1.5743 2.3037
Vref1 0.0043 0.0034 0.0036 0.0218 0.0051 0.0075
Vref2 0.0064 0.0050 0.0054 0.0327 0.0077 0.0112
Vref3 0.0131 0.0103 0.0110 0.0666 0.0157 0.0229
Table 4.7: Input-to-output participation factors for the 3-machine nine-bus system
(partial listing). The load at bus 5 is 4.4 pu.
69
0 0.5 1 1.5 2
10?6
10?5
10?4
10?3
10?2
10?1
Freq [Hz]
Power spectral density
Efd1Efd
2Efd
3
Figure 4.9: Power spectral densities of the states Efd1, Efd2 and Efd3. The load on
bus 5 was used as a bifurcation parameter. The load value is 4.0 pu. Band-limited
AWGN of zero mean and 0:05 power was added to Pm3, the mechanical power of
generator number 3.
70
0 0.5 1 1.5 2
10?5
10?4
10?3
10?2
10?1
Freq [Hz]
Power spectral density
Efd1Efd
2Efd
3
Figure 4.10: Power spectral densities of the states Efd1, Efd2 and Efd3. The load on
bus 5 was used as a bifurcation parameter. The load value is 4.4 pu. Band-limited
AWGN of zero mean and 0:05 power was added to Pm3, the mechanical power of
generator number 3.
71
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6
10?5
10?4
10?3
10?2
10?1
Freq [Hz]
Power spectral density of E
fd 1
PL
5
=4.0 pu
PL
5
=4.25 pu
PL
5
=4.4 pu
Figure 4.11: Power spectral density of Efd1 for three values of PL5 (the load on bus
5): PL5 = 4:0 pu (dash-dotted line), PL5 = 4:25 pu (dashed line) and PL5 = 4:4 pu
(solid line). Band-limited AWGN of zero mean and 0:05 power was added to Pm3,
the mechanical power of generator number 3.
72
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.610?6
10?5
10?4
10?3
10?2
Freq [Hz]
Power spectral density
Q27Q
14Efd
1
Figure 4.12: Power spectral densities of the state Efd1 and outputs Q14 and Q27.
The load on bus 5 was used as a bifurcation parameter. The load value is 4.4 pu.
Band-limited AWGN of zero mean and 0:05 power was added to Pm3, the mechanical
power of generator number 3.
73
Chapter 5
Participation Factor Assignment Using
Eigenvector Assignment-Based Feedback
Control
In this chapter, we develop state feedback control to assign modal participation fac-
tors for systems nearing instability using eigenvector assignment-based techniques.
We focus on assigning participation factors for state variables within a predeflned
group of state variables, where the group can be any subvector of the state vector,
or could consist of all system states. A procedure for computing the desired closed-
loop right eigenvector(s) (based on a given desired closed-loop participation factors)
is presented.
This chapter proceeds as follows. In Section 5.2, the algorithms for assign-
ing participation factors of system states in a mode associated with a single real
eigenvalue approaching zero are presented. In Section 5.3, the case of assigning
74
participation factors for systems with an open-loop mode associated with a pair of
conjugate complex eigenvalues approaching the imaginary axis is discussed.
5.1 Introduction
The goal of the control efiort is not to shift the eigenvalue(s) associated with a mode
that is approaching the imaginary axis (denoted as the critical mode). Rather, the
goal is to assign the participation factors of a predeflned group of states in the critical
mode. We would typically desire that some of these participation factors be small,
and possibly that others be relatively large so that the energy in the mode afiects
these physical states more than the \protected" ones. This could be viewed as using
some states, which could be associated with devices inserted into the system, as
vibration absorbers that tend to shield the most valued parts of the system from
the instability especially in its initial stages.
The feedback control will be designed to only change the closed-loop right
eigenvector(s) associated with the critical mode. The remaining closed-loop right
eigenvectors (associated with the remaining eigenvalues of the system) will be kept
the same as their open-loop values.
We address two cases, flrst the case of assigning participation factors of states
in a mode associated with a single real eigenvalue approaching zero, and second the
case of assigning participation factors of states in a mode associated with a pair of
conjugate complex eigenvalues approaching the imaginary axis.
75
Consider the closed-loop LTI system
_x(t) = (A+BF)x(t) (5.1)
associated with an open-loop system
_x(t) = Ax(t)+Bu(t) (5.2)
to which a state feedback u(t) = Fx(t) is applied.
In Chapter 4, it has been shown that noisy precursors can be used as a
warning signal indicating that a system is operating dangerously close to instability.
This warning signal can be used to trigger some control system as the one presented
in this chapter. The objective will be to move the system to a safer operating
condition in the sense that the \high-value" vulnerable parts of the system are more
\protected".
Minimizing the participation factors of a high-value group of states in the
critical mode implies that some other group of states will experience higher magni-
tudes of participation factors in the critical mode. Deciding which group of states
will be seen as \low-value" is based on system conditions and operators experience.
To avoid confusion, we will denote the group of high-value and low-value
states as the control group states. Also, we will arrange the system states such
that high-value states are followed by low-value states followed by the rest of the
system states. Based on the number of states in the control group and the number
of system inputs m, assigning participation factor associated with a critical mode
76
can be divided into two cases, (i) when nng ? m, and (ii) when nng > m, where
nng is the number of states in the control group. In the flrst case, we will show that
exact assignment of participation factors in the control group is possible. In the
second case, an optimal solution will be sought as the freedom of exact assignment
is limited by the number of system inputs (which is less than the number of states
in the control group).
5.2 Participation Factor Assignment: Case of
Real Eigenvalue Approaching Zero
In this section, we develop feedback control algorithms for assigning participation
factors of states in a mode associated with a single real eigenvalue approaching
zero. Such a mode is considered critical in the sense that it is associated with the
possibility of the system losing stability. Indeed, in such cases, the system would
typically undergo some type of stationary bifurcation at the critical parameter values
where the zero eigenvalue is actually attained. Next, we show how to compute the
desired closed-loop right eigenvector ri0d based on a given (desired) speciflcation of
the closed-loop participation factors. Note that the superscript (0) is used to denote
closed-loop quantities.
77
5.2.1 Computing desired closed-loop right eigenvector ri0d
The desired closed-loop i-th right eigenvector ri0d can be written as
ri0d :=
2
66
66
66
4
rs
rns
rnc
3
77
77
77
5
;
where, rs 2 Rns; rns 2 Rnns; rnc 2 Rn?(ns+nns) represent vectors of components
of ri0d corresponding to high-value, low-value and uncontrolled states, respectively.
fns; nns; nncgrepresentthenumberofhigh-value, low-valueanduncontrolledstates,
respectively. The high-value and low-value states form a control group. Let nng =
ns +nns be the number of control group states.
The case of zero closed-loop participation factors of high-value states in mode
i is straightforward and can be achieved by setting rs = 0 and rns = e 2Rnns where e
is a column vector of nonzero elements (to avoid zero participation). The remaining
(n?nng) components are left unspecifled (unconstrained).
Next, we consider the general case of nonzero closed-loop participation fac-
tors. Let the desired closed-loop participation factors of the control group states
fx2;:::;xnngg in mode i be given in terms of the closed-loop participation factor of
state x1 (flrst state of the group), i.e., let zk, k = 2;:::;nng be the desired ratio
between the closed-loop participation factor of state xk and that of state x1. Next,
we state a proposition that focuses on solving this participation factor assignment
problem in the case where we only desire to attain a certain ratio of closed-loop
78
participation factors for two states. The proposition shows how to compute the
ratio between components of the closed-loop right eigenvector ri0 (associated with
the real i-th eigenvalue ?i) to achieve a given (desired) ratio between the closed-loop
participation factors of two states fxk; xjg in mode i.
Proposition 5.1 Consider a linear time-invariant system (5.1{5.2) with the goal of
achieving a desired ratio z between the closed-loop participation factors of any two
states fxk; xjg in the i-th mode (viewed as the critical mode). It is possible to achieve
this desired ratio through state feedback without changing any right eigenvector except
for ri. Moreover, this is achieved by assigning the k-th and j-th components frik0; rij0g
of the closed-loop i-th right eigenvector such that their ratio is zli
k=l
i
j
, where lik; lij are
the k-th and j-th components of the open-loop i-th left eigenvector, respectively.
Proof: See Appendix A.1.
For the general case where more than two participation factors are to be
assigned, let zk, k = 2;:::;nng be the desired ratio between the closed-loop par-
ticipation factor of state xk and that of state x1. The flrst nng components of the
closed-loop i-th right eigenvector should be selected such that the ratio between the
k-th component and the flrst component is zkli
k=l
i1, for k = 2;:::;nng. The remaining
(n?nng) components are unspecifled (unconstrained).
79
5.2.2 Feedback control design: exactly assignable case
In this case, the number of states in the control group (formed of high-value and
low-value states) is less than or equal to the number of inputs m, i.e.,
nng = ns +nns ? m:
Recall the main assumption of the assignment problem, that is, the controller should
not change any open-loop eigenvalue or any open-loop right eigenvector except the
eigenvector associated with mode i. Saying that, observe that for all open-loop
pairs (?j; rj), fj = 1;2;:::;n; j 6= ig, conditions (1) and (2) of Proposition 2.1 are
trivially satisfled. Moreover, by deflnition, rj lies in the null space of (?jI?A), and
hence in the span of N?j (satisfying condition (3) of the Proposition). Therefore,
according to Preposition 2.1, there exists a feedback gain matrix F such that all
open-loop pairs (?j; rj), fj = 1;2;:::;n; j 6= ig can be kept invariant from their
original, uncontrolled values. Next we recall the steps to compute the state feedback
gain matrix F given in Section 2.3.1 with appropriate modiflcations to solve the
assignment problem beforehand.
Steps to compute F
1. Using the procedure given in Section 5.2.1, compute the desired vectors rs; rns
of components of the i-th closed-loop right eigenvector ri0 (to be exactly as-
signed) corresponding to high-level and low-level states, respectively.
80
2. Form a column vector e 2Rm?nng of nonzero elements.
3. Compute k to satisfy
2
66
66
66
4
rs
rns
e
3
77
77
77
5
= ~N?i k (5.3)
where ~N?i is the m?m matrix of the flrst m rows of N?i (deflned in (2.19))
associated with the flrst m states.1 The closed-loop right eigenvector ri0 is
then given by
ri0 = N?i k
4. Let
wi = ?M?i k;
where M?i is given by (2.19). Construct the real matrices W 2 Rm?n and
V 2Rn?n
W = [0 ::: 0 wi 0 ::: 0]
V = [:::Refrjg Imfrjg :::ri0:::rq :::]
where rj and rq represent open-loop right eigenvectors associated with complex
and real eigenvalues, respectively.
1In case ~N?
i is singular, a state from the control group has to be replaced by a state from
outside the group and test again for the singularity of the new ~N?i.
81
5. The feedback gain matrix F is given by
F = WV ?1:
Numerical example
Consider an example of a linear time-invariant system (5.2) ([23], pp. 108), where
A =
2
66
66
66
66
66
4
0 1 0 0
0 ?1:89 0:39 ?5:555
0 ?0:034 ?2:98 2:43
0:035 ?0:0011 ?0:99 ?0:21
3
77
77
77
77
77
5
and B =
2
66
66
66
66
66
4
1 0
0:76 ?1:6
?0:95 ?0:032
0:03 0
3
77
77
77
77
77
5
: (5.4)
The system?s four eigenvalues are f?0:1051; ?1:4811? 0:6239i; ?2:0127g. Since
the eigenvalue f?0:1051g is the closest to the origin, let us focus on its associated
participation factors for assignment through feedback. The open-loop participation
factors of system states fx1; x2; x3; x4g in the flrst mode are f1.1492, -0.0653,
0.0403, -0.1241g, respectively. The ratio between the participation factors offx1;x2g
is ?17:6, i.e., the participation of x1 in the flrst mode is much higher than the
participation of state x2. Suppose it is desired to reduce the ratio between the
participation factors of fx1;x2g to 0:5 (i.e., it is desired to push the participation of
the second state in the flrst mode to be twice that of the participation of the flrst
82
state). The right modal matrix R of the open-loop system is
R =
2
66
66
66
66
66
4
0:9934 ?0:4738+0:1996i ?0:4738?0:1996i 0:4447
?0:1044 0:8262 0:8262 ?0:8951
0:0314 ?0:1738+0:0881i ?0:1738?0:0881i ?0:0221
0:0357 ?0:0730+0:0990i ?0:0730?0:0990i ?0:0213
3
77
77
77
77
77
5
:
The open-loop left eigenvector components l11; l12 are f1.1568, 0.6258g, respectively.
Thus the design ratio for r110=r120 should be 0:5=(1:1568=:6258) = 0:2705. Applying
the steps above to compute the feedback gain matrix F yields
F =
2
66
4
0:5280 0:2857 0:5846 ?1:5851
0:1578 0:0854 0:1748 ?0:4738
3
77
5
The closed-loop system matrix Ac = A+BF is given by
Ac =
2
66
66
66
66
66
4
?0:5280 0:7143 ?0:5846 1:5851
?0:1488 ?1:9705 0:2253 ?5:10850
0:5067 0:2401 ?2:4190 0:90900
0:0192 ?0:0097 ?1:0075 ?0:1624
3
77
77
77
77
77
5
: (5.5)
The eigenvalues of Ac are the same as the open-loop eigenvalues (as expected), and
the closed-loop right modal matrix R0 is
R0 =
2
66
66
66
66
66
4
0:2457 ?0:4738+0:1996i ?0:4738?0:1996i 0:4447
0:9083 0:8262 0:8262 ?0:8951
0:0152 ?0:1738+0:0881i ?0:1738?0:0881i ?0:0221
?0:3382 ?0:0730+0:0990i ?0:0730?0:0990i ?0:0213
3
77
77
77
77
77
5
:
83
Theclosed-loopparticipationfactorsaref0.1389, 0.2778, 0.0095, 0.5738g. Notethat
the participation of x1 in the flrst mode is half that of x2 (as desired). Also note
that, although the desired goal was met, state x2 is not the state having the highest
participation in the flrst mode. State x4 is the one having the highest closed-loop
participation although it was not the one with the highest open-loop participation.
Thus, the exact assignment methodology can precisely assign participation factors
among the states in the control group but has no authority over the participation
factors for states outside the control group.
5.2.3 Feedback control design: not exactly assignable case
In this case,
nng = ns +nns > m:
Let the desired closed-loop participation of the control group states fx2;:::;xnngg
in mode i be given in terms of the closed-loop participation factor of state x1 (flrst
state of the group), i.e., let zk, k = 2;:::;nng be the desired ratio between the
closed-loop participation factor of state xk and that of state x1. Generally, the
freedom in specifying components of the closed-loop eigenvectors is limited by the
number of system inputs m as noted by Liu and Paton [23]. Thus, it is not possible
to exactly assign the desired ri0d (computed in Section 5.2.1) to achieve the specifled
ratios zk, k = 2;:::;nng. Instead, a best possible choice for an achievable eigenvector
is made. Recall the expression to compute the best achievable vector ri0a and the
84
corresponding coe?cient vector k (given in (2.26){(2.27)),
k = (NT?i N?i)?1 NT?i ri0d (5.6)
ri0a = N?i (NT?i N?i)?1 NT?i ri0d: (5.7)
Construction of the feedback gain matrix F is straightforward following steps (4)
and (5) in Section 5.2.2.
5.2.4 Feedback control design: \sacriflcial" case
In this case, no prior speciflcations on the closed-loop participation factors are given
except that it is required that the participation factors of the low-value states be
maximized with respect to the rest of the system states. Thus, the low-value states
will be \sacriflced" to \protect" the rest of the system states (including the high-
value states). A trivial assumption is
n?nns ? m;
otherwise, it would be possible to assign zeros to all participation factors of system
states except the low-value ones to achieve inflnite relative participation of low-
value states with respect to the rest of the system states. The proposed setting is
to maximize the ratio between the norms of two vectors. The entries of the two
vectors are those of the closed-loop participation factors of low-value states and the
rest of the system states, respectively. Next, we show the setup of the optimization
problem.
85
For simplicity of notation, the letter i used as a subscript with the critical
eigenvalue, or used as a superscript with the right/left eigenvector associated with
the critical mode, is omitted whenever it is clear from the context. The subscripts
ns and t are used to denote quantities that are a function of low-value and rest-of-
the-system states respectively. Also, the states are reordered such that low-value
states come flrst followed by the rest of the system states.
The closed-loop i-th right eigenvector can be written as:
r0 =
2
66
4
r0ns
r0t
3
77
5
=
2
66
4
N?ns
N?t
3
77
5k (5.8)
where, N?ns and N?t are the matrices formed of the rows of N? associated with
low-value and the rest of the system states, respectively. Deflne column vectors Pns
and Pt of closed-loop participation factors of low-value states and the rest of the
system states, respectively,
Pns := [p01 ::: p0ns]T
Pt := [p0ns+1 ::: p0n]T
where, the closed-loop participation factor of state q in the real eigenvalue i, p0qi =
liq0riq0, can be written in terms of open-loop quantities and riq0 as (see (A.2) in Ap-
pendix A)
p0qi = jRjjR0jliqriq0:
86
Deflne diagonal matrices Lns and Lt as
Lns =
2
66
66
66
66
66
4
li1
li2
...
lins
3
77
77
77
77
77
5
and Lt =
2
66
66
66
66
66
4
lins+1
lins+2
...
lin
3
77
77
77
77
77
5
;
respectively. Using (5.8), Pt and Pns can be written as
Pns = jRjjR0j Lns N?ns k; (5.9)
Pt = jRjjR0j Lt N?t k (5.10)
respectively. The solution to the optimization problem is to flnd the coe?cient
vector k that maximizes the ratio between the norms of the two vectors fPns; Ptg.
The statement of the problem is
max Q := Pns
TPns
PtTPt (5.11)
subject to kTk = 1
where (?)T indicates taking the transpose of a matrix or a vector. Using (5.9){(5.10),
the objective function Q in (5.11) can be written as
Q = k
T N?nsTLT
nsLnsN?ns k
kT N?tTLTt LtN?t k
Denote N?nsTLTnsLnsN?ns and N?tTLTt LtN?t as G and H respectively. The objective
function Q can be rewritten in terms of matrices G and H as
Q = k
T G k
kT H k
87
It is clear that G and H are symmetric matrices. If H is not positive deflnite, we
may select k to drive the denominator of this objective function to zero. This would
be the best case for our design; the participation factors of the low-value states in the
mode i will be inflnite relative to the participation factors of the rest of the system
states. More likely, however, is the case in which H is of full rank. In this case,
H will have a well-deflned matrix square root (H1=2 = LtN?t), and the maximum
value of the objective function Q in (5.11) is equal to the maximum eigenvalue of
E := H?1=2GH?1=2: (5.12)
The associated principal direction of Q, which we may denote as vmax, determines
k via
k = H?1=2vmax
Therefore, it is possible to analytically flnd a coe?cient vector k that solves the
optimization problem. Knowing the optimal coe?cient vector kopt, the closed-loop
right eigenvector ri0 can be computed using (5.8). Construction of the feedback gain
matrix F is straightforward following steps (4) and (5) in Section 5.2.2.
Numerical example
In this example, we revisit the example given in Section 5.2.2. Recall that the
system has an eigenvalue equal to f?0:1051g and it is the closest eigenvalue to the
origin. We focus on the case where states fx1; x2g are the high-value states while
88
states fx3; x4g are the low-value ones. No prior speciflcations on the closed-loop
participation factors are given except that it is required that the participation factors
of the low-value states be maximized with respect to the rest of the system states
(high-value states in this example). Since the number of high-value states is two
and the system has two inputs, one might think that a \good" solution is to exactly
assign zeros to the closed-loop participation factors of states fx1; x2g. This is not
possible as the rank of N?t (matrix formed of the rows of N? associated with the
high-value states) is two, i.e., no coe?cient vector k can be found to assign zeros to
the closed-loop participation factors of states fx1; x2g in the flrst mode. Thus, our
objective is to assign the closed-loop participation factors of the system states to
maximize the ratio of the participation factors of states fx3; x4g to states fx1; x2g.
Recall that open-loop participation factors of system states fx1; x2; x3; x4g
in mode one (?i = ?0:1051) aref1.1492, -0.0653, 0.0403, -0.1241g, respectively. The
ratio between the norms of the vectors of open-loop participation factors of states
fx3; x4g and fx1; x2g is 0:0129. Computing matrix E in (5.12) yields
E =
2
66
4
1:0437 1:8279
1:8279 3:2066
3
77
5
The eigenvalues of E are f0:0013; 4:249g, and therefore, the maximum attainable
value of the objective function Q is 4:249. The corresponding coe?cient vector kopt
that achieves that value is [0:5107 0:8598]T. Computing the feedback gain matrix
89
F yields
F =
2
66
4
0:5998 0:3245 0:6641 ?1:8007
0:1793 0:0970 0:1985 ?0:5383
3
77
5
The closed-loop participation factors of system states fx1; x2; x3; x4g in mode one
are f0.0015, 0.3244, 0.0053, 0.6687g, respectively. The ratio between the norms of
the vectors of closed-loop participation factors of states fx3; x4g and fx1; x2g is
4:249. This result shows a f4:249=0:0129g ? 329{fold increase between open-loop
and closed-loop cases.
5.3 Participation Factor Assignment: Case of
Eigenvalue Pair Approaching Imaginary Axis
In this section, we develop a state feedback control for assigning participation factors
for systems nearing another form of instability. The approaching instability studied
here is Hopf bifurcation, where a pair of complex conjugate eigenvalues cross the
imaginary axis. Let f?i, ?i+1 = ??ig be the pair of complex eigenvalues approaching
the imaginary axis. As in the case of real eigenvalue approaching imaginary axis
given in Section 5.2, the goal pursued by the control efiort is to achieve a prescribed
set of ratios between the norms of closed-loop participation factors of states in
a group denoted as the control group. This will be done without afiecting any
eigenvalues or any right eigenvectors except those right eigenvectors associated with
90
mode i. Next, we show how to compute the desired closed-loop right eigenvectors
fri0d; ri+10dg based on a given (desired) speciflcation of the closed-loop participation
factors.
5.3.1 Computing desired closed-loop right eigenvectors
fri0d; ri+10d = ri0d?g
The same notation used in Section 5.2.1 to specify difierent components of the
desired closed-loop right eigenvector ri0d will be adopted here, i.e.,
ri0d :=
2
66
66
66
4
rs
rns
rnc
3
77
77
77
5
;
where, rs 2 Cns; rns 2 Cnns; rnc 2 Cn?(ns+nns) represent vectors of components
of ri0d corresponding to high-value, low-value and uncontrolled states, respectively.
fns; nns; nncgrepresentthenumberofhigh-value, low-valueanduncontrolledstates,
respectively. The high-value and low-value states form a control group. Let nng =
ns +nns be the number of control group states.
The case of zero closed-loop participation factors of high-value states in mode
iisstraightforwardanditcanbeachievedbysettingrs = 0andrns = e 2Cnns where
e is a column vector of elements with nonzero norms (to avoid zero participation).
The remaining (n?nng) components are unspecifled (unconstrained).
For the general case of nonzero closed-loop participation factors, let the
91
desired norms of the closed-loop participation factors of the control group states
fx2;:::;xnngg in mode i be given in terms of the norm of the closed-loop partici-
pation factor of state x1 (flrst state of the group), i.e., let zk, k = 2;:::;nng be
the desired ratio between the norm of the closed-loop participation factor of state
xk and that of state x1. Next, we show how to write the closed-loop participation
factor of the j-th state xj in mode i in terms of the closed-loop right eigenvector ri0.
Recall that
p0ji = lij0rij0 (5.13)
Denote by R0 the matrix of closed-loop right eigenvectors
R0 =
2
66
66
66
66
66
66
66
66
4
r11 ::: ri10 ri10? ::: rn1
... ... ... ... ... ...
r1j ::: rij0 rij0? ::: rnj
... ... ... ... ... ...
r1n ::: rin0|{z}
ri0
rin0?|{z}
ri+10
::: rnn
3
77
77
77
77
77
77
77
77
5
where ri0 and ri+10 = ri0? are the only closed-loop eigenvectors that are difierent
from their corresponding open-loop values. The j-th component of the i-th closed-
loop left eigenvector lij0, which lies in the i-th row and j-th column of the matrix of
closed-loop left eigenvectors L0 (L0 = R0?1), is given by:
lij0 = (?1)
j+i
detfR0g detfMijg; (5.14)
92
where Mij (Mij 2 C(n?1)?(n?1)) is the matrix that results from removing the j-th
row and the i-th column of R0, i.e.,
Mij =
2
66
66
66
66
66
66
66
66
66
66
4
r11 ::: ri?11 ri10? ri+21 ::: rn1
... ... ... ... ... ... ...
r1j?1 ::: ri?1j?1 rij?10? ri+2j?1 ::: rnj?1
r1j+1 ::: ri?1j+1 rij+10? ri+2j+1 ::: rnj+1
... ... ... ... ... ... ...
r1n ::: ri?1n rin0?|{z}
ri+10
ri+2n ::: rnn
3
77
77
77
77
77
77
77
77
77
77
5
Let
Qj := [(?1)i+1jCi1j:::(?1)i+j?1jCi(j?1)j 0|{z}
j-th
(?1)i+jjCijj:::(?1)i+n?1jCi(n?1)j]T
where matrices Cim (Cim 2 C(n?2)?(n?2)) are formed by removing the m-th row
and i-th column of Mij, (m = 1;:::; n ? 1). The determinant of Mij can now
be obtained by cofactor expansion along the i-th column, yielding an expression in
terms of ri+10 = ri0? and Qj:
detfMijg = ri0HQj (5.15)
where (?)H indicates taking the Hermitian, or complex conjugate transpose, of a
matrix or vector. Using (5.13){(5.15), the closed-loop participation factor p0ji in
(5.13) can be written as
p0ji = f(?1)
j+i
detfR0gg?fr
i0HQjg?ri
j
0 (5.16)
93
Recall from (2.22) that ri0 can be written as
ri0 = N? k
for some vector k 2 Cm and N? deflned as given in (2.18). Using that, p0ji can be
expressed as
p0ji = f(?1)
j+i
detfR0gg?fk
HN?H Qjg?N?
j k (5.17)
where N?j is the j-th row of N?. To simplify the expression in (5.17), denote by ?j
the m?m matrix
?j := N?H Qj N?j: (5.18)
Equation (5.17) now becomes
p0ji = (?1)
j+i
detfR0g ?k
H?j k (5.19)
Note that square matrix ?j (?j 2C(m?m)) is singular and of rank one.
Next, we state a proposition that focuses on solving this participation fac-
tor assignment problem in the case where we only desire to attain a certain ratio
of closed-loop participation factors for two states. The proposition shows how to
compute the i-th closed-loop right eigenvector ri0 (associated with the complex i-th
eigenvalue ?i) to achieve a given (desired) ratio between the norms of closed-loop
participation factors of two states fxj; xlg in mode i.
Proposition 5.2 Consider a linear time-invariant system (5.1{5.2), with the goal
of achieving a desired ratio z between the norms of the closed-loop participation
94
factors of any two states fxj; xlg in the critical mode \i." Suppose that the matrices
f?j; ?lg don?t have the same null space. (Recall that matrices f?j; ?lg are given
by (5.18).) Then, it is possible to achieve the desired ratio z of participation factors
without afiecting any right eigenvectors except those associated with mode i.
Proof: See Appendix A.2.
Thus, accordingtoProposition5.2, itislikelytoflndvectorskj, j = 2;:::;nng
that solve for all desired ratios zj, j = 2;:::;nng. Despite that, it might be impos-
sible to flnd a single vector k that simultaneously solves for all desired ratios zj,
j = 2;:::;nng. A candidate coe?cient vector kcad could be a weighted sum of the
computed kj?s that solves for the zj?s. Those weights can be chosen based on the
most desired ratios. If for example, the ratios are equally desired, then the candidate
vector kcad is simply the average of the computed kj?s.
5.3.2 Feedback control design: exactly assignable case
Unlike exactly assignable cases discussed in Section 5.2.2 for systems with a real
eigenvalue approaching zero, the cases discussed in this section are limited. Two
cases are discussed, case (a) where we only desire zero closed-loop participation
of high-value states (of dimension less than number of inputs m) in mode i, and
case (b) where we only desire to attain certain ratio between norms of closed-loop
participation factors for two states only. The feedback control design follows a
similar approach to the one given for the real eigenvalue case except for some minor
95
changes.
Recall the main assumption of the assignment problem, that is, the controller
should not change any open-loop eigenvalue or any open-loop right eigenvector ex-
cept the eigenvectors (ri; ri+1 = ri?) associated with mode i. Saying that, the same
argument presented in the real eigenvalue case holds, i.e., according to Proposition
2.1, there exist a feedback gain matrix F such that all open-loop pairs (?j; rj),
fj = 1;2;:::;n; j 6= fi;i + 1gg can be kept invariant from their original, uncon-
trolled values. Next we recall the steps to compute the state feedback gain matrix
F given in Section 2.3.1 with the appropriate modiflcations to solve the assignment
problem beforehand.
Steps to compute F: case (a)
1. Set rs (vector of closed-loop right eigenvector components corresponding to
high-value states) to zeros.
2. Form a column vector e 2Cm?ns of elements of nonzero norms.
3. Compute k to satisfy
2
66
4
rs
e
3
77
5 = ~N?i k
where ~N?i is the m?m matrix of the flrst m rows of N?i (deflned in (2.18))
associated with the flrst m states.2 The closed-loop right eigenvector ri0 is
2In case ~N?
i is singular, see footnote in page 81.
96
then given by
ri0 = N?i k
4. Let
wi = ?M?i k:
Construct the real matrices W 2Rm?n and V 2Rn?n
W = [0 ::: 0 Refwig Imfwig 0 ::: 0]
V = [:::Refrjg Imfrjg :::Refri0g Imfri0g:::rq :::]
where rj and rq represent open-loop right eigenvectors associated with complex
and real eigenvalues, respectively.
5. The feedback gain matrix F is given by
F = WV ?1:
Steps to compute F: case (b)
1. Using the result stated in Proposition 5.2, compute the coe?cient vector k
that achieves the desired ratio of participation factors between two states.
The closed-loop right eigenvector ri0 is then given by
ri0 = N?i k
2. Follow steps (4) and (5) of the procedure to compute F given for case (a)
above.
97
5.3.3 Feedback control design: not exactly assignable case
This case involves the assignment of closed-loop participation factors of more than
two states. Assuming the existence of solution vectors kj, j = 2;:::;nng that
solve for all desired ratios zj, j = 2;:::;nng, a candidate coe?cient vector kcad
is a weighted sum of the computed kj?s that solves for the zj?s (as discussed in
Section 5.3.1). The weights are chosen based on system?s operation and operator?s
experience. The closed-loop right eigenvector ri0 is given by
ri0 = N?i kcad:
Construction of the feedback gain matrix F is straightforward following steps (4) and
(5) of the procedure to compute F given for case (a) of the exactly assignable case.
Next, we present a numerical example of participation factor assignment involving
more than two states.
Numerical example
The example in Section 5.2.2 is revisited but with a slight modiflcation. The system
considered here has an extra input, i.e., the input matrix B has an extra column,
B =
2
66
66
66
66
66
4
1 0 2
0:76 ?1:6 ?5
?0:95 ?0:032 3
0:03 0 2:8
3
77
77
77
77
77
5
:
98
As noted in the example in Section 5.2.2, The system?s four eigenvalues are
f?0:1051;?1:4811 ? 0:6239i;?2:0127g. For the mode associated with the pair of
complex conjugate eigenvalues (eigenvalues 2 and 3), let us focus on its associ-
ated participation factors for assignment through feedback. The open-loop par-
ticipation factors of system states fx1; x2; x3; x4g in that mode are f?0:0446 ?
0:1059i; 0:0129 + 0:1872i; 0:3868 ? 1:1663i; 0:6449 + 1:0850ig, respectively. The
norms of the open-loop participation factors are f0:1149; 0:1877; 1:2287; 1:2621g.
Suppose that three states fx1; x2; x4g are selected to form the control group with
state x4 chosen as the high-value state and states fx1; x2g are the low-value ones.
The open-loop participation ratios are jp12jjp42j = 0:0910 and jp22jjp42j = 0:1487. Note that
the open-loop participation of the high-value state x4 is higher than the participation
of the low-value states fx1; x2; x4g in the second (and third) mode.
Suppose it is desired to assign the participation of states fx1; x2; x4g to
achieve the following desired ratios between norms of closed-loop participations:
zd1 = jp012jjp0
42j
= 1:2 and zd2 = jp022jjp0
42j
= 2:24. Using the procedure given in the proof of
Proposition 5.2 (see Appendix A), the coe?cient vectors fk1; k2g that separately
solve for the desired ratios fzd1; zd2g, respectively, are
k1 = [0:3299 ?0:1856+0:3824i 0:8452+0:3967i]T;
k2 = [0:3531 ?0:1558+0:3824i 0:7623+0:3967i]T
Choosing kcad to be the average of fk1; k2g yields,
kcad = [0:3415 ?0:1707+0:3824i 0:8037+0:3967i]T
99
Table 5.1 shows the difierent possibly achievable closed-loop ratios fz1; z2g for
difierent coe?cient vectors fk1; k2g and for the candidate vector kcad.
z1 = jp012jjp0
42j
z2 = jp022jjp0
42j
k1 1.20 1.84
k2 1.45 2.24
kcad 1.31 2.02
Table 5.1: Closed-loop ratios fz1; z2g for coe?cient vectors fk1; k2; kcadg. The
desired closed-loop ratios are fzd1 = 1:2; zd2 = 2:24g.
100
Chapter 6
Application to Electric Power Networks
In this chapter, we apply the results obtained in Chapter 5 on participation factor
assignment for LTI systems to an example from electric power networks. The system
(at equilibrium) is operating close to instability where a pair of complex eigenvalues
is close to the imaginary axis. Modal analysis of linearized model of that network
is performed. Generators highly afiected by the critical mode are identifled. State
feedback control is applied to assign the participation factors of a group of the state
variables to protect high-value generating units.
This chapter proceeds as follows. In Section 6.2, Western System Coordinat-
ing Council (WSCC) 3{generator, 9{bus power system model, is examined. This
system has a critical electrical mode associated with the system exciters. A state
feedback control is designed for assigning closed-loop participation factors. The goal
fulfllled is to trade the risk of participation in the critical mode between a high-value
generator and the other generators in the system denoted as low-value generators.
101
In Section 6.3, remarks on the use of state feedback in power networks are given.
6.1 Introduction
For small disturbance analysis in electric power networks, the participation factor
of a mode in a state can be interpreted, from the view point of probability, as the
participation \risk". Larger norm of participation factor of a state in a particular
mode means it has more efiect on and will be more afiected by this mode. It may
be desirable that the participation factor of a particular state variable be reduced
in some critical system modes. This desire is based on the fact that the higher the
magnitude of participation factor the higher the tendency of such state to be afiected
in a serious way by an impending (or true) instability. In power systems operation,
nuclear units, for example, are regarded as high value subsystems. Reducing par-
ticipation factors of these units in some critical mechanical and/or electrical mode
would reduce the risk of premature tripping, as well as alleviate the efiect on these
units of a troublesome mode being excited.
6.2 Three-Generator Nine-Bus Power System
In this section, we revisit the Western System Coordinating Council (WSCC) 3-
generator, 9-bus power system model, studied for monitoring impending instabilities
in Section 4.2.3 (see Figure 4.8 for the line diagram). In this model, a subcritical
102
Hopf bifurcation occurs as the load on bus 5 is increased beyond 4.5 pu [33]. It
is observed that the critical modes associated with the pair of complex conjugate
eigenvalues closest to the imaginary axis are the electrical ones associated with the
exciter. The system is composed of three synchronous generators equipped with
three identical IEEE-type I exciters. A fourth order model is used to capture the
dynamics of the synchronous generator:
_? = ?b(! ?1)
_! = (Tm ?Te ?D(! ?1))=M
_e0q = (?fs(e0q)?(xd ?x0d)id +v?f)=T0d0
_e0d = (?e0d +(xq ?x0q)iq)=T0q0;
where the subscripts d and q refer to the direct and quadrature axes, respectively.
fs(:) captures the saturation characteristic of the synchronous generator. The elec-
trical power is
Pe = (vq +raiq)iq +(vd +raid)id
and current link is described by the equations:
0 = vq +raiq ?e0q +(x0d ?x1)id
0 = vd +raid ?e0d +(x0q ?x1)iq
103
The exciter is depicted in Figure 6.1 and described by the following equations:
_vm = (V ?vm)=T
_vr1 = (Ka(Vref ?vm ?vr2 ? KfT
f
vf)?vr1)=Ta (6.1)
vr =
8
>>>
>>><
>>>>
>>:
vr1 if vrmin ? vr1 ? vrmax ? vr1
vrmax ifvr1 > vrmax;
vrmin ifvr1 < vrmin:
_vr2 = ?(KfT
f
vf +vr2)=Tf
_vf = ?(vf(1+Se(vf))?vr)=Te
The generators and exciters data are given in Section 4.2.3. Details of loads and
 
Figure 6.1: Exciter for the three-generator, nine-bus system.
transmission lines can be found in [33]. The system has 24 states arranged as follows:
x(t) = [?1 !1 e0q1 e0d1 vm1 vr11 vr21 vf1 :::]T:
104
When the load at bus 5 is 4.4 pu, the linearization of the system at the operating
pointhasacomplex conjugatepair ofeigenvalues?fi;i+1g = ?0:04602?j2:1151. This
isthepairofcomplexconjugateeigenvaluesclosesttotheimaginaryaxis. Denotethe
mode associated with that pair of eigenvalues as the critical mode. Computing the
participation factors of system states in that mode reveals that states representing q-
axis transient voltage, e0q, are the states with the highest participation in the critical
mode. Table 6.1, shows the states having the highest open-loop participation factors
in the critical mode. Note that, the norms of the participation factors were divided
by the sum of the norms of the participation factors to get more sense out of the
numbers for comparison.
State jriqj \riq jpqij = jliqjjriqjx50n
k=1jlikjjrikj
e0q3 0.0108 66.802 0.19887
e0q1 0.0189 34.437 0.14365
e0q2 0.0230 43.401 0.09643
e0d1 0.0339 9.702 0.07081
e0d2 0.0323 5.327 0.05138
vr21 0.0481 -59.185 0.05207
vr22 0.0433 -55.0678 0.07031
Table 6.1: The states with the highest norms of open-loop participation factors in
the critical mode associated with ?fi;i+1g = ?0:04602?j2:1151. The load at bus 5
is 4.4 pu.
Observe that e0q3, q-axis transient voltage of the third generator, is the state
with the highest participation in the critical mode indicating that the third gen-
105
erator is the most afiected generator by the critical mode. Assume that the third
generator is a high-value generator (nuclear unit for example). We would like to
lower the participation of this generator in the critical mode to alleviate any unde-
sirable efiects that may possibly occur if the critical mode is excited. The critical
mode of the system can be revealed from the time domain simulation by applying a
2% perturbation to the initial condition of the q-axis transient voltage of generator
2 (e00q2) while keeping the rest of the states at their equilibrium values. The resulting
q-axis transient voltage of generators 1; 2, and 3 are shown in Figure 6.2. From
Figure 6.2, it is obvious that after short-time transients of the three generators, a
poorly damped mode with frequency around 0:33 Hz can be detected (note the sus-
tained oscillations). This is consistent with that the critical mode is associated with
a pair of eigenvalues with frequency of 0:3337 Hz. Note also how the 3 states reacts
difierently to the perturbation. e0q3, for example, is the most afiected followed by e0q2
followed by e0q1. This is expected from the values of the norms of the participation
factors of the q-axis transient voltage states in the critical mode (see Table 6.1).
The open-loop ratios between the norms of the participation factors of states
fe0q1, e0q2g and state e0q3 are f0:7223; 0:4849g, respectively. Suppose it is desired
that the closed-loop ratios between the norms of the participation factors of states
fe0q1, e0q2g and state e0q3 are fz13 = 2:92; z23 = 2:03g, respectively. Particularly, we
need to trade the risk of participation in the critical mode between the third genera-
tor (high-value generator) and the other two generators (low-value generators) of the
106
system. Next, we design a state feedback control to assign the system participation
factors in the critical mode to fulfll the goal mentioned above.
0 10 20 30 40 50 60 70
1.142
1.143
1.144
1.145(a)
e q 3?
0 10 20 30 40 50 60 70
0.86
0.865
0.87
0.875
(b)
e q 2?
0 10 20 30 40 50 60 700.92
0.922
0.924
0.926
Time (sec)
e q 1?
(c)
Figure 6.2: Open-loop response for a perturbation of 2% in the initial condition of
the q-axis transient voltage of generator 2 (e00q2) of the WSCC 3-gen., 9-bus power
system. (a) e0q3, (b) e0q2 and (c) e0q1.
6.2.1 Feedback control design
Recall the procedure for the control design given in Section 5.3.3. We have two
desired ratios fz13; z23g relating the closed-loop norms of the participation factors
of q-axis transient voltage states, e0qk; k = 1; 2; 3, to achieve. The control action
will take place through the the reference input of the exciter connected to each
generator. The input matrix B is a matrix of dimension (24?3) of zeros except at
107
entries (14;1); (18;2) and (22;3) corresponding to input vref of generators 1; 2 and
3, respectively (see (6.1)). The nonzero entries are identical and equal to KaTa = 100.
Using the procedure given in the proof of Proposition 5.2 (see Appendix A), the
coe?cient vectors fk1; k2g that separately solve for the desired ratios fz13; z23g,
respectively, are
k1 = [0:2753+0:1238i 0:5075?0:4070i 0:1508?0:5426i]T;
k2 = [0:2347?0:0288i 0:6427?0:1942i ?0:0374?0:5528i]T:
Choosing kcad to be the average of fk1; k2g yields,
kcad = [0:2550+0:0475i 0:5751?0:3006i 0:0567?0:5477i]T:
Table 6.2 shows the difierent possibly achievable closed-loop ratios for difierent co-
e?cient vectors fk1; k2g and for the candidate vector kcad.
z13 z23
k1 2.92 1.227
k2 1.9425 2.03
kcad 2.4314 1.6341
Table 6.2: Closed-loop ratios fz13; z23g for coe?cient vectors fk1; k2; kcadg. The
desired closed-loop ratios are fz13 = 2:92; z23 = 2:03g.
Table 6.3 shows the states having the highest closed-loop participation factors
in the critical mode. As expected, the participation of the q-axis transient voltage
state of generator 3 in the critical mode is reduced, and it is not the state with the
108
State jpqij = jliqjjriqjx50n
k=1jlikjjrikj
e0q3 0.0864
e0q1 0.20754
e0q2 0.13934
e0d1 0.07036
e0d2 0.04882
vr12 0.05019
vr21 0.07605
Table 6.3: The states with the highest norms of closed-loop participation factors in
the critical mode associated with ?fi;i+1g = ?0:04602?j2:1151. The load at bus 5
is 4.4 pu.
highest participation in the critical mode. Time domain simulations of this model
were conducted using PSAT [16]. Figures 6.3-6.5 show the difierent time response
of the q-axis transient voltage states, fe0q1; e0q2; e0q3g, for a 2% perturbation to the
initial condition of the q-axis transient voltage of generator 2 (e00q2) before and after
control. The Figures show the improvement in the response of the high-value state
e0q3 of the third machine and the efiects of increased participation in the critical
mode of the low-value states fe0q1; e0q2g of the flrst and second machine, respectively.
109
10 20 30 40 50 60 70
0.918
0.919
0.92
0.921
0.922
0.923
0.924
0.925
0.926
0.927
Time (sec)
e q 1?
Before controlAfter control
Figure 6.3: Plot of q-axis transient voltage of generator 1 (e0q1) versus time for a
perturbation of 2% in the initial condition of the q-axis transient voltage of generator
2 (e00q2) before and after control.
110
10 20 30 40 50 60 70
0.855
0.86
0.865
0.87
0.875
Time (sec)
e q 2?
After controlBefore control
Figure 6.4: Plot of q-axis transient voltage of generator 2 (e0q2) versus time for a
perturbation of 2% in the initial condition of the q-axis transient voltage of generator
2 (e00q2) before and after control.
111
10 20 30 40 50 60 701.141
1.1415
1.142
1.1425
1.143
1.1435
1.144
1.1445
1.145
1.1455
Time (sec)
e q 3?
Before controlAfter control
Figure 6.5: Plot of q-axis transient voltage of generator 3 (e0q3) versus time for a
perturbation of 2% in the initial condition of the q-axis transient voltage of generator
2 (e00q2) before and after control.
6.3 Remarks on Use of State Feedback in Power
Networks
The work in this chapter assumes access to the system state vector for purposes
of feedback control. There are several ways to deal with this assumption from a
practical perspective. One approach is to use a dynamic observer to calculate the
state vector from available output measurements. Another approach is to calculate
a state feedback control and then approximate it by an output feedback controller.
The third approach is to actually make real-time state measurements available, and
this is a real possibility in light of recent technological developments in real-time
112
phasor measurements [45]. As noted by Phadke [45], synchronized phasor measure-
ments provide the flrst real possibility of providing a dynamic state estimator. Using
such measurements, it is possible for the substations to provide a continuous stream
of phasor data to the control center. In this way, a state vector estimate that follows
the system dynamics can be constructed. With normal dedicated communication
circuits, a continuous data stream of one phasor measurement every 2-5 cycles (33.3
- 83.33 msec) can be sustained. Considering that typical power system dynamic
phenomena fall in the frequency range of 0-2 Hz, it is possible to observe power
system dynamic phenomena in real-time with high fldelity at the control center.
113
Chapter 7
Conclusions and Future Directions
7.1 Conclusions
In this dissertation, we have developed monitoring and control systems to improve
the performance of systems that are sometimes required to operate at the edge of
their stability envelopes. Our main contributions are as follows:
? We used the concept of modal participation factors extensively in this work,
and generalized the concept in several directions. In particular, we introduced
output participation factors and input-to-output participation factors.
? A signal-based approach for real-time detection of impending instability in
nonlinear systems was considered. The main idea pursued involves using a
small additive white Gaussian noise as a probe signal and monitoring the spec-
tral density of one or more measured states or outputs for certain signatures
of impending instability. Input-to-state and input-to-output participation fac-
114
tors were introduced as tools to aid in selection of locations for probe inputs
and states or outputs to be monitored, respectively. Since these participation
factors are model-based, the work presented combines signal-based and model-
based ideas toward achieving a robust methodology for instability monitoring.
We also introduced a proximity index to help assess closeness to instability
based on power spectral density measurement.
? We also introduced the problem of participation factor assignment by feed-
back, and made important progress on this problem. Feedback algorithms
were developed for assigning modal participation factors using eigenvector
assignment-based techniques. The goal is to reduce the interaction between a
selected group of states (high-value group) and an undesirable mode. In partic-
ular, we addressed two cases, flrst the case of a mode associated with a single
real eigenvalue approaching zero (Stationary bifurcation), and second the case
of a mode associated with a pair of conjugate complex eigenvalues approach-
ing the imaginary axis (Hopf bifurcation). A novel procedure for computing
the desired closed-loop right eigenvector(s) associated with the critical mode
(based on given constraints on the desired closed-loop participation factors)
was presented
? We applied our results to a well-studied electric power network, the WSCC
3-generator, 9-bus network. As expected, participation factor assignment has
a clear efiect on the relative time domain responses of system states.
115
7.2 Future Directions
Among the many open problems of interest are the following:
? More work is needed on generalized participation factors and on participation
factor assignment by state and output feedback. Particularly, the problem
of flnding a deflnition for participation of outputs in modes that reduces to
the original state participation when the output is simply one of the states
is yet to be solved. Unlike state participation problem that has identical
expressions for participation of modes in states and states in modes, the output
participation problem is quite difierent and an expression for participation of
outputs in modes is not expected to be the same as the expression proposed
for participation of modes in outputs in Chapter 3.
? A particularly challenging problem involves detection not only nearness to in-
stability, but also detecting the severity of the impending instability from the
point of view of nonlinear system behavior. For example, an oscillatory insta-
bility can be of the hunting type, in which small amplitude oscillations occur,
or it can be divergent, resulting in complete loss of operation. Although this
can be determined using analytical models using known methods of bifurca-
tion analysis, it is not known how this can be achieved using a signal-based
approach.
? Another direction involves studying use of other probe signals besides AWGN.
116
Examples include periodic signals, chaotic signals covering an appropriate fre-
quency range, and colored noise signals. The relative advantages and disad-
vantages of the various probe signals should be considered. In this regard,
connections to past work in real-time probing of power systems and aircraft
dynamics should be studied. In research aircraft, for example, it is common to
use \chirp" signals to probe the aircraft for its stability properties in various
parts of its  ight envelope.
117
Appendix A
Proof of Propositions 5.1 and 5.2
A.1 Proof of Proposition 5.1
The matrix of closed-loop right eigenvectors R0 is given by
R0 = [r1 ::: ri0 ::: rn];
where ri0 (associated with the real eigenvalue ?i) is the only closed-loop right eigen-
vector afiected by the control law. The matrix of closed-loop left eigenvectors
L0 = R0?1 is given by
L0 = [l10T ::: li0T ::: ln0T]T:
The closed-loop i-th left eigenvector li0 is given by
li0 = 1detfR0g[ci1 ::: cik ::: cin];
where the cik?s, k = 1;2;:::;n, represent matrix cofactors, i.e.,
cik = (?1)k+i ?detfMikg;
118
where Mik is the matrix of dimension (n?1) that is obtained by eliminating the k-th
row and the i-th column of R0. Since the only change to the open-loop right modal
matrix R is to the components of its i-th column ri, the cik?s, k = 1;2;:::;n are only
function in open-loop left eigenvectors lq?s, q = 1;:::; n; q 6= i. Therefore, the k-th
component of the closed-loop i-th left eigenvector li0 can be written in terms of the
k-th component of the open-loop i-th left eigenvector lik as
lik0 = detfRgdetfR0glik: (A.1)
Using (A.1), the closed-loop participation factor of state k in mode i is reduced to
p0ki = lik0rik0
= detfRgdetfR0glikrik0 (A.2)
Next, we use (A.2) to write the desired ratio z between the closed-loop participation
factors of states fxk;xjg in mode i as
z = p
0
ki
p0ji =
lik0rik0
lij0rij0 (A.3)
= l
i
kr
i
k
0
lijrij0; (A.4)
Therefore, to achieve the desired ratio z, the closed-loop right eigenvector compo-
nents rik0; rij0 should be chosen such that rik
0
rij0 =
z
lik=lij. ?
119
A.2 Proof of Proposition 5.2
The closed-loop participation factor of the j-th state xj in mode i in terms of the
closed-loop right eigenvector ri0 is given by (see (5.13)):
p0ji = (?1)
j+i
detfR0g ?k
H?j k (A.5)
The ratio between the norms of the closed-loop participation factors of two states
fxj; xlg in the critical mode \i" is:
jp0jij
jp0lij =
jlij0jjrij0j
jlil0jjril0j (A.6)
= jk
H?j kj
jkH?l kj; (A.7)
where matrices f?j; ?lg are given by (5.18). Matrices f?j; ?lg (2 Cm?m) are
singular and of rank one. The goal is to flnd a coe?cient vector k (2 Cm) that
achieves a desired ratio z between the norms of the closed-loop participation factors
of the two states fxj; xlg in the critical mode \i," i,e., flnd k that solves:
jp0jij
jp0lij =
jkH?j kj
jkH?l kj (A.8)
= z:
In general, if there are any solutions other than 0, one would expect the
solutions to form a variety of codimension 1 inCm. If we can flnd some k0 for which
the right side of (A.8) is less than z and some k1 for which it is greater than z, then
there is a solution k,
k = t k0 +(1?t) k1; 0 < t < 1: (A.9)
120
The squares of the norms in (A.8) are quartic polynomials in t. Solving
these quartic polynomials for t, we can flnd k that solves for the desired ratio. One
possible choice of k0 is a vector that lies in the null space of ?j and not in the
null space of ?l (in this case, the right side of (A.8) is equal to zero). One possible
choice of k1 is a vector that lies in the null space of ?l and not in the null space
of ?j (in this case, the right side of (A.8) is equal to inflnity). Note that this is
always possible providing that the assumption on matrices f?j; ?lg of Proposition
5.2 holds. Knowing t, the coe?cient vector k that solves for the desired ratio z can
be computed using (A.9). ?
121
Appendix B
Table of Symbols
Symbol Meaning
detfAg or jAj determinant of A
Refvg Real part of v
Imfvg Imaginary part of v
AT; xT transpose of matrix A, vector x
A?; x? conjugate of matrix A, vector x
AH; xH conjugate transpose of matrix A, vector x
jaj absolute value (magnitude, if a is complex) of a
kxk norm of vector x
ej j-th unit vector (j-th entry is 1; 0 otherwise)
Rn(Cn) real (complex) n-dim. space
Rm?n(Cm?n) real (complex) m?n matrices
122
Bibliography
[1] F. L. Pagola, I. J. Perez-Arriaga and G. C. Verghese, \On sensitivities, residues
and participations: Applications to oscillatory stability analysis and control,"
IEEE Transactions on Power Systems, Vol. 4, No. 1, pp. 278{285, Feb. 1989.
[2] I. J. Perez-Arriaga, G. C. Verghese and F. C. Schweppe, \Selective modal analy-
sis with applications to electric power systems, Part I: Heuristic introduction,"
IEEE Transactions on Power Apparatus and Systems, Vol. 101, No. 9, pp.
3117{3125, Sep. 1982.
[3] P. Kundur, Power System Stability and Control. McGraw Hill, New York, 1994.
[4] C. L. DeMarco, \Design of predatory generation control in electric power sys-
tems," Thirty-First Annual Hawaii International Conference on System Sci-
ences HICSS, Vol. 3, pp. 32{38, Jan. 1998.
[5] E. Abed, A. Tesi and H. Wang, \Control of bifurcation and chaos," in The
Control Handbook, (CRC Press), 1996.
123
[6] E. H. Abed and J. H. FU, \Local feedback stabilization and bifurcation control,
I. Hopf bifurcation," Systems and Control Letters, Vol. 7, pp. 11{17, February
1986.
[7] E. H. Abed and J. H. FU, \Local feedback stabilization and bifurcation control,
II. stationary bifurcation," Systems and Control Letters, Vol. 8, pp. 467{473,
May 1987.
[8] A.N. Andry, E.Y. Shapiro and J.C. Chung, \Eigenstructure Assignment for
Linear Systems,"IEEE Trans. on Aerospace and Electronic Systems, Vol. 19,
No. 5, pp. 711{729, Sept. 1983.
[9] L. C. P. da Silva, Y. Wang, V. F. da Costa and W. Xu, \Assessment of generator
impact on system power transfer capability using modal participation factors,"
IEE proceedings on Generation, Transmission and Distribution," Vol. 149, No.
5, pp. 564{570, Sept. 2002.
[10] Y. Y. Hsu, C. L. Chen, \Identiflcation of Optimum Location of Stabilizer Ap-
plications Using Participation Factors," IEE Proceedings on Generation, Trans-
mission and Distribution," Vol. 134, pp. 238{44, May 1987.
[11] I. J. Perez-Arriaga, L. Rouco, F. L. Pagola and J. L. Sancha, \The role of
participation factors in reduced order eigenanalysis of large power systems,"
124
IEEE International Symposium on Circuits and Systems, Vol. 927, pp. 923{
927, June 1988.
[12] H. O. Wang, E. H. Abed, and M. A. Hamdan, \Bifurcation, chaos, and crises in
voltage collapse of a model power system," IEEE transactions on Circuits and
Systems- I: Fundamental Theory and Applications, Vol. 41, No. 4, pp. 294{302,
March 1994.
[13] P. Huang and Y. Hsu, \Eigenstructure assignment in a longitudinal power
system via excitation control," IEEE transactions on Power Systems, Vol. 5,
No. 1, pp. 96{102, Feb. 1990.
[14] E. H. Abed and P. P. Varaiya, \Nonlinear oscillations in power systems," In-
ternational J. Electrical Power and Energy Systems, Vol. 6, No. 1, pp. 37{43,
Jan. 1984.
[15] E. H. Abed, D. Lindsay and W. A. Hashlamoun, \On participation factors for
linear systems," Automatica, Vol. 36, No. 10, pp. 1489{1496, Oct. 2000.
[16] F. Milano, \PowerSystem Analysis ToolboxDocumentationfor PSAT," version
1.3.0, 2004 (available online at http://www.power.uwaterloo.ca/~fmilano/).
[17] B. C. Moore, \On the  exibility ofiered by state feedback in multivariable
systems beyond closed loop eigenvalue assignment," IEEE Trans. on Automatic
125
Control, pp. 689{692, Oct. 1976.
[18] M. A. Hassouneh, H. Yaghoobi and E. H. Abed, \Monitoring and control of
bifurcations using probe signals," in Dynamics, Bifurcations and Control, Lec-
ture Notes in Control and Information Sciences, Vol. 273, F. Colonius and L.
Gruene, Eds., Berlin: Springer Verlag, pp. 51{65, 2002.
[19] J. F. Hauer and M. J. Beshir, \Dynamic performance validation in the western
power system," Association of Power Exchanges, Kananaskis, Alberta, Oct.
2000.
[20] C. G. Verghese, I. J. Perez-Arriaga and F. C. Schweppe, \Measuring state
variable participation for selective modal analysis," IFAC symposium on digital
control, New Delhi, India.
[21] M.A. Hassouneah, M. S. Saad and E. H. Abed, \Signal-Based Instability Mon-
itoring of Electric Power Systems," 16th IFAC World Congress, Prague, Czech
Republic, 2005.
[22] E. H. Abed, M.A. Hassouneah and M. S. Saad, \Instability monitoring and con-
trol of power systems," Applied Mathematics for Restructured Electric Power
Systems, J. H. Chow, F. F. Wu and J. A. Momoh, Eds., Springer, pp. 159{178,
2005.
126
[23] G. P. Liu and R. J. Patton, Eigenstructure Assignment for Control System
Design. Wiley, New York, 1998.
[24] C. W. Helstrom, Probability and Stochastic Processes for Engineers, Second
Edition, Macmillan, 1991.
[25] T. Kim and E. H. Abed, \Closed-loop monitoring systems for detecting impend-
ing instability," IEEE Trans. Circuits and Systems{I: Fundamental Theory and
Applications, Vol. 47, No. 10, pp. 1479{1493, Oct. 2000.
[26] N. Martins and L. T. G. Lima, \Determination of suitable locations for power
system stablizers and static VAR compensators for damping electromechan-
ical oscillations in large scale power systems," IEEE Transactions on Power
Systems, Vol. 5, No. 4, pp. 1455{1469, Nov 1990.
[27] E. H. Abed, H. O. Wang and R. C. Chen, \Stabilization of period-doubling
bifurcations and implications for control of chaos," Physica D, Vol. 70, pp.
154{164, Jan 1994.
[28] C. A. Canizares, S. Corsi and M. Pozzi, \Modeling and implementation of
TCR and VSI based FACTS controllers," in Technical Report AT-UCR 99/595,
ENEL Ricerca, Milan, Italy, Dec 1999.
[29] Z. Feng, V. Ajjarapu and D. J. Maratukulam, \A comprehensive approach for
127
preventive and corrective control to mitigate voltage collapse," IEEE Trans. on
Power Systems, Vol. 15, No. 2, pp. 791{797, May 2000.
[30] T. Kim, Noisy Precursors for Nonlinear System Instability with Application to
Axial Flow Compressors. Ph.D. thesis, University of Maryland, College Park,
1997.
[31] I. Genc, H. Schattler and J. Zaborszky, \Clustering the bulk power system with
applications towards hopf bifurcation related oscillatory instability," Electric
Power Components and Systems, Vol. 33, No. 2, pp. 181{198, Feb 2005.
[32] V. M. I. Genc, \Hopf Bifurcation Related Coherent Oscillations in Electric
Power Systems with a Clusterd Texture," Ph.D. Dissertation, Washington Uni-
versity, Saint Louis, MO, Dec 2001.
[33] P. W. Sauer and M. P. Pai, Power System Dynamics and Stability, Prentice
Hall, New Jersey, 1998.
[34] G. C. Verghese, I. J. Perez-Arriaga and F. C. Schweppe, \Selective modal analy-
sis with applications to electric power systems, Part II: The dynamic stability
problem," IEEE Trans. Power Apparatus and Systems, Vol. 101, No. 9, pp.
3126{3134, Sep. 1982.
[35] N. Mithulananthan, C.A. Canizares, J. Reeve and G.J. Rogers, \Comparison
128
of PSS, SVC, and STATCOM Controllers for Damping Power System Oscil-
lations," IEEE Transactions on Power Systems, Vol. 18, No. 2, pp. 786{792,
May 2003.
[36] H. Ghasemi, C. A. Canizares and A. Moshref, \Oscillatory stability limit pre-
diction using stochastic subspace identiflcation," IEEE Transactions on Power
Systems, Vol. 21, No. 2, pp. 736{745, May 2006.
[37] K. Wiesenfeld, \Noisy precursors of nonlinear instabilities," Journal of Statis-
tical Physics, Vol. 38, No. 5-6, pp. 1071{1097, 1985.
[38] K. Wiesenfeld, \Virtual Hopf phenomenon: A new precursor of period doubling
bifurcations," Physical Review A, Vol. 32, Sep. 1985.
[39] C. Jefiries and K. Wiesenfeld, \Observation of noisy precursors of dynamical
instabilities,"Physical Review A, Vol. 31, pp. 1077{1084, 1985.
[40] H. Yaghoobi and E. H. Abed, \Optimal actuator and sensor placement for
modal and stability monitoring," Proceedings of the American Control Confer-
ence, pp. 3702{3707, San Diego, CA, June 1999.
[41] H. Yaghoobi, Stability Monitoring and Control of Nonlinear Systems Susceptible
to Oscillatory Bifurcations. Ph.D. thesis, University of Maryland, College Park,
2000.
129
[42] J. Lu, H. D. Chiang and J. S. Thorp, \Eigenstructure assignment by decentral-
ized feedback control," IEEE Transactions on Automatic Control, Vol. 38, No.
4, pp. 587{594, April 1993.
[43] J. F. Hauer and J. G. DeSteese,\A Tutorial on Detection and Characteriza-
tion of Special Behavior in Large Electric Power Systems,", Paciflc Northwest
National Laboratory Richland, Washington 99352, July 2004.
[44] C. Canizares, \Voltage stability assessment: Concepts, practices and tools,"
Technical Report SP101PSS, Institute for Electrical and Electronic Engineers,
2002.
[45] A. G. Phadke, \Synchronized phasor measurements in power systems," IEEE
Computer App. in Power, Vol. 6, pp. 10{15, April 1993.
130