A Methodology for Estimating the Parameters of Steam Turbine Generator Shaft Systems  
for Subsynchronous Resonance Studies 
 
by 
 
Krishna Sambarapu  
 
 
 
 
A dissertation submitted to the Graduate Faculty of 
Auburn University 
in partial fulfillment of the 
requirements for the Degree of 
Doctor of Philosophy 
 
Auburn, Alabama 
August 4, 2012 
 
 
 
 
Keywords: Steam Turbine Generator, Parameter Estimation,   
Prony?s Method, Subsynchronous Resonance, SSR, Electrical Power Systems 
 
 
 
 
 
 
Approved by 
 
Mark S. Halpin, Chair, Alabama Power Company Distinguished Professor 
R. Mark Nelms, Professor and Chair of Electrical & Computer Engineering Department 
Charles A. Gross, Professor Emeritus, Square D Power Professor 
University Reader, P. K. Raju, Thomas Walter Distinguished Professor of Mechanical  
Engineering 
   
 
 
ii 
 
 
 
 
 
 
Abstract 
 
 
 The increase of coal and nuclear power steam turbines over the past few decades 
combined with transmission line series capacitors creates a potential drawback known as 
subsynchronous resonance (SSR) by increasing the risk of interaction between electrical power 
systems and turbine-generators? rotor torsional systems. Current work presented here adds to a 
growing body of solutions and preventive measures developed by power industry professionals 
to address this problem by proposing a new methodology to estimate the shaft spring constants 
non-intrusively using measurement signals from the shaft systems. 
A significant contribution of this work has been in determining, after a thorough 
investigation of eigen properties of several steam turbine generator shaft systems, which of the 
non-intrusively available measurement signals, is the best candidate to be used for extracting the 
shaft torsional modes. The finding that the use of the high-pressure rotor section speed signal is 
significantly better than the instantaneous power signal in extracting the torsional modes through 
modal decomposition has been verified using measurement signals obtained from simulations, as 
well as real measurements. The modal components, which were extracted using Prony?s method, 
were later used along with rotor section moment-of-inertia values to successfully estimate the 
shaft spring constants of steam turbine generator shaft systems that were simulated in ATP. In 
addition, an algorithm has been developed to reorder the output modes of Prony?s method to 
identify the power system modes. 
 
 
iii 
 
 
 
 
 
 
 
Acknowledgments 
 
 
 I would like to express my deepest gratitude to my advisor and mentor, Dr. Mark S. 
Halpin, without whom this dissertation would not have been possible. I have greatly benefited 
from his expert guidance, support and owe him immensely for shaping my professional identity. 
I would also like to thank my committee members Dr. Mark Nelms, Dr. Charles Gross and Dr. P. 
K. Raju for their inestimable advice and guidance. I am also greatly indebted to Dr. John Hung & 
Late. Dr. Scott Hodel for their invaluable tutelage in state space control systems.  
Most of all, I would like to thank my parents Mr. Sambarapu Narayana and Mrs. 
Yadamma for their unconditional love, support and countless sacrifices. It is my life?s 
profoundest blessing to have been born as their child. I would also like to express  my gratitude 
to my grandparents, Late. Mr. Madanala Balakistaiah and Mrs Pochamma and Late. Mr. 
Sambarapu Chennaiah & Mrs. Laxmi. I would also like to thank my brother-in-law & sister, Mr. 
Botla Maheshwar & Mrs. Anuradha for their constant support and encouragement. I am also 
greatly indebted to my brothers & sisters-in-law: Mr. Narender Kumar & Mrs. Vani and Mr. 
Chandra Mohan & Mrs. Pallavi. I would also like to acknowledge my nieces & nephews for 
making my life richer: Rohith & Haritha, Shashaank & Ojo, Karthik & Nishaant. I would also 
like to thank all my friends at Auburn University & Mississippi State University and all my 
relatives for their wonderful support. 
 
 
 
iv 
 
 
 
 
 
 
 
Table of Contents 
 
 
Abstract ......................................................................................................................................... ii 
Acknowledgments........................................................................................................................ iii 
List of Tables ............................................................................................................................. viii 
List of Figures .............................................................................................................................. xi 
List of Abbreviations ................................................................................................................. xvi 
Chapter 1   Introduction ................................................................................................................ 1 
Chapter 2  Subsynchronous Resonance Review ........................................................................... 5 
 2.1. Types of SSR Phenomena .......................................................................................... 6 
 2.2. SSR Analysis Methods .............................................................................................. 8 
 2.3. Steam Turbine Description  ..................................................................................... 10 
 2.4. TG Shaft System Modeling ..................................................................................... 11 
 2.5. Conclusion ............................................................................................................... 12 
Chapter 3  Review of Prony?s Method ....................................................................................... 13 
 3.1 History and Literature Review .................................................................................. 13 
 3.2 Prony?s Method Description ..................................................................................... 14 
 3.3 Example .................................................................................................................... 21 
 3.4 Conclusions ............................................................................................................... 25 
Chapter 4 Reordering algorithm ................................................................................................. 26 
 4.1 Proposed Reordering Algorithm ............................................................................... 26 
v 
 
 4.2 Example .................................................................................................................... 29 
 4.3 Conclusions ............................................................................................................... 32 
Chapter 5 Effect Of Prony Parameters ........................................................................................ 34 
 5.1 Prony Driver GUI Program Description ................................................................... 34 
 5.2 Effect of Modal Order ............................................................................................... 36 
 5.3 Effect of Sliding Window ......................................................................................... 44 
 5.4 Effect of Sampling Rate ............................................................................................ 48 
 5.5 Guidelines ................................................................................................................. 54 
 5.6 Conclusions ............................................................................................................... 55 
Chapter 6  Eigen Analysis of STG Shaft Systems ...................................................................... 56 
 6.1 State Space Theory Review ...................................................................................... 57 
 6.2 Eigen Values ............................................................................................................. 59 
 6.3 Right Eigen Vectors .................................................................................................. 59 
 6.4 Left Eigen Vectors .................................................................................................... 60 
 6.5 Participation Factors ................................................................................................. 62 
 6.6 Derivation of State Equations for Kundur 1039 System .......................................... 62 
 6.7 Conclusions ............................................................................................................... 73  
Chapter 7  Application of Eigen Analysis On Test Cases .......................................................... 74 
 7.1 Remarks on Terminology Used ................................................................................ 75 
 7.2 System 1: 555 MVA STG System from Kundur pp 1036 ........................................ 76 
      7.2.1 Right Eigen Vector Analysis  ........................................................................... 77 
      7.2.2 Participation Factors ......................................................................................... 79 
      7.2.3 Impulse Response Analysis .............................................................................. 81 
vi 
 
      7.2.4 Conclusions for Kundur 1036 System .............................................................. 83 
 7.3 System 2: 960 MVA STG System from Kundur pp 1038 ........................................ 84 
      7.3.1 Right Eigen Vector Analysis ............................................................................ 85 
      7.3.2 Participation Factors ......................................................................................... 86 
      7.3.3 Impulse Response Analysis .............................................................................. 87 
      7.3.4 Conclusions for Kundur 1038 System .............................................................. 89 
 7.4 System 3: 191 MVA STG System from Kundur pp 1039 ........................................ 89 
      7.4.1 Right Eigen Vector Analysis ............................................................................ 89 
      7.4.2 Participation Factors ......................................................................................... 90 
      7.4.3 Impulse Response Analysis .............................................................................. 91 
 7.5 System 4: 892.4 MVA STG System from IEEE-FBM ............................................. 92 
      7.5.1 Right Eigen Vector Analysis ............................................................................ 93 
      7.5.2 Participation Factors ......................................................................................... 95 
      7.5.3 Impulse Response Analysis .............................................................................. 97 
 7.6 General Conclusions From Right Eigen Vector Analysis ........................................ 98 
 7.7 Overall Analysis and Conclusions using Participation Factors .............................. 101 
 7.8 Overall Analysis and Conclusions using Impulse Response Analysis ................... 103 
 7.9 Conclusions ............................................................................................................. 104  
Chapter 8  IEEE-FBM Full Test System Application and Analysis ......................................... 106 
 8.1 Test Case 1: IEEE-FBM ......................................................................................... 107 
      8.1.1 Simulation in ATP/ATPDraw Software ......................................................... 110 
      8.1.2 Analysis of IEEE-FBM Simulation Results ................................................... 112 
 8.2 Test Case 2: Kundur 1036 System .......................................................................... 118 
vii 
 
      8.2.1 Simulation in ATP/ATPDraw Software ......................................................... 120 
      8.2.2 Analysis of Kundur 1036 System Simulation Results .................................... 120 
 8.3 Conclusions ............................................................................................................. 125  
Chapter 9  Real Measurement Analysis .................................................................................... 126 
 9.1 Data Collection ....................................................................................................... 126 
 9.2 Analysis of Data Set 1............................................................................................. 127 
 9.3 Analysis of Data Set 2............................................................................................. 138 
 9.4 Comparison of Results and Analysis ...................................................................... 142 
 9.5 Conclusions ............................................................................................................. 144 
Chapter 10 Conclusions ............................................................................................................ 146 
Bibliography   ........................................................................................................................... 149 
 
 
  
 
 
  
  
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
viii 
 
 
 
 
 
 
List of Tables 
 
 
Table 3-1 A 6-Sample Record to Illustrate Prony?s Method ...................................................... 21 
Table 4-1 Output of Prony?s Method .......................................................................................... 31 
Table 4-2 Reordered List of Prony?s Method Output Modes ..................................................... 32  
Table 5-1 The LSE Values for Successive Windows ................................................................. 46 
Table 5-2 Effect of Sliding Window on Frequency of Mode1 ................................................... 46  
Table 5-3 Effect of Sliding Window on Frequency of Mode 2 .................................................. 46 
Table 5-4 Effect of Sliding Window on Frequency of Mode 3 .................................................. 46  
Table 5-5  Sampling Rate (nSkipSample) vs LSE - A................................................................ 51 
Table 5-6 Sampling rate (nSkipSample) vs LSE - B .................................................................. 52  
Table 6-1  Rotor Inertias ............................................................................................................. 63 
Table 6-2  Shaft spring constants ................................................................................................ 63  
Table 6-3 S.I. and p.u  Values of Kundur 1039 STG Rotor Section Moment of Inertias........... 67 
Table 6-4 S.I. and p.u. Values of Kundur 1039 STG Shaft Spring Constants ............................ 67  
Table 6-5 Right Eigen Vector Matrix of Kundur 1039 STG Shaft System ................................ 68 
Table 6-6 Three Right Eigen Vectors of Modes 52 Hz, 29 Hz, and 22 Hz Mode ...................... 69  
Table 6-7 Participation Factors of Kundur-pp-1039 System ...................................................... 70 
Table 6-8 Impulse Response Modal Magnitudes of Kundur-pp-1039 System........................... 73  
Table 7-1 Rotor Inertias .............................................................................................................. 76 
Table 7-2 Shaft Spring Constants ............................................................................................... 77  
ix 
 
Table 7-3  Right Eigen Vectors for Kundur-pp-1036 System .................................................... 77 
Table 7-4 Participation Factors of Kundur-pp-1036 System ...................................................... 80  
Table 7-5 Impulse Response Modal Magnitudes of Kundur-pp-1036 system ........................... 82 
Table 7-6 Rotor Section Inertias ................................................................................................. 84  
Table 7-7 Shaft Spring Constants ............................................................................................... 85 
Table 7-8  Right Eigen Vectors for Kundur-pp-1038 System .................................................... 85  
Table 7-9 Participation Factors of Kundur-pp-1038 System ...................................................... 87 
Table 7-10 Impulse Response Modal Magnitudes of Kundur-pp-1038 System......................... 88  
Table 7-11 Rotor Spring Mass Parameters ................................................................................. 93 
Table 7-12  Right Eigen Vectors for IEEE-FBM System........................................................... 94  
Table 7-13 Participation Factors of IEEE-FBM Machine .......................................................... 95 
Table 7-14 Impulse Response Modal Magnitudes of the IEEE-FBM System ........................... 97  
Table 8-1 Generator Impedances, Time Constants and Operating Conditions ........................ 108 
Table 8-2 Rotor Circuit Parameters .......................................................................................... 109  
Table 8-3 Rotor Spring Mass Parameters ................................................................................. 109 
Table 8-4 List of SSR frequencies ............................................................................................ 110  
Table 8-5a Significant Constituent Modal Components in the FBM?s Power Signal  ............. 111  
Table 8-5b Shaft Torsional Modal Components in the FBM?s Power Signal  ......................... 114 
Table 8-6 Comparison Between Estimated and FBM?s Actual Frequencies in p(t) signal ...... 114  
Table 8-7a Significant Constituent Modal Components in FBM?s HP Speed Signal .............. 116  
Table 8-7b Shaft Torsional Modal Components in FBM?s HP Speed Signal .......................... 116 
Table 8-8 Comparison Between Estimated vs. Actual Shaft Spring Constant Values (p.u.) ... 118  
Table 8-9 Specifications for Kundur 1036 Shaft System Parameters....................................... 119 
x 
 
Table 8-10 Kundur 1036 System Generator?s Electrical Data ................................................. 119  
Table 8-11a Significant Constituent Modes in Kundur 1036?s p(t) ......................................... 121  
Table 8-11b Shaft Torsional Modes in the Kundur1036?s p(t)................................................. 121  
Table 8-12a Significant Constituent Modes in Kundur 1036?s HP Speed Signal .................... 123  
Table 8-12b Shaft Torsional Modes in the Kundur1036?s HP Speed Signal ........................... 123 
Table 8-11b Results of Modal Decomposition on the HP Rotor Section Velocity Signal ....... 118  
Table 8-13 Comparison Between Estimated and Actual Spring Constant Values ................... 124 
Table 9-1 Prony?s Method?s Output List of Constituent Modal Components for Data Set 1 .. 131  
Table 9-2 Selected List of Constituent Modes in the Signal Detected by Prony?s Method ..... 133 
Table 9-3 List of Constituent Modes in the Signal Detected by PM with Algorithm .............. 134  
Table 9-4 List of Constituent Modes in the Signal Detected by PM with Algorithm .............. 134 
Table 9-5 Combined List Modal Components ......................................................................... 135  
Table 9-6 Prony?s Method?s Output List of Constituent Modal Components for Data Set 2 .. 140 
Table  9-7 Selected List of Constituent Modes in the Signal Detected by Prony?s Method .... 140  
Table  9-8 Results of Modal Decomposition Using Prony?s Analysis on the Two Data Sets .. 143 
 
 
 
 
 
 
 
 
xi 
 
 
 
 
 
 
List of Figures 
 
 
Figure 3-1 Lumped Spring-Mass Model of the Turbine-Generator Shaft  ................................. 15 
Figure 4-1 Flowchart for Reordering Algorithm ........................................................................ 29  
Figure 4-2 MATLAB Code for Creating the Sample Waveform ............................................... 29  
Figure 4-3 Plot of the Sample Waveform. .................................................................................. 30  
Figure 4-4 MATLAB Code for Generating the LSE Norms ...................................................... 31  
Figure 5-1 GUI Built in MATLAB to Allow User Selection ..................................................... 35  
Figure 5-2 Sample Waveform Consisting Of Three Sinusoids and Noise. ................................ 36  
Figure 5-3 Comparison of Actual Signal and PM?s Output: Stem Plot Showing the LSE    
          Magnitudes of Constituent Sinusoids. .............................................................................. 37  
Figure 5-4a Comparison of Actual Signal and PM?s Output along with Stem Plot of the  
          Constituent Modes? LSE Magnitudes with numModes Parameter = 10 .......................... 38  
Figure 5-4b Comparison of Actual Signal and PM?s Output along with Stem Plot of the 
          Constituent Modes? LSE Magnitudes with numModes Parameter = 20 .......................... 39  
Figure 5-4c Comparison of Actual Signal and PM?s Output along with Stem Plot of the 
          Constituent Modes? LSE Magnitudes with numModes Parameter = 30 .......................... 39  
Figure 5-4d Comparison of Actual Signal and PM?s Output along with Stem Plot of the 
          Constituent Modes? LSE Magnitudes with numModes Parameter = 60 .......................... 40  
Figure 5-4e Comparison of Actual Signal and PM?s Output along with Stem Plot of the  
          Constituent Modes? LSE Magnitudes with numModes Parameter = 70 .......................... 40  
Figure 5-4f Comparison of Actual Signal and PM?s Output along with Stem Plot of the  
          Constituent Modes? LSE Magnitudes with numModes Parameter = 80 .......................... 41  
Figure 5-4g Comparison of Actual Signal and PM?s Output along with Stem Plot of the  
          Constituent Modes? LSE Magnitudes with numModes Parameter = 90 .......................... 41  
xii 
 
Figure 5-5 Plot of LS Error as a Function of numModes ........................................................... 42  
Figure 5-6 Total LSE vs numModes with Noise as Varying Parameter for a 3-Sinusoid .......... 43  
Figure 5-7 Total LSE vs numModes with Noise as a Varying Parameter for a 1-Sinusoid ....... 44  
Figure 5-8 Comparison of given and estimated waveforms ....................................................... 45  
Figure 5-9a Comparison of Estimated Vs. Actual Magnitudes of  2.5 Hz Mode ....................... 47  
Figure 5-9b Comparison of Estimated Vs. Actual Magnitudes of  15 Hz Mode ........................ 47  
Figure 5-9c Comparison of Estimated Vs. Actual Magnitudes of  180 Hz Mode ...................... 47  
Figure 5-10a Comparison of Actual Signal and PM?s Output Along With Stem Plot of the  
          Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 1 ......................... 49  
Figure 5-10b Comparison of Actual Signal and PM?s Output Along With Stem Plot of the  
          Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 2 ......................... 49  
Figure 5-10c Comparison of Actual Signal and PM?s Output Along With Stem Plot of the  
          Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 3 ......................... 50  
Figure 5-10d Comparison of Actual Signal and PM?s Output Along With Stem Plot of the  
          Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 4 ......................... 50  
Figure 5-10e Comparison of Actual Signal and PM?s Output Along With Stem Plot of the  
          Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 5 ......................... 51  
Figure 5-11 Effect of Sampling Rate on the Relationship Between numModes and LSE ......... 53  
Figure 6-1 Example of a 4 State System Consisting of Inductors and Capacitors ..................... 58  
Figure 6-2 Schematic Layout of Rotors and Shafts of Kundur pp. 1039 Problem ..................... 63  
Figure 6-3 Electrical Equivalent of Shaft system ....................................................................... 64  
Figure 6-4 ATPDraw Schematic of Simulation of Impulse Response of Kundur 1039 Shaft  
          System  .............................................................................................................................. 71  
Figure 6-5 Plot of HP Rotor Section Velocity ............................................................................ 72  
Figure 7-1 Schematic Layout of Rotor and Shaft Sections of Kundur pp. 1036 system ............ 76  
Figure 7-2 Kundur 1036 Right Eigen Vector (mode shape) Absolute Magnitudes Across Rotor  
          Sections ............................................................................................................................. 78  
xiii 
 
Figure 7-3 Kundur 1036 System?s Participation Factors Across All Rotor Sections ................. 80  
Figure 7-4 Kundur 1036 System Impulse Response Results ...................................................... 83  
Figure 7-5 Schematic layout of STG shaft System of Kundur pp. 1038 problem ...................... 84  
Figure 7-6 Kundur 1038 Right Eigen Vector Absolute Magnitudes Across Rotors Sections .... 85  
Figure 7-7 Kundur 1038 System?s Participation Factors Across All Rotors Sections ............... 87  
Figure 7-8. Kundur 1038 System Impulse Response Results ..................................................... 88  
Figure 7-9  STG Shaft System Schematic for Kundur 1039 System .......................................... 90  
Figure 7-10 Kundur 1039 Right Eigen Vector Absolute Magnitudes Across Rotor sections .... 89  
Figure 7-11 Kundur 1039 System?s Participation Factors Across All Rotor Sections ............... 91  
Figure 7-12 Kundur 1039 System Impulse Response Results .................................................... 92  
Figure 7-13 STG Shaft System for IEEE-FBM .......................................................................... 93  
Figure 7-14 IEEE-FBM Right Eigen Vector Absolute Magnitudes Across Rotors Sections .... 94  
Figure 7-15 IEEE-FBM Participation Factors Across All Rotor Sections ................................. 95  
Figure 7-16 IEEE-FBM Participation Factors Across Generator and HP Rotor Sections ......... 96  
Figure 7-17 IEEE-FBM Impulse Response Results ................................................................... 97  
Figure 7-18 Comparison of Mode Shapes of Highest Frequency Modes ................................... 98  
Figure 7-19 Comparison of Mode Shapes of Highest Frequency Modes ................................... 99  
Figure 7-20 Comparison of Mode Shapes of Highest Frequency Modes ................................. 100  
Figure 7-21 Comparison of Mode Shapes of Highest Frequency Modes ................................. 101  
Figure 7-22 Comparison of Participation Factors of Highest Frequency Modes ..................... 102  
Figure 7-23 Comparison of Mode Shapes of Highest Frequency Modes ................................. 103  
Figure 7-24 Improvement of Using the HPRSV Signal Instead of the GNRSV Signal ........... 104  
Figure  8-1 FBM Test Case?s Network ..................................................................................... 107  
xiv 
 
Figure  8-2 Lumped Spring Mass Model of the Steam Turbine ............................................... 110  
Figure 8-3 ATPDraw Schematic Showing the IEEE FBM System .......................................... 111  
Figure 8-4 Instantaneous Power p(t) ......................................................................................... 112  
Figure 8-4b Comparison of Actual Vs Estimated FBM?s p(t) Signals ..................................... 113  
Figure  8-5 Speed Signal of the HP Rotor Section ................................................................... 115  
Figure  8-6 Waveform Comparison of Actual vs. Estimated HP signal ................................... 117  
Figure 8-7 Kundur 1036 Shaft System Configuration .............................................................. 119  
Figure  8-8 ATPDraw Schematic of the Kundur 1036 System ................................................ 120  
Figure  8-8b Comparison of Actual Vs Estimated Waveforms of Kundur 1036 System?s p(t) 122  
Figure  8-8c Comparison of Actual Vs Estimated Waveforms of Kundur 1036 System?s HP  
          Rotor Section Speed Signal ............................................................................................ 124  
Figure 9-1 Plot of Actual p(t) Measurement Data .................................................................... 128  
Figure  9-2 Plot of Filtered p(t) Signal from Data Set 1 ........................................................... 129  
Figure  9-3 Zoomed-in Plot of the Disturbance Event in the p(t) Signal .................................. 129  
Figure  9-4 Plot of Proper Candidate Data Set for Prony?s Method (Data Set 1) ..................... 130  
Figure 9-5 A Plot of Comparison Between the Actual and Reconstructed Signal Using the Output  
          Modes of Prony?s Method for Data Set 1 ....................................................................... 132  
Figure  9-6 Scatter Plot of Constituent Modes of Data Set 1 .................................................... 136  
Figure 9-7 Scatter Plot of Constituent Modal Components with Major Clusters Demarcated to  
          Show the Four Dominant Modes .................................................................................... 137  
Figure 9-8 Histogram of Modal Components for Data Set 1 with Bin Width = 1 Hz .............. 137  
Figure 9-9 Plot of Filter p(t) signal from Data Set 2................................................................. 138  
Figure 9-10 Plot of Proper Candidate Data Set for Prony?s Method (Data Set 2) .................... 139  
Figure 9-11 A Plot of Comparison Between the Actual and Reconstructed Signal Using the  
          Output Modes of Prony?s Method for Data Set 2 ........................................................... 139  
xv 
 
Figure 9-12 Scatter plot of Modal Components with Four Major Clusters Demarcated.......... 141  
Figure 9-13 Histogram of Modal Components in Data Set 2 with Bin Width = 1 Hz ............. 142  
Figure 9-14 Steam Turbine Generator Shaft System Configuration of the Studied System .... 144 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
xvi 
 
 
 
 
List of Abbreviations 
 
 
SSR Subsynchronous resonance 
LTI Linear time invariant 
IEEE Institute of Electrical and Electronics Engineers 
FBM First Benchmark Model 
SBM Second Benchmark Model 
STG Steam turbine generator  
STGS Steam turbine generator shaft 
RSV Rotor section velocity 
PM Prony?s Method 
GUI Graphical user interface 
HP High pressure 
 
IP Intermediate pressure   
 
LPA Low pressure-A 
 
LPB Low pressure-B 
 
GENR Generator 
 
 
 
 
 
 
CHAPTER 1 
INTRODUCTION  
 
Series capacitor compensation is a scheme used in power systems to alleviate problems 
that arise because of high reactance of long transmission lines. By reducing effective reactance, it 
increases transient stability and improves the load carrying capability of these lines. By 
controlling this reactance it allows better control over load sharing between parallel transmission 
lines. Because of the advantages it provides economically, it has been used widely in power 
systems. However, with the increase of coal and nuclear power steam turbines over the past few 
decades, series capacitors have become a potential drawback by increasing the risk of interaction 
between electrical power systems and turbine-generators? rotor torsional systems- a phenomenon 
known as subsynchronous resonance (SSR) or subsynchronous oscillation. ?Subsynchronous 
oscillation is an electric power system condition where the electric network exchanges 
significant energy with a turbine-generator at one or more of the natural frequencies of the 
combined system below the synchronous frequency of the system following a disturbance from 
equilibrium [1].? When these energy interchanges manifest themselves in rotors as pulsating 
torques, they have the potential to cause irreversible damage to rotor shaft systems such as cracks 
and fractures. Two such severe incidents occurred in the early 1970s in Nevada's Mohave 
Generating Station resulting in damaged rotor shafts [2]. As a result, power industry 
professionals began building a framework of solutions and preventive measures. They developed 
2 
 
definitions [1,3] of the SSR problems, identified the different aspects of SSR, developed the 
IEEE First Benchmark Model (FBM) [4] and the IEEE Second Benchmark Model (SBM) [5] 
and performed tests on some turbo-generating (TG) units.  
Building preventive measures around SSR frequencies using test simulations requires 
accurate input data. This data consists of steam turbine rotor shaft system parameters such as 
moments of inertia, shaft spring constants and the effective dampings due to mechanical 
hysteresis and system loads. In most cases, the damping is negligible. While moments of inertia 
can be determined with considerable accuracy, the spring constants are not always precisely 
known. Traditionally, turbo-generator manufacturers do not provide this information. When 
requested, the manufacturers estimate or compute these parameters by performing extensive 
simulations based on the rotor?s material composition and physical dimensions and other design 
parameters. These simulations and tests can also be performed by a third party. Because of the 
extent of effort and time required, determining spring constants is usually an expensive task. The 
main focus of this research is to develop a measurement-based methodology to estimate these 
spring constants without necessitating expensive and intrusive tests. A key advantage of this 
online measurement method is that it does not require the unit to be taken offline for testing 
procedures. 
This methodology consists of a two-step approach:  
? Monitor and record measurements of the instantaneous power at the generator terminals 
(or available rotor velocity signals) whenever a fault or disturbance occurs. When a large 
disturbance occurs, it tends to excite all the torsional modes resulting in sustained rotor 
shaft oscillations at these frequencies for a few seconds. These oscillations can be 
detected in the measured signals. This recorded ring-down signal is fed to a Prony based 
3 
 
algorithm to extract the modal information - amplitude, phase, frequency and exponential 
damping - from the signal.  
? Using the extracted frequencies and known TG system rotor masses, the necessary shaft 
spring constants are estimated using a non-linear least squares algorithm. 
The different types of SSR phenomena and the types of analysis available to detect these 
phenomena are reviewed in Chapter 2. Various techniques available to model the entire shaft 
system will be examined after reviewing the mechanical features of the STG shaft system. 
One of the main features of this work is the use of Prony's method for modal decomposition 
instead of the FFT because of its various advantages over the FFT. In Chapter 3 is provided a 
thorough review, analysis and implementation of Prony's method. 
During the course of this work, it was observed that there tend to exist certain general 
characteristics of STG torsional modes when measured in rotor speed signals or 
instantaneous power signals. This observation led to the development of a new algorithm, 
which reorders the output modes of Prony's method and enables the user to differentiate these 
torsional modes from extraneous signals such as noise etc. This algorithm is presented in 
Chapter 4. 
There are many advantages of Prony's method over the FFT. For instance, it provides 
modal damping information whereas the FFT does not. However, there are some 
disadvantages in using Prony's method. The accuracy of its results is compromised by a 
sensitivity to noise. Moreover, as with any system identification tool, its parameters have to 
be chosen carefully. These practical issues and their solutions are examined in Chapter 5.  
 
4 
 
Before the proposed methodology is tested using field data and data obtained from 
simulations of fully connected test systems consisting of the STG shaft systems and electrical 
machine models of generators and connected electrical networks, it is very beneficial to study the 
properties of isolated STG shaft systems. Treating the isolated lumped spring mass model of the 
STG shaft system as an LTI system, the following three analysis are performed on it: 
? Right eigen vector analysis, 
? Participation factor analysis, and 
? Analysis of impulse response determined by time-domain simulation. 
These three analyses are illustrated with the aid of an example shaft system in Chapter 6. The 
results of a thorough investigation of these three analyses on four STG shaft systems whose data 
was provided in the literature are presented in Chapter 7. All these STG shaft systems have very 
similar configuration and overall general characteristics with respect to rotor section moments-
of-inertia.  This has revealed some extremely useful information that has direct implication to the 
methodology presented. 
In Chapter 8, two test cases will be utilized to illustrate and validate the proposed 
methodology. These test cases have been built and simulated in ATP and ATPDraw software, 
which are used by power industry professionals to perform time-domain simulations. 
In Chapter 9, the analyses of real field data and the results of modal extraction performed 
using Prony?s method are presented. The real field data consists of instantaneous three-phase 
power measurements that were obtained from a steam turbine generator unit whose shaft system 
consists of six rotor sections. It is important to note that this was the only type of field data 
available for this dissertation work. Rotor section speed signals were not available. 
 
 
 
 
 
 
CHAPTER 2 
SUBSYNCHRONOUS RESONANCE REVIEW 
 
Although subsynchronous resonance was first mentioned in the literature in 1937 while 
discussing series capacitors in transmission lines and studying negative damping of electrical 
machinery [6][7], there were no issues/problems regarding SSR, so it was ignored until 1971. In 
that year, SSR caused the shaft failure of a turbine-generator shaft in the Mohave Generating 
Station in Nevada. 
Tacoma Narrows bridge (Washington, US) collapsed in 1940 because of resonance 
induced by wind, which coincided with the bridge's natural frequency, causing extremely large 
amplitude swings and ultimately the collapse of the bridge. SSR is a similar dynamic 
phenomenon which involves the TG mechanical shaft system and the electrical network that 
contains series capacitors coupled electromagnetically to the generator described by Park's 
equations. Electrical quantities in any electrical circuit can be classified into two components, (1) 
forced and (2) natural, as follows: 
)c o s ()c o s ()( 222111 2 ???? ? ???? ? teAtAti t 
While the ?1 frequency is determined solely by the network's driving or active sources, the ?2 
frequency is dependent completely on the network parameters. In a transmission system, these 
network parameters include the transmission line reactance, capacitance and the sub transient 
reactance of the generator. These currents flowing through the generator armature create slip 
6 
 
frequency torques in the TG rotor because of Park's equations. When the frequencies of these 
slip frequency torques nearly coincide with any of the natural frequencies of the TG shaft 
system, it has the potential to cause serious damage, including loss of fatigue life or even cracks 
in the rotor. The Mohave incident refers to a TG shaft system failure that occurred in the early 
1970s because of self-excited torsional oscillations in the TG shaft system.  It occurred in the 
500 KV coal-fired Mohave Generating Station in Nevada, U.S.A., which consisted of two 798 
MW steam turbine generating units. The generating station was radially connected to the load, 
and the fault occurred when one of the two radial long transmission lines was opened with the 
other line containing series capacitors in place. The opening of the line induced heavy currents in 
the other line, which, along with the generator, transmission line and series capacitor, formed an 
LC tank circuit and thus induced 29.5 Hz current in the generator armature. These currents in 
turn induced 30.5 Hz undamped torques in the TG shaft. The resulting damage included 
circumferential cracks one inch deep, melting or burning of the TG shaft and burning of collector 
rings. One of the immediate remedies included removal of series capacitors and installing a TEX 
type relay, which tripped the unit when SSR current flowed in the armature.  But, because of the 
loss of transient stability, the removal of the series capacitors was not deemed correct, and a 
permanent solution of installing resistive filters was proposed.   
 
2.1 Types of SSR Phenomena 
There are many types of subsynchronous resonance phenomena. Three of them are 
considered important, namely: 
1. Induction generator effect, 
2. Torsional interaction, and  
7 
 
3. Transient torque. 
 
As mentioned before, these are SSR phenomenon, i.e. they include series capacitors in lines. 
The first two can be considered small signal phenomena (during the initial stages of a fault) and 
hence can be analyzed using linear analysis tools like Eigen analysis. The third is a non-linear 
phenomenon which will require a more detailed study. The first one does not include TG modal 
interaction, i.e. it occurs when the rotor is running at synchronous speed. The latter two involve 
TG modal interaction.  
 
Induction generator effect  
Any disturbance causes transient currents to flow through the electrical network and generator 
armature at the combined network?s natural frequencies, both above and below the synchronous 
frequency. When the subsynchronous frequency currents flow through the generator armature, 
they view the synchronously rotating rotor?s circuit as negative resistance. If this negative 
resistance is greater than the positive resistance encountered in the network and generator 
armature, the net negative resistance will give rise to high amplitude sustained subsynchronous 
currents. 
 
Torsional Interaction 
This phenomenon involves torsional modal interaction. As mentioned previously, any electrical 
disturbance causes subsynchronous fault currents to flow through generator armature. When 
these induce subsynchronous torques in the rotor which coincide with the rotor?s natural modes, 
they induce rotor oscillations at these subsynchronous frequencies. These induce voltages in the 
8 
 
generator armature at the subsynchronous frequencies and are sustained when they are in phase. 
If the magnitudes of these self-excited subsynchronous torques are large enough to overcome the 
overall mechanical damping involved, they will have the potential to cause large amplitude 
oscillations that can damage the rotor. This phenomenon was responsible for the Mohave 
incident. 
 
Transient Torque 
This non-linear phenomenon is caused by large short circuit currents that oscillate at one of the 
system?s natural frequencies (resulting from system disturbances). Because resultant rotor 
torques induced are directly proportional to these disturbance currents, the resulting peak torques 
can be quite large and have the potential to cause irreversible damage to the rotor.             
 
2.2 SSR Analysis Methods 
Frequency Scanning Technique 
The frequency scanning technique is a small signal analysis tool that computes the reactance and 
resistance as a function of frequency looking into the network from behind the stator winding. It 
is especially suitable for detecting IGE conditions. When, at a particular frequency, the computed 
reactance is nearly equal to zero and the resistance is less than zero, i.e. negative, then IGE 
induced stator currents can be expected. It can also aid in detecting torsional interaction and 
transient torque phenomenon, which could occur when there is series resonance or a reactance 
minimum. Because of its simplicity, it is widely used in preliminary studies in North America. 
 
 
9 
 
 
Eigen Analysis 
Another popular method of studying SSR problems is eigen analysis. In this method, the entire 
system is mathematically modeled using a single set of linear differential equations. These linear 
equations are solved using numerical analysis tools like singular value decomposition (SVD) to 
yield the eigen values of all system components involved. This has the advantage over frequency 
scanning techniques in that the eigen values provide damping as well as frequency values; any 
positive real values of eigen values will immediately inform the engineer of a potential unstable 
mode in the system. 
 
EMTP/ATP/ATPDraw analysis 
The EMTP/ATP simulation program is a very powerful software tool for analyzing SSR 
conditions. Because it uses stepwise numerical integration of differential equations, it allows the 
user to enter into the program detailed attributes of each component in the system. Instead of 
using linear approximations (like eigen analysis and frequency scan), this software models all the 
non-linear attributes of complex components involved in the system study. Also, in comparison 
to other stability programs, which use per phase or positive sequence models, it uses three-phase 
information allowing simulation/modeling of unbalanced events. Because it can model non-
linear phenomenon it is most suitable for transient torque analysis. The resulting plots of the 
simulation results yield amplitude information of electrical quantities, in addition to frequency 
and damping. 
 
 
10 
 
 
2.3 Steam Turbine Description 
The invention of steam turbine is attributed to Carl Gustav de Laval of Sweden and 
Charles Parson of Great Britain during the 1880s independently. Over the decades, steam 
turbines have grown in size and complexity. Whereas, the typical ratings for large steam turbines 
prior to the 1950s was less than 200 MW, modern large ST ratings can exceed 1200 MW. In 
contrast to older units, most modern STs rarely contain gear boxes and so their total length can 
exceed 150 ft (in the case of tandem units) and their weight can be in the range of several 
hundred tons  [8]. The rotors are usually made of Cr.Mo.V steel and are either forged or welded. 
The entire rotor system can contain, in addition to the rotor, a multitude of other parts such as 
turbine blades, discs, retaining rings, etc.. As a result, the rotor system might contain an infinite 
number of modes of torsional vibration  [8]. However, most of these modes interact with higher 
frequency electrical torques and can be safely ignored for SSR studies. Because of the extremely 
large power output (>1000 MW), the system is divided into multiple stages of turbine cylinders 
and shafts, each producing a part of the total output. There are generally two such types of 
compound units (1) Tandem-compound and (2) Cross compound units. In tandem compound 
units, all the cylinders are on a single shaft. Consequently, the total length of the shaft can be 
very large, since it includes many stages (sometimes exceeding 10 stages) in addition to the 
generator and exciter. The cross-compound unit contains two shafts which allow for lesser rating 
and reduced total shaft length. However, they have a disadvantage in that they need two separate 
sets of generator and exciter controllers, one for each shaft.  
 
 
11 
 
 
2.4 TG Shaft System Modeling 
  SSR is caused by series capacitors in transmission lines. This compensation is such that 
the total capacitive reactance is less than the total line reactance (usually less than 70%). This 
leads to SSR frequencies less than 60 Hz. This, in turn, greatly helps simplify the modeling 
requirements of the TG shaft. 
 
There are three types of TG system models. They are: 
1. Lumped spring mass models, 
2. Advanced continuum models, and 
3. Finite element analysis (FEA) models. 
 
While the first two are theoretical analytical models, the third one is a 
computational/numerical model. All three models lead to methods which involve assembling sets 
of linear differential equations and applying equations of equilibrium based on Hooke?s law and 
Newton?s laws of dynamics to determine eigen values and eigen vectors (mode shapes).   
In the lumped spring mass model, each stage of the turbine (HP, IP etc) as a whole is 
modeled as a rigid rotating mass, which is connected by massless springs representing the shaft 
sections. In the advanced continuum model, each rotor disk in each turbine stage is modeled as a 
rigid mass. Thus, the entire model might contain hundreds of linear differential equations. This 
fact is further confirmed in FEA analysis. 
Finite element analysis is a computational technique generally used by manufacturers for 
mechanical design purposes. Here, all the rotor?s segments, such as turbine blades, retaining 
12 
 
rings, etc, are modeled as discrete regions linked in a geometric model representing the entire 
turbine generator shaft. The results of these models reasonably match those of the basic lumped 
spring mass model. This allows us to use the lumped mass model in SSR analysis methods. 
 
2.5 Conclusion  
In conclusion, the history, the theory of SSR and the different types of analysis tools 
commonly used in power industry are reviewed in this chapter. Eigen analysis and time-domain 
simulation methods using ATP/EMTP have been used extensively to examine the characteristics 
of the shaft torsional and electrical systems involved. While FEA and advanced continuum 
models are useful, this work, like most SSR analysis in power industry, uses lumped spring-mass 
models, which yields sufficient information to aid in analysis and building of SSR preventive 
measures. Prony's method, reviewed in the next chapter, is then used for modal decomposition of 
signals obtained from time-domain simulations of lumped spring-mass models. 
 
 
 
 
 
 
 
 
CHAPTER 3 
REVIEW OF PRONY?S METHOD 
 
Prony?s method used in this work is an exponential parameter estimation technique. It is 
used to fit a linear combination of exponential terms to a uniformly sampled signal, thereby 
extracting the four parameters for each mode from the signal: amplitude, phase, frequency and 
exponential damping [9][10][11]. This technique was developed in the 1790s by Baron de Prony. 
However, because it was computationally expensive, it was not applied until the advent of 
modern fast and efficient computers. Ever since, it has been a practical tool in the field of 
modern spectrum analysis [11]. However, it is only during the last two decades that it found 
practical use in power system applications, such as extracting modal information from transient 
stability programs and system field tests [10]. Because this method also extracts modal damping 
coefficients, it has an advantage over Fourier analysis for the analysis of electromechanical 
oscillations that routinely involve exponentially decreasing (or increasing) waveforms.  
 
3.1 History and Literature Review 
Prony's Method was invented by Gaspard Clair Fran?ois Marie Riche de Prony, who was 
a French mathematician at ?cole Polytechnique in the 18th century. He used the method as 
described in [12], which fitted a sum of exponentials to a set of equally spaced data points, in his 
work on expansion of gases. After the original paper by Prony, [9] is one of the earliest texts on 
14 
 
the subject in the 1900s. Hildebrand wrote this in the 1950s as part of describing other 
approximation procedures in this numerical analysis book. It is presented as an exponential 
approximation technique, which simplifies the solution of a set of non-linear equations into a 
much simpler three-step process involving two least-squares solutions and one polynomial root- 
finding process. After that, Prony?s analysis (PA), was considered a good tool among the 11 
tools used for spectrum analysis, along with Blackman-Tukey, ARMA etc.  
Different interpretations were provided for PA. It was considered that the exponential 
approximation is the homogenous solution of a non-homogenous constant coefficient linear 
difference equation. It was found that pre-filtering vastly improved the accuracy of PA results 
because PA is found to be highly sensitive to noise. In one application it was used after applying 
a moving average smoothing filter to the extraction of modal parameters in South African tie-line 
power swings and German faulty neutral-to-ground resonant voltage measurements.   
Kumaresan showed that it was possible to create a pre-filter out of the data samples 
themselves and used it to pre-filter the data samples. In [10] Prony's method was introduced to 
the field of power systems. In it, the authors introduce Prony?s method to the field of power 
system computations as a complement to eigen analysis and Fourier analysis. Also in the paper, 
the authors describe the method mathematically, list the strategic issues and shortcomings, and 
partially validate, with test cases, its viability as a modal analysis tool for transient stability 
program output. 
 
3.2 Prony?s Method Description 
Prony?s method is an exponential approximation technique that is used to fit a function 
consisting of a linear combination of exponentially damped sinusoids to a given linearly spaced 
15 
 
data sample set [9,10,11]. Extensive developments during the past few decades in linear 
prediction theory, especially in speech processing, have provided new interpretations to the 
original method [13,14,11].  
Prony's method can be best explained using the example of an LTI system such as the 
lumped spring-mass model of the turbo-generator rotor shaft shown in Figure 3-1.  Its dynamics 
are represented by the differential equation for the jth rotor mass shown in (3.1) 
                     ? ? ? ?
kjjkjiijjjj KKTdtdJ ????? ???????????
 
? ? ? ?ijjikjjkjj DDD ????? ??????????                    (3.1) 
In (3.1), Ji and Kij represent the moment of inertia of the ith rotor section and the spring constant 
of the shaft connecting the ith and jth rotor sections, respectively. Di and Dij represent the damping 
coefficients. ??i, ??i and ?Ti represent the incremental angular velocity, angular displacement 
and input torque of the ith rotor mass, respectively. 
 
 
 
 
 
 
 
 
 
Figure 3-1. Lumped Spring-Mass Model of the Turbine-Generator Shaft [15] 
 
u(t) y(t) 
Kij Kjk 
??k ??j ??i 
?Ti ?Tj ?Tk 
Jk Jj Ji 
Dij Djk 
Di Dj Dk 
16 
 
These dynamics can be rewritten in state space form as given in (3.2). 
                           )(~][)(~ txAdttxd ?                                                             (3.2) 
In (3.2), the )(~tx  are the states of the system consisting of n states (??i & ??i) and [A] is the 
state transition matrix. Assuming the states of the system are brought to an initial state )(~0tx  and 
any further disturbances or inputs (?Ti)  are removed, then the states of the system can be 
described by (3.3).  
? ?? ? ? ? ? ?? ??? ?? ?? ni tini tiTi ii etxReptxqtx 1 01 0 ~~~~)(~ ??                                        (3.3) 
In (3.3), ?i, ip~  and iq~  are the ith eigen value, right eigen vector and left eigen vector, respectively  
[10].  
The output of this system y(t) is given by (3.4).  
? ? ? ? ? ?? ? ? ? ? ? ? ?? ??? ?? ??? ni tini tiTi ii etxRCeptxqCtxCty 1 01 0 ~~~~)(~)( ??                             (3.4) 
In (3.4), [C] and [Ri] are the output matrix of the LTI system and its ith residue matrix. 
 Let y(k) be the measured signal for y(t) consisting of N equally spaced samples at intervals 
k=0,1,2,?N-1. The y(t) equation shows this exponential nature of individual components in the 
output. Prony?s method tries to address this nature by fitting a similar structure to the measured 
output as shown in (3.5). 
  ?
??
p
i
ti ieBty
1)(?
?                                                              (3.5) 
This exponential structure, which can also be presented as a sum of exponentially damped 
sinusoids shown in (3.6), is thus suitable to estimate the individual properties like frequencies 
and dampings of modes in the turbo-generator shaft?s output signals. 
17 
 
?? ?? Qi iiti tfeAty i1 )2c o s ()(? ???                                           (3.6) 
In (3.6), Ai, ?i, fi, ?i are the amplitude, exponential damping, frequency and phase values of the ith 
modal component.  
The discrete-time version )(?ky  for )(?ty can be written in a convenient form as shown in (3.7). 
     ?
??
P
i
kii zBky
1)(?
                                                             (3.7) 
In (3.7), zi is the exponential part of the modal component as shown in (3.8). 
ti iez ?? ?                                                                  (3.8) 
Expanding (3.7) for sample instants k=0,1,2? N-1gives (3.9). 
                                                 nBBBBy ????? ?321)0(  
nn zBzBzBzBy ????? ?332211)1(                                       (3.9) 
     2233222211)2( nn zBzBzBzBy ????? ?  
     ?  
                                          1133122111)1( ???? ?????? NnnNNN zBzBzBzBNy ?   
Rewriting (3.9) in matrix form gives (3.10). 
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
??
?
?
?
?
?
?
??
?
?
?
?
?
??
?
?
?
?
?
??
??
??
?
?
?
?
?
?
?
?
?
?
?
?
?
? ????? nNnNN
n
n
N
nn
N
nn
nn
B
B
B
zzz
zzz
zzz
zBzB
zBzB
zBzB
Ny
y
y
?
?
????
?
?
?
?
?
?
?
2
1
11
2
1
1
11
2
1
1
00
2
0
1
11
11
11
11
00
11
)1(
)1(
)0(
                (3.10) 
[Y]=[Z][B]                                                                  (3.11) 
In (3.11), Y, B and Z represent matrices as shown in (3.12). 
 
 
18 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)1(
)1(
)0(
Ny
y
y
Y ? , 
?
?
?
?
?
?
?
?
?
?
?
?
?
nB
B
B
B ?2
1
, 
??
?
?
?
?
?
??
?
?
?
?
?
?
??? 11
2
1
1
11
2
1
1
00
2
0
1
N
n
NN
n
n
zzz
zzz
zzz
Z
?
????
?
?
                       (3.12) 
The procedure for Prony?s method formulated in linear prediction theory to solve the linear 
equations (3.10) for Bi?s and zi?s can be summarized in the following steps  [10]:  
(1) Build a linear prediction model (LPM) for y(k) where the nth sample is a function of n past 
samples and unknown prediction coefficients {?1, ?2, ??n}  
)0()2()1()( 21 ynynyny n??? ?????? ?                                 (3.13) 
Apply this model for all successive samples as shown in (3.14). 
)1()1()()1( 21 ynynyny n??? ?????? ? 
?                                                                    (3.14) 
)1()3()2()1( 21 ????????? nNyNyNyNy n??? ? 
Rearranging (3.14) in matrix form gives (3.15) 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
????
??
??
??
)1(
)2(
)1(
)0(
)1()3()2(
)2()0()1(
)1()1()0(
)0()2()1(
2
1
Ny
ny
ny
ny
nNyNyNy
ynyny
ynyny
ynyny
n
?
?
?
????
?
?
?
?
?
?
                       (3.15) 
(2) Solve the LPM shown in (3.15) to find the predictor coefficients {?1, ?2, ??n}. 
The ??s in the linear equations shown in (3.15) can be solved directly if  N=2n or by the method 
of least squares if N>2n 
(3) Construct the LPM?s associated characteristic polynomial P(z) using the predictor 
coefficients as shown in (3.16)  
               nnnn zzzzP ??? ????? ?? 2211)(                                         (3.16) 
19 
 
(4) Determine the roots of the polynomial P(z). These roots are the required complex frequencies 
zi  [9]  [10]. 
(5) Substitute zi and y(k) in (3.10) and solve the resulting linear equations by using the least 
squares method to determine Bi. 
 (6) Substitute zi in (3.8) to find ?i. 
 
MATLAB Implementation of Prony Analysis 
The basic Prony?s method has been implemented in MATLAB the following way.    
 
 | 1 |  function [ vrMagn_n_PhaseRad, vrFreqRad_n_Decay, Dc_offset ] = PronySimple3( t, fn_data, 
timestep, numModes ) 
 | 2 |   
 | 3 |  N = length(fn_data);  
 | 4 |  init_f = fn_data(1); 
 | 5 |  fn_data = fn_data(2:N)-fn_data(1:N-1); 
 | 6 |  numDataPoints = length(fn_data); 
 | 7 |  
 | 8 |  % A : #cols = number of modes 
 | 9 |  % A : #rows = number of data points - number of modes 
 |10 |      
 |11 |  A = zeros( numDataPoints-numModes, numModes); 
 |12 |       
 |13 |  for r = 1:numDataPoints-numModes 
 |14 |      for c = 1:numModes 
 |15 |          A(r,numModes-c+1) = fn_data(c+r-1);         
 |16 |      end 
 |17 |  end 
 |18 |         
 |19 |  % B is always a column vector (only one column) 
 |20 |  % #rows = number of data points - number of modes 
 |21 |  B = zeros( numDataPoints-numModes, 1); 
 |22 |   
 |23 |  for r = 1:numDataPoints-numModes         
 |24 |      B(r,1) = -fn_data(numModes+r); 
 |25 |  end 
 |26 |     
 |27 |  X = A\B; 
 |28 |   
 |29 |  mu_ =roots([ 1 transpose(X)]); 
 |30 |   
 |31 |   
 |32 |  vrFreqRad_n_Decay = log(mu_)/timestep; 
 |33 |   
 |34 |  % M : #cols = number of modes 
 |35 |  % M : #rows = number of data points  
 |36 |  M = zeros(numDataPoints, numModes); 
 |37 |   
 |38 |  for r = 1:numDataPoints 
 |39 |      for c = 1:numModes 
 |40 |          M(r,c) = (mu_(c)-1)  *  (mu_(c) )^(r-1);  
 |41 |      end 
 |42 |  end 
 |43 |   
 |44 |  % F is always a column vector (only one column) 
 |45 |  % #rows = number of data points  
 |46 |   
20 
 
 |47 |  vrMagn_n_PhaseRad = M \ transpose(fn_data)  ; 
 |48 |   
 |49 |   
 |50 |  Dc_offset = init_f - sum(vrMagn_n_PhaseRad); 
 |51 |   
 |52 |   
 |53 |  return 
 |54 |   
 |55 | 
 
 
 
 
The original PA fits a sum of complex sinusoids to a given f(t). But, when a dc offset is present 
in the signal the algorithm needs some modifications. These modifications are done in the lines 
4, 5 and 50 as follows: 
 
 
 | 4 |  init_f = fn_data(1); 
 | 5 |  fn_data = fn_data(2:N)-fn_data(1:N-1); 
 ?  
 |50 |  Dc_offset = init_f - sum(vrMagn_n_PhaseRad); 
 |51 |   
 
 
 
The overall structure of the code can be be divided into mainly 3 parts:  
1: Solution of N>n linear equations in the form of AX=B to determine ??s in (3.15); 
2: Solving for the roots of an nth degree polynomial; and 
3: Solution of N>n linear equations in the form of AX=B to determine magnitudes and phase 
angles. 
In both the cases of solving the linear equations, it is necessary to use the SVD method ( svd() 
command in MATLAB) to accurately compute the solution for these equations. The main reason 
for that is because the matrices are essentially Vandermonde matrices and hence they are highly 
singular relative to a computer?s working precision. These changes are shown as follows: 
 
 |11 |  A = zeros( numDataPoints-numModes, numModes); 
 |12 |       
 |13 |  for r = 1:numDataPoints-numModes 
 |14 |      for c = 1:numModes 
 |15 |          A(r,numModes-c+1) = fn_data(c+r-1);         
 |16 |      end 
 |17 |  end 
 |18 |         
 |19 |  % B is always a column vector (only one column) 
21 
 
 |20 |  % #rows = number of data points - number of modes 
 |21 |  B = zeros( numDataPoints-numModes, 1); 
 |22 |   
 |23 |  for r = 1:numDataPoints-numModes         
 |24 |      B(r,1) = -fn_data(numModes+r); 
 |25 |  end 
 |26 |     
 |27 |  X = A\B; 
 |28 | 
 
 
 |34 |  % M : #cols = number of modes 
 |35 |  % M : #rows = number of data points  
 |36 |  M = zeros(numDataPoints, numModes); 
 |37 |   
 |38 |  for r = 1:numDataPoints 
 |39 |      for c = 1:numModes 
 |40 |          M(r,c) = (mu_(c)-1)  *  (mu_(c) )^(r-1);  
 |41 |      end 
 |42 |  end 
 |43 |   
 |44 |  % F is always a column vector (only one column) 
 |45 |  % #rows = number of data points  
 |46 |   
 |47 |  vrMagn_n_PhaseRad = M \ transpose(fn_data)  ; 
 
The only non-linear process is the finding of the roots of the polynomial to determine the zi.  
 
 
 
 
3.3 Example 
 
In this section, Prony?s method will be illustrated with the aid of a numerical example. 
Assume that it is required to fit to the record given in Table 3-1 a function of the form (3.17). 
The data samples are recorded at equally spaced intervals of 1 second beginning at 0 seconds. 
The actual form and the sinusoidal parameters of the function f(x) are shown in (3.18)   
 
xaxa eCeCxf 21 21)( ??                                                     (3.17)  
 
 
Table 3-1 A 6-Sample Record To Illustrate Prony?s Method 
 
x f(x) 
0 1.085950012 
1 -0.017837743 
2 -0.244420304 
3 -0.024960803 
4 0.051582021 
5 0.011677152 
22 
 
 
 
)16.0)23.0(2c o s (1.1)( 75.0 ?? ? xexf x ?                                        (3.18) 
 
 
 
It is more convenient to represent these computations if the quantities are converted using the 
identities given in (3.19a) & (3.19b) to that shown in (3.20). 
 
           11 ??ae                                                             (3.19a) 
               22 ??ae                                                            (3.19b) 
 
(3.18)  can now be written in the form of (3.20)   
 
                        xx CCxf 2211)( ?? ??                                                     (3.20) 
 
Substituting the data samples for x = 0,1,2,3,4,5 gives the following equations (3.21a) through 
(3.21f) 
 
                            085950012.1)0(21 ??? fCC                                        (3.21a) 
 
         017837743.0)1(2211 ???? fCC ??                                      (3.21b) 
 
244420304.0)2(222211 ???? fCC ??                                      (3.21c) 
 
024960803.0)3(322311 ???? fCC ??                                      (3.21d) 
 
051582021.0)4(422411 ??? fCC ??                                         (3.21e) 
 
011677152.0)5(522511 ??? fCC ??                                         (3.21f) 
 
Rewriting these equations in matrix form gives (3.22) 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)5(
)4(
)3(
)2(
)1(
)0(11
2
1
5
2
5
1
4
2
4
1
3
2
3
1
2
2
2
1
21
f
f
f
f
f
f
C
C
??
??
??
??
??
                                                            (3.22)
 
23 
 
 
Before C1 and C2 are determined using (3.22), ?1 and ?2 need to be determined first. Let us 
assume that ?1 and ?2 are the roots of a polynomial shown in (3.23)  
 
0212 ??? ????                                                              (3.23) 
or 
? ?? ? 021 ??? ????                                                          (3.24) 
 
Multiplying (3.21a) by ?2, (3.21b) by ?1 and adding the result to (3.21c) gives (3.25)    
       ? ? 22212 085950012.1)0( ??? ??? fCC  
? ? ? ? 1122111 017837743.01 ????? ???? fCC  
244420304.0)2(222211 ???? fCC ??  
 
? ? ? ? )0()1()2( 21222122121121 fffCC ?????????? ????????               (3.25) 
 
Because ?1 and ?2 are roots of the polynomial, the left hand side of the equations computes to 
zero. Repeating the same process using (3.21b), (3.21c) and (3.21d) gives (3.26b). Repeating 
again gives (3.26c) and (3.26d).     
 
? ? ? ? 0085950012.1017837743.0244420304.0 21 ????? ??                          (3.26a) 
? ? ? ? 0017837743.0244420304.0024960803.0 21 ?????? ??                         (3.26b) 
? ? ? ? 0244420304.0024960803.0051582021.0 21 ????? ??                         (3.26c) 
? ? ? ? 0024960803.0051582021.0011677152.0 21 ???? ??                        (3.26d) 
 
Rearranging these equations in matrix form yields (3.27) 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
????
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
??
?
011677152.0
051582021.0
024960803.0
244420304.0
024960803.0051582021.0
244420304.0024960803.0
017837743.0244420304.0
085950012.1017837743.0
2
1
?
?      (3.27) 
 
This is in the form of AX=B. The solution, X, can be found using the least squares method as 
shown in (3.28). However, because of the nature of the numerical values involved it was 
24 
 
observed that using singular value decomposition (SVD) resulted in considerably better results. 
The solution vector, X, found in Matlab using SVD (3.29) is shown in (3.29).   
 
? ? BAAAX TT 1??                                                              (3.28) 
 
 
X=A\B                                                                   (3.29) 
 
Solving (3.29) gives ?1 and  ?2  as shown in (3.29b)  
 
 
????????????? 0148430 .2 2 3 1 3 0 1 6 4965310 .1 1 8 4 0 6 4 5-21??                                        (3.29b) 
 
 
Substituting ?1 and ?2 in (3.23) and determining its roots, yields ?1 and ?2 values as:  
 
?1 = 0.05920322748265 + 0.46864180138360i and 
?2 = 0.05920322748265 - 0.46864180138360i. 
 
The exponential parameters of (3.18), a1 and a2, are then determined as:  
 
 
a1 = log(?1) = -0.75 + 1.44513262065131i  and 
a2 = log(?2) = -0.75 - 1.44513262065131i. 
 
Substituting ?1 and ?2 in (3.22) gives (3.30). 
 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
0.01167
0.05158
0.02496-
0.24442-
0.01783-
1.08595
) 0.019026i - (0. 01 38 23) 0.019026i + (0. 01 38 23
0.023985i )  0.043628 (0.023985i ) - 0.043628 (
0.097997i )  (-0. 03 88 00 0.097997i ) - (-0. 03 88 00
0.055490i ) - (-0. 21 61 200.055490i ) + (-0. 21 61 20
0.468641i ) - 0.059203 (0.468641i ) + 0.059203 (
11
2
1
C
C
                  (3.30)
 
 
The overly determined least square problem, (3.30), is solved in Matlab using SVD to obtain C1 
and C2 as: 
 
C1 = (0.542975 + 0.087625i) and 
C2 = (0.542975 - 0.087625i). 
25 
 
 
Substituting a1, a2, C1 and C2 in (3.18) gives (3.31). 
 
?? xexf )1 . 4 4 5 1 i + - 0 . 7 5 0 0(0 . 0 8 7 6 2 5 i ) + ( 0 . 5 4 2 9 7 5)(  
                                 xe )1 . 4 4 5 1 i - - 0 . 7 5 0 0(0 . 0 8 7 6 2 5 i ) - 5. . ( 0 . 5 4 2 9 7.    (3.31) 
 
Using the complex exponential identities for cosine, (3.31) is converted to its sinusoidal form as 
shown in (3.32). The parameters of this sinusoid are equal to those of the specified function in 
(3.17). Thus, we have seen that by using Prony?s method on a record consisting of six samples, 
we have been able to determine the original parameters of the function used in constructing the 
record. 
 
)16.0)23.0(2c o s (1.1)( 75.0 ?? ? xexf x ?                                    (3.22) 
 
 
 
3.4 Conclusion 
 
The mathematical analysis of Prony's method and its implementation in MATLAB are 
reviewed in this chapter. Prony?s method, an exponential parameter estimation technique that is 
used widely in communications as a signal analysis tool, has been chosen over the FFT in this 
work because of it many advantages. It is ideally suited to extract modal information from 
transient stability programs? outputs and STG shaft signal measured data. 
While Prony?s method faithfully fits a specified number of modal components to a given 
signal, it is necessary to reorder the output resultant modes of Prony's method to distinguish the 
signal modes from the noise modes. The next chapter presents the details of a reordering 
algorithm developed during the course of this work. The practical issues such as choosing the 
proper parameters of Prony?s method to obtain better results are examined in Chapter 5. 
 
 
 
 
 
 
 
CHAPTER 4 
REORDERING ALGORITHM 
 
In this work, Prony?s analysis is used to perform modal decomposition of a given 
function. When the noise in the signal is low or negligible, the specified modal order can be set 
equal to the true order. However, usually there is noise or nonlinearities in the signal. In those 
cases, the modal order should be specified many times larger than the true internal modal order.  
Prony?s analysis is then performed to extract the modes in the given function. The resulting 
output of Prony?s analysis is in no particular order. Therefore, it is hard to distinguish the true 
modes from noise. A separate sorting algorithm is required for that purpose. One could employ 
any kind of sorting algorithm based on the constituent parameters: (1) magnitude, (2) frequency, 
or (3) exponential decay. However, these have their own limitations. Therefore, for torsional 
modes, it was observed during the course of this work, that a least-squares-error norm based 
sorting algorithm gives much better results. A weighting factor is computed for each constituent 
sinusoidal mode. These weighting factors are sorted in descending order to give the final list of 
modes in which the torsional modes are accumulated at the top of the list. As with Prony?s 
analysis, this reordering algorithm was also coded in MATLAB. 
 
 
 
 
 
 
 
27 
 
4.1 Proposed Reordering Algorithm 
Let fT(k) be the given linearly sampled signal. Prony?s method can be used to extract the 
constituent modal components by fitting the given signal to a sum of n complex exponentials as 
shown in (4.1) by computing the complex values Bi and zi. 
 
knnkkT zBzBzBkf ???? 2211)(                                             (4.1) 
 
The output of the Prony algorithm is a list of n complex value pairs where each pair represents an 
exponentially decaying sinusoid. The complex conjugate pairs can be converted to parameters of 
sets of exponentially damped sinusoids as shown in (4.2). 
 
? ? ? ?? ? ? ?iiiiiiii AzBzB ??? ,,,,,, ** ? 
or 
? ?
? ?
? ??
?
?
?
?
?
?
?
?
?
?
?
NNNNA
A
A
???
???
???
,,,
,,,
,,,
2222
1111
?                                                                   (4.2) 
 
Let the ith mode be reconstructed to give fi(k) as shown in (4.3). 
 
)c o s (:)( iitii teAkf i ??? ?                                                          (4.3) 
 
The given signal fT(k) can thus be represented as a sum of individual exponentially damped 
sinusoids and the error as shown in (4.4) and (4.5).  
28 
 
 
ekfkfkfkfkf NT ?????? )()()()()( 321 ?                                        (4.4) 
ekfkf RT ?? )()(                                                       (4.5) 
? ?iiiii Akf ??? ,,,)( ?                                                     (4.6) 
 
In (4.6), fi(k)  is the signal built using the ith set of parameters.  Let us compute the error norm for 
the first sinusoid f1(k). This indicates how much error is generated if this mode were not present 
in the reconstructed waveform fR(k). 
? ?)()()()( 321 kfkfkfkfe NT ????? ?                                      (4.7)
  
ekfe
or
kfkfe
N
i
i
T i
??
?? ?
?
?
)(
)()(
11
1
1
1
                                                        
(4.8)
 
ekfe
or
kfkfe
jj
N
i
ji
Tj i
??
?? ?
?
?
)(
)()(
1
                                                         
(4.9)
 
 
In this way the error norms for each sinusoid ? ?)(,),(),( 21 kfkfkf N? are computed resulting in the list 
L given in (4.10). 
L:? ?Neee ,,, 21 ?                                                             (4.10) 
29 
 
 
Sorting the list L in descending order according to magnitude yields the required dominant 
modes list. The entire algorithm can be summarized in a flowchart as shown in Figure 4-1. 
 
 
 
 
 
 
 
Figure 4-1. Flowchart for Reordering Algorithm 
4.2 Example 
This method can be illustrated with the aid of the following example. The sampled 
waveform in this example is constructed out of four modal components (three exponentially 
decaying sinusoid modes and one exponentially decaying dc mode) and noise (w(t)) as shown in 
(4.11).                    
)(17)5.22c o s (5.9)1802c o s (1.0)152c o s (5.3)( 05.002.08.01.0 twetetetetf tttt ????? ???? ???   (4.11) 
The MATLAB code for constructing this signal is listed in Figure 4-2.  
 
 
 
 
Figure 4-2. MATLAB Code for Creating the Sample Waveform 
Prony?s Method 
fT(k) 
Bi ,zi 
 
Reordering algorithm 
Bj ,zj 
 
clearall;clc; 
t = 0:0.0001:1; 
noise_ = (rand(size(f))-0.5) * 1e-5; 
f =  3.5 * exp(-0.1 *t) .* cos(2*pi* 15*t)  +   
0.1 * exp(-0.8 *t) .* cos(2*pi*180*t)   + ... 
     9.5 * exp(-0.02*t) .* cos(2*pi*2.5*t)  +  
17   * exp(-0.05*t)  + noise_ ; 
plot( t, f ); 
 
30 
 
A plot of f(t) vs time is shown in Figure 4-3. 
 
Figure 4-3. Plot of the Sample Waveform. 
Prony?s method can now be applied on this signal to extract the underlying modal components. 
Note that while applying the method, it is required to specify the method the number of modes in 
the signal. Generally speaking, whenever there is noise or any kind of non-linearities in the 
signal, a modal order much higher than the true estimated modal order needs to be specified. For 
example, in this case, the modal order is specified as 200 instead of the true modal order, four. 
Prony?s method fits these extraneous modes to noise in the signal. Experience, or trial and error, 
is usually used to estimate this modal order. 
The output, a list of modes, of the MATLAB function implementing Prony?s method is given 
in Table 4-1. Note that these modes are in no particular order. It is difficult to distinguish the real 
signal modes from the noise modes. Because the signals that are of significance in the current 
work are slowly decaying lower frequency power system signals, which span the range of 1 Hz 
to a few 100 Hz, the described method can be applied to separate the signal modes from the 
noise modes. The MATLAB code shown in Figure 4-4 implements this procedure. 
 
31 
 
 
 
Table 4-1. Output of Prony?s Method 
 
Mode    Magnitude  decay f(hz)    Phase  
1 17.0672167 -0.050 0.000 0.00 
2 2.78E-08 -722.469 0.000 -3.14 
3 3.248E-07 -6.565 497.942 -0.94 
4 6.276E-07 -3.990 491.587 -2.88 
5 8.207E-07 -14.330 485.991 -2.96 
6 6.989E-07 -9.757 482.354 2.10 
7 6.564E-07 -5.200 477.294 -0.45 
8 3.9346E-06 -18.105 471.655 3.11 
9 2.7019E-06 -13.365 467.724 1.24 
13 1.02827E-05 -17.819 445.327 -2.78 
14 9.9657E-06 -14.832 443.801 0.73 
 ?   ?   ?   ?   ?  
94 0.0000015 -7.454 267.251 -0.69 
95 2.0571E-06 -5.710 246.345 2.58 
96 1.0887E-06 -7.345 254.272 -2.59 
97 1.7877E-06 -8.856 242.250 2.07 
98 1.2883E-06 -15.795 259.649 -1.94 
99 1.7721E-06 -25.196 262.943 3.14 
100 1.63115E-05 -121.743 247.006 -2.58 
101 0.000032574 -362.771 283.476 -0.12 
 
 
 
 
 
 
 
 
Figure 4-4. MATLAB Code for Generating the LSE Norms 
 
The MATLAB function create_vrLSE() shown in Figure 4-4 is used to reorder the list of modes 
in Table 4-1. Shown in Table 4-2 is the output of this MATLAB function (the reordered list of 
modes). 
function [ vrLSE ] = create_vrLSE(   t, 
vrMagnFinal, vrDecayFinal, vrFreqHzFinal, 
vrPhaseRadFinal, fdiff) 
vrLSE = zeros(size(vrMagnFinal));  
for k=1:length(vrMagnFinal)   
fndata_estk = RebuildFn(t, 0, vrMagnFinal, 
vrDecayFinal, vrFreqHzFinal, 
vrPhaseRadFinal, k); 
lse_k = norm(fdiff+fndata_estk); 
vrLSE(k) = lse_k;  
% LSE without mode k 
end 
return 
 
32 
 
 
 
Table 4-2. Reordered List of Prony?s Method Output Modes 
 
Mode #  LSE  Magnitude decay f (hz) phi 
1 5.27E+02 1.71E+01 -0.05 0 0 
2 2.10E+02 9.50E+00 -0.02 2.5 0 
3 7.45E+01 3.50E+00 -0.1 15 0 
4 1.58E+00 1.00E-01 -0.8 180 0.1 
5 2.65E-04 1.00E-05 -14.832 443.8 0.7 
6 2.65E-04 1.03E-05 -17.819 445.33 -2.8 
7 2.64E-04 9.30E-06 -16.161 361.72 -2.7 
8 2.64E-04 3.26E-05 -362.771 283.48 -0.1 
9 2.64E-04 8.10E-06 -17.014 358.89 1.6 
 ?   ?   ?   ?   ?   ?  
97 2.62E-04 3.00E-07 -6.072 294.37 0.6 
98 2.62E-04 3.00E-07 -6.521 97.48 1.9 
99 2.62E-04 3.00E-07 -6.565 497.94 -0.9 
100 2.62E-04 1.00E-07 -2.333 326.05 -1 
101 2.62E-04 0.00E+00 -722.469 0 -3.1 
 
 
As seen in Table 4-2, the power signal modes are all sorted at the top because of the specific 
error norm based reordering method. 
4.3 Conclusions  
An algorithm to reorder the output modes of Prony's method, which was developed 
during the course of this work, is presented in this chapter. This reordering is performed to help 
distinguish the signal modes from noise modes in the given signal. A MATLAB implementation 
is presented along with an example that illustrates the method.  In this work, this methodology 
has been applied every time after modal decomposition is performed using Prony analysis as 
shown in Chapter 5, 6, 7 and 8. However, this algorithm was not applied to real measurements in 
Chapter 9 because of reasons discussed in that chapter. 
33 
 
While the issues regarding the output of Prony's method have been dealt with in this 
chapter, the next chapter deals with the practical issues of choosing the input parameters of 
Prony's method that aid in getting better results. 
 
 
 
 
 
 
 
CHAPTER 5 
EFFECT OF PRONY PARAMETERS 
 
The main focus of the previous chapter was the reordering algorithm for the output modes of 
Prony?s analysis. The focus of this chapter will be the practical issues of applying Prony 
analysis. These discussions will focus primarily on the effects of varying Prony analysis 
parameters on the output of Prony analysis. Basically there are three parameters: 
? Modal Order ( numModes ), 
? Sampling frequency ( nSkipSample ), and 
? Window parameters - length & position.  
These parameters need to be chosen with extreme consideration just like those of other 
identification methods. The user would set these parameters before running the Prony algorithm. 
A Matlab GUI has been developed to facilitate changing and setting these parameters in an easy 
and intuitive way. This careful selection of parameters is especially needed when the data is not 
explicitly linear or contains significant noise. 
 
5.1 Prony Driver GUI Program Description 
This section contains a brief explanation of the Matlab tool that has been developed to run 
Prony?s method. This tool has a GUI shown in Figure 5-1 built using the Matlab GUIDE 
program. It facilitates setting and modification of Prony parameters and data window properties. 
35 
 
 
 
 
 
Figure 5-1. GUI Built in Matlab to Allow User Selection of Data Window Properties and Prony 
Parameters 
The user would enter ?PronyDriver? at the ?>>? prompt in Matlab command window. A file 
open dialog pops up and allows the user to select a data file. This data file contains two columns 
which represent time and signal data. Once the program reads the file data, it displays the signal 
data in a plot to allow the user to select data window properties such as tmin and tmax.  
36 
 
The GUI will then pop up and the user can set the required Prony parameters and data window 
properties. The program will then run Prony?s method and display its output in the Matlab 
command window. 
 
5.2 Effect of Modal Order 
The process can be explained with the aid of an example. Shown in Figure 5-2 is the plot of a 
signal built in Matlab. This signal has three sinusoidal components and noise, i.e., a linear part 
and a non-linear part. 
)()5.112c o s (15.0)5.52c o s (5.0)5.12c o s (5.2)( 01.02.02.0 twtetetetf ttt ??????? ??? ???      (5.1) 
In (5.1) w(t) represents noise. 
 
Figure 5-2 Sample Waveform Consisting Of Three Exponentially Decaying Sinusoids and Noise. 
 
Prony?s method (PM) is then applied on this signal for modal decomposition. An input parameter 
to the PM, n, determines the modal order of the constructed model. PM will try to fit this model, 
a waveform consisting of a linear combination of n sinusoids, to the given signal. Shown in 
Figure 5-3 is the output of PM. Shown in Figure 5-3a are both the given signal and the 
37 
 
constructed model. The error between these two signals is shown in Figure 5-3b. Shown in 
Figure 5-3c is the stem plot, which represents the magnitudes of LSE of each individual 
component in the model. LSE of a sinusoidal component is a kind of weighting factor which 
represents its 2-norm. Because the signal consists of three sinusoidal components, the stem plot 
should display only three non-zero stems with the rest being zero magnitudes. However, because 
of noise in the signal, PM is unable to fit the model to the given data perfectly resulting in error 
shown in Figure 5-3b. This results in non-zero LSE magnitudes of noise modes. 
 
Figure 5-3 Comparison of Actual Signal and PM?s Output: Stem Plot Showing the LSE 
Magnitudes of Constituent Sinusoids. 
38 
 
Repeating this procedure with different values of numModes while keeping all other parameters 
constant, gives Figures 5-4a through 5-4g. In these figures there are 7 sets of PM results with 
numModes = {10, 20, 30, 60, 70, 80, 90}.  
 
Figure 5-4a Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 10 
 
39 
 
 
Figure 5-4b Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 20  
 
 
 
 
Figure 5-4c Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 30 
 
40 
 
 
Figure 5-4d Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 60 
 
 
 
Figure 5-4e Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 70 
41 
 
 
Figure 5-4f Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 80 
 
 
 
Figure 5-4g Comparison of Actual Signal and PM?s Output Along with Stem Plot of the 
Constituent Modes? LSE Magnitudes with numModes Parameter = 90  
 
42 
 
From these figures, it is evident that as numModes increases the model waveform increasingly 
matches the given signal which results in decreasing error magnitude. In other words, the 2-norm 
of the error signal decreases as numModes increases. Also, the true signal modes? LSEs in the 
stem plot becomes increasingly obvious, i.e., they stand out from other noise stems more 
distinctly because as numModes increases the LSE magnitudes of the noise modes will gradually 
diminish. Comparing the LSE magnitudes of noise modes in set (1) and set (8), the magnitude 
has decreased from 75 to less than 0.1. This property has been used in creating an LSE based 
reordering algorithm, discussed in Chapter 4, which reorders the modes based on LSE 
magnitudes. Thereafter, another algorithm is used to separate the signal modes from the noise 
modes.  
Increasing the numModes usually results in decreased error. However, there is a 
particular characteristic to this decrease.  This characteristic is captured and illustrated in Figure 
5-5. It can be observed that the error initially decreases only slightly. However, after a threshold 
it decreases exponentially. Note that the y-axis is logarithmic. 
 
Figure 5-5 Plot of LS Error as a Function of numModes 
 
43 
 
Next, the same experiment is repeated with the magnitude of noise in the signal as the varying 
parameter. The relationship between these LSEs and numModes for various noise magnitudes 
are plotted as shown in Figure 5-6. It can be observed from this figure that as the noise in the 
signal increases, it takes more numModes for the LSE to fall below a certain specified level.  
 
 
Figure 5-6 Total LSE vs numModes with Noise as Varying Parameter for a 3-Sinusoid Signal 
 
Next, the previous experiment is repeated with a different modal composition in the given signal. 
This time a signal composed of only one sinusoid instead of three sinusoids is used. The 
corresponding relationship between LSE and numModes is shown in Figure 5-7. 
 
44 
 
 
Figure 5-7 Total LSE vs numModes with Noise as a Varying Parameter for a 1-Sinusoid Signal 
 
5.3 Effect of Sliding Window Parameters on Prony?s Method?s Output 
As mentioned previously, the output of the Prony analysis for a given signal depends on 
three parameters that are chosen by the user. These are the selected modal order, the sampling 
rate and the data window?s size and position. The effects of the window parameters ? length and 
position ? on the output of Prony analysis will be discussed in this section. For convenience sake, 
a new parameter, named nSkipSample, will be used instead of sampling rate. The nSkipSample 
parameter downsamples the given record by sampling only every nth point. The window size and 
position are specified by two parameters tmin, which specifies the beginning point of the data 
window, and tmax, which specifies the ending point of the data window. 
Consider the same signal that was used in the previous chapter as shown in (5.2).  
)()5.22c o s (5.9)1802c o s (1.0)152c o s (5.3)( 02.08.01.0 twtetetetf ttt ??????? ??? ???    (5.2)  
This data record is fed to Prony?s method to extract its constituent modes using the following 
input parameters:  
45 
 
Modal Order (numModes) : 50 
nSkipSample : 4 
Window size and position: tmin = 0 and tmax = 2.0 
Shown in Figure 5-8 is the output of the Prony?s method. As seen in this Figure, the Prony model 
matches the actual given data very closely. 
 
 
Figure 5-8 Comparison of Given and Estimated Waveforms 
 
The error norm, LSE, for this window is 2.41E-5. Next, the same experiment is repeated while 
sliding the window 0.4 seconds every time. The results of these tests are shown in Figures 5-9a 
through 5-9c and Tables 5-1 through 5-4. Here, as the very small LSE numbers indicate, the 
model matches the given signal very closely.  
 
 
46 
 
Table 5-1. The LSE Values for Successive Windows 
Window #  tmin tmax LSE 
1 0 2 2.41E-05 
2 0.4 2.4 7.52E-06 
3 0.8 2.8 1.91E-05 
4 1.2 3.2 3.67E-05 
5 1.6 3.6 4.39E-05 
6 2 4 7.92E-06 
 
Table 5-2 Effect of Sliding Window on Frequency of Mode1 
Window # Magnitude Decay Frequency (Hz) Phase 
1 9.4999807 -0.02 2.5 0.0015708 
2 9.4242842 -0.02 2.5 0.0015708 
3 9.3491911 -0.02 2.5 0.0015708 
4 9.2746952 -0.0199999 2.5 0.0015708 
5 9.2007935 -0.0199999 2.5 0.0015708 
6 9.1274815 -0.02 2.5 0.0015708 
 
Table 5-3 Effect of Sliding Window on Frequency of Mode 2 
Window # Magnitude Decay Frequency (Hz) Phase 
1 3.499965 -0.1 15 0.0094248 
2 3.3627294 -0.1 15 0.0094248 
3 3.2308749 -0.1 15 0.0094248 
4 3.1041905 -0.1 15 0.0094248 
5 2.9824734 -0.1 15 0.0094248 
6 2.865529 -0.1 15 0.0094248 
 
Table 5-4 Effect of Sliding Window on Frequency of Mode 3 
Window # Magnitude Decay Frequency (Hz) Phase 
1 0.099992 -0.8 180 0.1130973 
2 0.0726091 -0.8 180 0.1130973 
3 0.052725 -0.8 180 0.1130973 
4 0.0382862 -0.8 180 0.1130973 
5 0.0278015 -0.8 180 0.1130973 
6 0.020188 -0.8 180 0.1130973 
47 
 
 
 
Figure 5-9a Comparison of Estimated and Actual Magnitudes of Mode 1 (2.5 Hz) mode 
 
 
Figure 5-9b Comparison of Estimated and Actual Magnitudes of Mode 2 (15 Hz) mode 
 
 
Figure 5-9c Comparison of Estimated and Actual Magnitudes of Mode 3 (180 Hz) mode 
48 
 
 
5.4 Effect of varying sampling rate parameter on Prony?s Method?s Output 
              In this section, the effects of varying the sampling-rate parameter on the PM?s output, 
and consequently its LSE, are observed. As mentioned previously, a new parameter, named 
nSkipSample, will be used instead of sampling rate. The nSkipSample parameter downsamples 
the given record by sampling only every nth point.  Consider the same data set that was used in 
previous section. It is given in (5.3). 
)()5.22c o s (5.9)1802c o s (1.0)152c o s (5.3)( 02.08.01.0 twtetetetf ttt ??????? ??? ???      (5.3) 
Because of noise in the signal, the numModes parameter was set to 100 and was kept constant 
throughout all five test cases i.e., nSkipSample = {1,2,3,4,5}. This data record is originally 
sampled at 10 kHz and contains two seconds of data. Therefore the maximum number of modes 
that we could specify to PM as a parameter is 20,000/2 = 10,000. However, because the highest 
frequency is 180 Hz, we could specify a much lower numModes parameter, usually less than 
500. As mentioned previously, experience is often used in this choice. As can be seen in Figure 
5-10a through Figure 5-10e,  as the nSkipSample parameter increases, the PM fits the model to 
given data sample more and more completely, i.e., the PM model increasingly matches the given 
data sample. Consequently, the LSE decreases as seen in Table 5-5.  
 
49 
 
 
Figure 5-10a Comparison of Actual Signal and PM?s Output Along With Stem Plot of the 
Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 1 
 
 
 
 
 
Figure 5-10b Comparison of Actual Signal and PM?s Output Along With Stem Plot of the 
Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 2 
50 
 
 
Figure 5-10c Comparison of Actual Signal and PM?s Output Along With Stem Plot of the 
Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 3 
 
 
 
 
Figure 5-10d Comparison of Actual Signal and PM?s Output Along With Stem Plot of the 
Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 4 
51 
 
 
Figure 5-10e Comparison of Actual Signal and PM?s Output Along With Stem Plot of the 
Constituent Modes? LSE Magnitudes with nSkipSample Parameter = 5 
 
 
Table 5-5  Sampling Rate (nSkipSample) vs LSE - A 
nSkipSample LSE-error 
1 1317.27 
2 21.1131 
3 0.52178 
4 0.06744 
5 0.05099 
 
It is to be noted that as the nSkipSample parameter rises the resolution decreases, and hence the 
PM might miss detecting some higher frequencies. In this case, because the highest frequency is 
180 Hz and the lowest sampling rate (case nSkipSample = 5) is 2 kHz, there is no such problem. 
52 
 
From Table 5-6 it can be inferred that although the total LSE is decreasing quickly, the 
individual parameters, i.e., magnitudes, frequencies, etc, converge to their correct values slowly. 
It requires nSkipSample of greater than or equal to four for the magnitude of the frequency 2.5 
Hz mode to converge to its correct value, 9.5. 
Table 5-6 Sampling rate (nSkipSample) vs LSE - B 
  Mode #  Magnitude ExpDecay Frequency(Hz) Phase(rad) 
nSkipSample:   1  1 19.7376133 -1.882675 0 0 
  2 3.9084642 -0.2312709 14.9710909 0.1716462 
  3 1.1300583 -144.9632392 0 3.1415927 
  4 0.1001143 -0.8003881 179.9999787 0.1125296 
nSkipSample:   2  1 9.695243 -0.0408697 2.5001401 -0.0012846 
  2 3.4998881 -0.0999648 14.9999748 0.0096094 
  3 0.1000002 -0.8001528 180.0000277 0.113005 
nSkipSample:   3  1 9.4940861 -0.0193665 2.4999749 0.0017882 
  2 3.4998786 -0.0999756 15.0000111 0.0093589 
  3 0.0999947 -0.8000499 180.0000022 0.1130859 
nSkipSample:   4  1 9.5007977 -0.0200874 2.5000007 0.0015574 
  2 3.4999868 -0.1000066 15.0000012 0.0094178 
  3 0.0999847 -0.7998804 180.000008 0.1130771 
nSkipSample:   5  1 9.5007017 -0.0200781 2.5000061 0.0015255 
  2 3.4999951 -0.1000082 14.9999978 0.0094376 
  3 0.0999979 -0.8001168 179.999992 0.1131259 
 
In the previous example, the numModes parameter was kept constant, while the sampling 
rate varied. Now the sampling rate is kept constant (nSkipSample = 1) and numModes is 
increased while keeping the window size and position constant. Next, the LSE for each 
numMode is recorded. Plotting the corresponding LSE vs. numMode parameter gives us data for 
one curve. Repeating the same example with varying sampling rate   (nSkipSample = 2, 3) would 
give the family of curves shown in Figure 5-11. 
53 
 
 
 
Figure 5-11 Effect of Sampling Rate on the Relationship Between numModes and LSE 
 
It can be observed from these families of curves that while it took about 150 
modes/numModes for obtaining sufficiently less final LSE with the SkipSample2 parameter, it 
took only about 80 numModes for obtaining the final LSE with the SkipSample4 parameter. One 
inference: in some cases increasing the sampling rate (thereby the number of samples) does not 
necessarily give us better results. It will take more numModes, and hence, become more 
computationally expensive and time consuming. In other words, the sampling rate can be 
54 
 
decreased to a certain extent, provided it is much higher than the highest frequency. 
Additionally, it can also be inferred that, as long as the sampling rate is much higher than the 
highest signal frequency, we can afford to reduce the sampling frequency to obtain the same 
result with the same LSE, but much quicker and less computationally expensive. 
  
5.5 Guidelines 
Based on the observations working with Prony?s method, some general guidelines have 
been developed for setting these parameters. These mostly apply for extracting power system 
torsional modal information. The guidelines for real signal measurements are slightly different 
from those of transient stability program (TSP) output signals, because of the presence of 
measurement noise and other nonlinearities in real signal measurements. 
 
Setting the modal order parameter (numModes) : 
As mentioned earlier, the real signal measurements contain measurement noise and other 
nonlinearities. To compensate for this, it is required to choose a larger numModes as an input 
parameter for Prony's method. It is difficult to choose the correct value in a single step [10]. It is 
usually an iterative process where the user sets the numModes parameter and executes Prony's 
method. The output modes are used to reconstruct the signal and if it is not a good fit, a higher 
numModes parameter is set and the process is repeated. However, it was observed that too high a 
modal order, while resulting in a good fit, could result in incorrect modal results probably 
because the matrices used in Prony?s method?s becoming increasingly singular. Therefore, 
numModes should be increased gradually, until a good fit for TSP output signals or a reasonable 
enough fit for real measurement signals is obtained. 
55 
 
Choosing the window location and size parameters: 
When choosing the window position for a given TSP output record, it is advisable to choose later 
because these torsional modes do not decay quickly in comparison to other short duration modes, 
such as electrical stator transients. For real measurement records, it is advisable to choose earlier 
because as the length of the record increases, the probability of non-linear phenomena -  
switching events or other influences because of corrective measures, such as governors and 
exciters - distorting the LTI nature of the overall signal increases. Finally, the window size 
should include at least two cycles of the lowest modal frequency in most cases. 
 
 Setting the sampling rate parameter: 
Increasing this parameter is another way of modifying numModes and the total record length. 
The sampling rate should be such that each mode in the signal would at least have four samples 
per cycle similar to the Nyquist criteria. Increasing the sampling rate parameter thereafter does 
result in a better fit. However, increasing it too much does not necessarily increase fit quality, 
because the matrices involved in Prony?s analysis tend to become singular and give incorrect 
results. 
 
5.6 Conclusions 
The process of selecting the various parameters for Prony?s method is complex. It 
requires careful consideration as there are many parameters that affect directly the results and 
their accuracy. Unlike Fourier analysis, each data record must be treated individually. While 
there are some empirical rules, as discussed in the previous section, the user does need some 
training/practice before obtaining good results. The most important parameter is numModes, 
which is used to fit not only the signal modes, but also any noise or non-linearities in the signal. 
 
 
 
 
 
 
 
 
 
 
CHAPTER 6  
 
EIGEN ANALYSIS OF STEAM TURBINE GENERATOR SHAFT SYSTEMS 
 
The main aim of this work, as mentioned earlier, is to estimate the shaft spring constants 
so that they could be used to build better models for SSR analysis. These shaft spring constants 
are estimated using rotor mass values and measured modal frequency values obtained from either 
the rotor speed signals or the power signal p(t). The accuracy of the estimated shaft spring 
constants is therefore directly dependent on the accuracy of the measured frequencies. This 
accuracy is in turn dependent on the signal-to-noise ratio of the frequencies in the signal. In other 
words, it is highly preferable to have a signal where the magnitudes of modal frequencies are 
large enough that they could be measured reliably. One of the major factors that affects the 
magnitudes of modal components in a signal is the eigen properties of the concerned shaft 
system.  In this and the next chapter, these aspects, eigen properties of the STG shaft systems, 
will be investigated to determine their effects on this work.  In order to apply eigen analysis, the 
isolated lumped-spring-mass mechanical model of the STG turbine will be treated as an LTI 
system. There are certain unique properties of STG shaft systems because of the way the rotor 
masses and the interconnecting shafts of the lumped-spring-mass model are connected. These 
properties will be examined with the aid of a few shaft systems whose data is provided in the 
literature [4][16]. 
 
 
57 
 
 
6.1 State-Space Theory Review [17] 
One of the main reasons eigen analysis is used here is that it makes it easier to analyze 
state systems of different sizes in a similar unified manner. Using eigen analysis allows the use 
of linear algebra principles and well-established computer routines to easily analyze the system 
characteristic features such as eigen values, eigen vectors and participation factors. A brief 
review of the relevant state space theory used in this work is provided in the following section. 
Let the dynamics of an nth order non-linear system be described by the set of first order 
differential equations given by (6.1) and its outputs be described by (6.2). 
? ?)(),(),(),()(),()( 2121 tutututxtxtxfdt tdx mnii ???                                             (6.1) 
i = 1, 2, ? n  
? ?)(),(),(),()(),()( 2121 tutututxtxtxgty mnii ???                                      (6.2) 
In (6.1) and (6.2) fi and gi are non-linear functions of variables, also known as system state 
variables x1(t), x2(t) ? xn(t),  and inputs u1(t), u2(t) ? um(t). 
The states, xi, are a minimal set of linearly independent system variables, either real, 
physical (such as capacitor voltage and inductor current) or abstract mathematical variables 
associated with the system equations. Knowledge of these state variables at t=t0 along with the 
inputs for t>t0 is sufficient to completely determine the system behavior and all other system 
variables. The selection of state variables is not unique. Any selected set would yield the same 
information about the state except in a different coordinate system in the dynamic system's n-
dimensional euclidean state space. If the number of state variables is chosen greater than the true 
58 
 
order of the system, then not all of them will be linearly independent as shown in the following 
example.  
There are at least two sets of state variables that can be used to represent the system 
shown in Figure 6-1: (1) v1, v2, i1, i2 and (2) v1, v2, v3, i1, i2. Only the first set one is the correct 
choice because it is a minimal set, whereas set (2) is not. Set (2) is not a minimal set because any 
one of the three voltage variables (v1, v2, v3,) can be represented as a linear combination of the 
other variables. 
 
 
 
 
 
 
 
 
 
 
 
Figure 6-1. Example of a Four State System Consisting of Inductors and Capacitors 
In the case of a linear time-invariant system, the non-linear functions fi and gi can be replaced by 
constant coefficient matrices and these equations can be written in vector-matrix form as shown 
in (6.3a) and (6.3b). 
 
)()()( tButAxdt tdx ??                                                        (6.3a) 
C1 C2 C3 
L1 L2 v1 v2 v3 
i1 i2 
59 
 
)()()( tDutCxty ??                                                          (6.3b) 
 
The matrix A is commonly referred to as the system?s state matrix. B, C and D are all constant 
coefficient matrices. Using the traditional transfer function approach, the transfer function of this 
system can be derived after taking the Laplace transform of the equations and setting the initial 
conditions to 0. The resulting transfer function is shown in (6.4). 
 
? ?AsI DAsIBAsIadjCsG ? ???? )()(
                                                     (6.4) 
Setting the denominator of G(s), also referred to as the characteristic polynomial, to 0 gives the 
characteristic equation (6.5).  
 
0??AsI                                                                   (6.5) 
 
6.2 Eigen values 
The roots of the characteristic equation shown in (6.5) are called the eigen values, ?i, of the 
system. When (6.3a) and (6.3b) are used to describe the dynamics of an isolated linear time-
invariant TG shaft system, these eigen values correspond to its natural frequencies or torsional 
mode frequencies.  
 
6.3 Right eigen vectors 
60 
 
A right eigen vector of matrix A associated with ?i is any non-zero vector that satisfies equation 
(6.6). 
iiiA ??? ?                                                                       (6.6) 
i = 1, 2, ? n  
Rearranging terms in (6.6) gives (6.7). 
? ? 0?? iiIA ??                                                                (6.7) 
Because equation (6.7) is a homogenous equation, the eigen vector can only be determined  
within a scalar multiplier.  Generally, this vector is scaled such that it is of unit length. In this 
current work, for displaying purposes only, these vectors will be occasionally scaled such that 
the absolute value of the largest element in the vector is 1. In the context of TG rotor system 
analysis, these right eigen vectors, also called mode-shapes, indicate the fixed relative amplitude 
and phase of turbine rotor sections, when they are oscillating under steady state conditions at a 
given natural frequency. 
 
6.4 Left eigen vectors 
The left eigen vector, ?i, of the state matrix A associated with ?i, is any non-zero vector that 
satisfies equation (6.8).  
iii A ??? ?                                                                    (6.8) 
i  = 1, 2, ? n  
Just like the right eigen vectors, the left eigen vectors, ?i, can only be determined within a scalar 
multiplier and are usually scaled to have unit length.  
1?ii??                                                                            (6.9) 
61 
 
Let us define matrix ?, shown in (6.10),  as a diagonal matrix consisting of eigen values, ?i, of  
the state matrix, A. 
?
?
?
?
?
?
?
?
?
?
?
?
??
n?
?
?
?
????
?
?
00
00
00
2
1
                                                              (6.10)
 
Let us define matrices ? and ?, shown in (6.11) and (6.12), as consisting of right and left eigen 
vectors, ?i and ?i, of  the state matrix, A. 
? ?n?? ?1??                                                                       (6.11) 
? ?n?? ?1??                                                                     (6.12) 
Using (6.6) and (6.11) gives (6.13) and (6.14). 
????A                                                                             (6.13) 
????? A1                                                                          (6.14) 
Using (6.9), (6.11) and (6.12) gives (6.15) and (6.16). 
I???                                                                               (6.15) 
1????                                                                               (6.16)
 
The current work on eigen analysis involves computing these matrices in Matlab. First, 
the right eigen vector matrix, ?, is computed from the state matrix, A. Then the left eigen vector 
matrix, ?, is computed by inverting the right eigen vector matrix as shown in (6.16). 
The formulae discussed in the previous sections have been used to study the eigen 
properties of the linear time-invariant TG shaft system. Because of the configuration, physical 
properties and the differential equations governing the system, the isolated undamped TG shaft 
62 
 
system?s state matrix, A, consists of only real elements. Furthermore, all its eigen values occur in 
complex conjugate pairs with real parts equal to 0, with each of the complex conjugate pairs 
corresponding to one undamped oscillatory mode of the form given in (6.17).  
)cos( iii ta ?? ?                                                                         (6.17) 
 
6.5 Participation Factors [16] 
An additional important technique of quantifying how strongly the modes and states are 
associated is participation factors. These combine both the right and left eigen vectors of a mode 
to create a new vector that shows how this mode is associated with the system states as shown in 
(6.18).  
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
inni
ii
ii
ni
i
i
i
p
p
p
p
??
??
??
22
11
2
1
?                                                                (6.18) 
Collecting the participation vectors of all the modes gives us the participation matrix as shown in 
(6.19). 
? ?npppP ?21?                                                          (6.19) 
 
Because we will be dealing with a lot of TG shaft system analyses, let us first take a look at how 
they are derived with the aid of one STG shaft system example provided in [16].  
 
6.6 Derivation of State Equations for Kundur 1039 System  
These derivations are better illustrated with the aid of a test system whose rotor layout 
schematic is shown in Figure 6-2 and its electrical equivalent in Figure 6-3. This data belongs to 
63 
 
a 191 MVA 3,600 rpm STG unit. In this equivalent scheme, the rotor section masses are 
represented by capacitors and shaft spring constants are represented by inverse inductances. The 
torques transmitted in shaft springs and the rotor section velocities are represented by inductor 
currents and capacitor voltages respectively. 
 
 
 
 
 
 
 
 
 
Figure 6-2 Schematic Layout of Rotors and Shafts of Kundur pp. 1039 Problem 
 
 
 
Table 6-1  Rotor Inertias 
Rotor H (MW-s/MVA) 
HP 0.099 
IP 0.337 
LP 3.68 
G 0.946 
 
Table 6-2  Shaft Spring Constants 
Shaft K (p.u. torque/rad) 
HP-IP 37.95 
IP-LP 81.91 
LP-G 82.74 
 
 
 
 
 
 
 
 
 
 
 
 
HP IP LP G 
   
64 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
Figure 6-3 Electrical Equivalent of Shaft System 
 
The following differential equations for voltage across an inductor and current through a 
capacitor will be used to represent the dynamics of the electrical equivalent circuit shown in 
Figure 6-3.   
)()(
)()(
tvdttdiL
tidttdvC
?
?
      
 
Applying Kirchoff?s laws on the electrical equivalent circuit shown in Figure 6-3 and using these 
equations yields equations (6.20) through (6.24). 
 
 
)()( tidt tdvC IPHPHPHP ???                                                     (6.20) 
 
)()()( titidt tdvC LPIPIPHPIPIP ?? ??                                       (6.21) 
 
)()()( titidt tdvC GLPLPIPLPLP ?? ??                                       (6.22) 
CHP CIP CLP CG 
LHP-IP LIP-LP LLP-G vHP vIP vLP vG 
iHP-IP iIP-LP iLP-G 
65 
 
 
)()( tidt tdvC GLPGG ??                                                      (6.23) 
 
)()()( tvtvdt tdiL IPHPIPHPIPHP ????                                          (6.24a) 
 
)()()( tvtvdt tdiL LPIPLPIPLPIP ????                                           (6.24b) 
 
)()()( tvtvdt tdiL GLPGLPGLP ????                                            (6.24c) 
 
 
Using equations (6.20), (6.21), (6.22) and (6.23) yields (6.25). 
.  
 
dt tdvCdt tdvCdt tdvCdt tdvC LPLPIPIPHPHPGG )()()()( ????  
 
KtvCtvCtvCtvC LPLPIPIPHPHPGG ????? )()()()(                                         (6.25) 
 
 
Because the main concern in this work is eigen values, vectors and participation factors, i.e., 
dynamics only, the non-transient part, K, can be safely ignored as shown in (6.26). 
)()()()( tvCtvCtvCtvC LPLPIPIPHPHPGG ????                                                      (6.26) 
 
Substituting (6.26) in (6.24c) yields (6.27).  
 
)(1)()()( tvCCtvCCtvCCdt tdiL LPGLPIPGIPHPGHPGLPGLP ???????? ??????????????????????                        (6.27) 
 
 
Rewriting the differential equations, (6.20), (6.21), (6.22), (6.24a), (6.24b) and (6.27)  in matrix 
form yields (6.28). 
 
66 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
????
??
??
?
?
?
)(
)(
)(
)(
)(
)(
000
11
0
0000
11
00000
1
1
000
11
0000
0
11
000
)(
)(
)(
)(
)(
)(
tv
tv
tv
ti
ti
ti
CC
CC
C
CL
C
LCL
C
CL
C
LL
LL
dt
tdv
dt
tdv
dt
tdv
dt
tdi
dt
tdi
dt
tdi
LP
IP
HP
GLP
LPIP
IPHP
LPLP
IPIP
HP
GGLP
LP
GLPGGLP
IP
GGLP
HP
LPIPLPIP
IPHPIPHP
LP
IP
HP
GLP
LPIP
IPHP
   (6.28) 
 
The state matrix of the STG shaft system is now given by (6.29). 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
? ?
?
?
?
?
?
?
?
?
?
????
??
??
000
11
0
0000
11
00000
1
1
000
11
0000
0
11
000
LPLP
IPIP
HP
GGLP
LP
GLPGGLP
IP
GGLP
HP
LPIPLPIP
IPHPIPHP
CC
CC
C
CL
C
LCL
C
CL
C
LL
LL
A
          (6.29)
 
 
 
Using the generator ratings of 191 MVA and 3600 rpm as base values and choosing a time base 
value of 1 second establishes the base values for the rest of the relevant parameters as shown: 
?m-base = 377 mech-rad/sec, 
VA-base = 191,000,000 VA, 
67 
 
H-base =  ? ?
2
*2
basem
baseVA?   = 2687.6992 seconds, 
Torque-base = 506631.2997 N.m, 
J-base = 1343.849601 kg.m^2, and 
K-base = 1343.849601 N.m/mech-rad. 
 
Provided in Tables 6-1 and 6-2 are the rotor sections? moment of inertia values and shaft spring 
values in seconds and p.u. torque/rad. Using the base values computed, these are converted to S.I 
values and p.u. values as shown in Tables 6-3 and 6-4. 
 
Table 6-3 S.I. and p.u  Values of Kundur 1039 STG Rotor Section Moment of Inertias 
J (kg-m2) J-p.u C-p.u (C=J) 
266.0822211 0.198 0.198 
905.7546314 0.674 0.674 
9890.733066 7.36 7.36 
2542.563446 1.892 1.892 
 
Table 6-4 S.I. and p.u. Values of Kundur 1039 STG Shaft Spring Constants 
K (Nm/m.rad) K-p.u. L p.u. ( L = 1/K ) 
19226657.82 14307.15 6.98951E-05 
41498169.76 30880.07 3.23833E-05 
41918673.74 31192.98 3.20585E-05 
 
 These tables also provide the equivalent capacitances and inverse-inductances. Plugging these 
values into the state matrix, A, gives (6.30).  
 
68 
 
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
0000.1358-0.13580
00001.4836-1.4836
000005.05-
152535.611112.083264.38000
30880.1-30880.10000
014307.15-14307.15000
A
  (6.30)
 
 
The matrix of eigen values (W) and matrix of right eigen vectors (V), shown in Table 6-5, are 
computed using the Matlab command eig(A) on the state matrix, A. The matrix V consists of the 
right eigen vectors and the diagonal matrix W matrix gives the corresponding eigen values. The 
eigen values are 52.74, 29.66 and 22.33.   
 
 
Table 6-5 Right Eigen Vector Matrix of Kundur 1039 STG Shaft System 
0.792379 0.792379 0.326785 0.326785 0.068212 0.068212 
-0.609188 -0.609188 0.904531 0.904531 0.237142 0.237142 
0.028983 0.028983 -0.273741 -0.273741 0.969072 0.969072 
0.012077 i -0.012077 i 0.008856 i  -0.008856 i   0.002455 i   -0.002455 i  
-0.006275 i 0.006275 i 0.004600 i  -0.004600 i   0.001786 i   -0.001786 i  
0.000262 i  -0.000262 i -0.000859 i 0.000859 i   0.000709 i   -0.000709 i  
 
 
Each column in the matrix V gives the eigen vector of one eigen value in a diagonal position of 
the W matrix. Because we are interested in the rotor section velocities, these are selected and 
displayed in Table 6-6. Note that in the first column of Table 6-5, which represents the V matrix, 
there are six elements. These are the row elements in the eigen vector for the 52.74 Hz mode and 
correspond to the six state variables, iHP-IP(t), iIP-LP(t), iLP-G(t), vHP(t), vIP(t) and vLP(t), 
respectively. The row element for vG(t) is obtained using equation (6.26) and the row elements 
69 
 
for the remaining state variables are obtained as in the previous step. From these seven eigen 
vector elements, the four corresponding to the voltages are selected and normalized such that the 
maximum absolute value is 1. The resulting eigen vectors for the three modes are shown in Table 
6-6. 
Table 6-6 Three Right Eigen Vectors of Modes 52 Hz, 29 Hz, and 22 Hz Mode 
eigen vector for   52.74 hz :  
v1 HP 1 
v2 IP -0.51962 
v3 LP 0.021667 
v4 GENR -0.00383 
   
   
eigen vector for   29.66 hz :  
v1 HP 1 
v2 IP 0.519373 
v3 LP -0.097 
v4 GENR 0.087664 
   
   
eigen vector for   22.33 hz :  
v1 HP -0.6726 
v2 IP -0.48934 
v3 LP -0.19416 
v4 GENR 1 
 
 
The inverse of matrix V shown in Table 6-5 gives the matrix consisting of left eigen vectors, 
where each row is the left eigen vector of the corresponding eigen value in the W matrix. 
Multiplying the ith mode?s right eigen vector and left eigen vector gives the column of 
participation factors for the ith mode and the corresponding states. The participation factors for 
all the modes are given in Table 6-7. 
 
 
70 
 
Table 6-7. Participation Factors of Kundur-pp-1039 System 
Rotor 
52.73 Hz 
Mode 
29.66 Hz 
Mode 
22.33 Hz 
Mode 
HP 0.258171 0.213546 0.018504 
IP 0.237287 0.196085 0.03334 
LP 0.004505 0.074687 0.057315 
GENR 3.61E-05 0.015682 0.390841 
 
Impulse Response 
To obtain the impulse response magnitudes of the modes in individual rotor section 
velocities, the Kundur shaft system was built and simulated in ATPDraw/ATP software. The 
schematic for this system is shown in Figure 6-4. In this equivalent electrical circuit, the 
capacitors represent the rotor sections? moments of inertia and the inductances represent the 
inverse values of the shaft spring constants. 
71 
 
 
 
Figure 6-4 ATPDraw Schematic of Simulation of Impulse Response of Kundur 1039 Shaft 
System 
All the electrical quantities in this circuit are zero when the simulations begin. The 
impulse signal is applied as a short duration torque into the generator rotor section, which is 
represented by the capacitor connected to the node labeled XX0036. This causes the quantities to 
rapidly rise (for the duration of the impulse) and then ring-down back to zero. The 
ATPDraw/ATP software writes the capacitor voltage signals into a binary output.pl4 file. 
Another program is used read to the .pl4 output file and display the saved waveforms. Shown in 
72 
 
Figure 6-5 is the plot of the HP rotor section velocity signal. Note that the HP rotor section was 
represented by the capacitor connected to the node XX0020 in Figure 6-4. 
 
 
Figure 6-5 Plot of HP Rotor Section Velocity 
 
The software also allows us to convert these waveforms to ascii format and store them in 
a .txt file. Prony?s method and the reordering algorithm presented in Chapter 4 are then used to 
extract the constituent modes in each rotor section velocity signal. The results are given in Table 
6-8. Each column in this table gives the magnitude of the particular mode in a given signal. For 
73 
 
instance, shown in Column 1 (titled HP) are the magnitudes of the 52.7 Hz, 29.7 Hz and 22.3 Hz 
modes in the HP rotor section velocity signal: 0.0000104, 0.0019394 and 0.0038962 
respectively. The purpose of performing this test is to further illustrate the association of 
different modes and the states of the shaft system. This serves as a good complement to the 
theoretical analysis involving right eigen vectors and participation factors discussed in the 
previous sections. 
Table 6-8. Impulse Response Modal Magnitudes of Kundur-pp-1039 system 
  HP IP LP GENR 
Frequency(Hz) Magnitude Magnitude Magnitude Magnitude 
52.7 0.0000104 0.0000054 0.0000002 4.87E-08 
29.7 0.0019394 0.0010071 0.0001881 0.0001703 
22.3 0.0038962 0.0028337 0.001123 0.0057919 
 
6.7 Conclusions 
The basic concepts of state space theory are revisited and reviewed in this chapter in order to 
aid in analyzing the eigen characteristics of STG shaft systems. An STG system whose data was 
provided in [16] was used to illustrate the way in which these eigen properties, such as eigen 
values, right eigen vectors and participation factors are computed. ATPDraw/ATP software was 
used to simulate the impulse response of this shaft system. The purpose of impulse response 
analysis is to verify, to a certain degree, the theoretical approximations that were obtained using 
eigen analysis. These methods are used in the next chapter, in which four STG shaft systems are 
analyzed, to help reveal general characteristics of all these shaft systems. 
 
 
 
 
 
 
 
CHAPTER 7 
APPLICATION OF EIGEN ANALYSIS ON TEST CASES 
 
 Before the proposed methodology is applied to test a complete system (in the next 
chapter) consisting  of both the mechanical shaft system and the connected electrical power 
system, it is beneficial to study the behavior of a few typical STG shaft systems whose data is 
presented in the literature [4,16]. It is hoped that applying eigen analyses on these test cases will 
help reveal crucial information about how typical STG shaft systems aspects affect the proposed 
methodology?s capabilities in determining the parameters of STG systems. 
In this chapter, the concepts of eigen analysis introduced in the last chapter will be applied to 
four separate systems. These systems will be separately analyzed to help draw conclusions 
regarding the effect of these eigen properties on the measured magnitudes of rotor speed signals. 
The four shaft systems presented are:  
1. 555 MVA, 3600 rpm STG shaft system (Kundur pp 1036) [16], 
2. 960 MVA, 1800 rpm STG shaft system (Kundur pp 1038) [16],  
3. 191 MVA, 3600 rpm  STG shaft system (Kundur pp 1039) [16], and 
4. 892.4 MVA, STG shaft system, IEEE-FBM [4]. 
 
 
 
75 
 
 
7.1 Remarks on Terminology Used 
   Eigen analysis has proved to be an extremely valuable analysis tool in this work. It has helped 
reveal some fundamental aspects of the STG shaft systems that directly affect the proposed 
methodology?s performance. While applying these analyses the actual speed signals of various 
rotor sections will be referred to as state variables in the sense that these are all purely 
mathematical variables.  
   However, in real life there are some important constraints. First is that not every state variable 
in the STG shaft system is measureable. We only have access to speed signals on either end of 
the shaft system and the instantaneous power p(t) at the generator terminals. In most cases for 
large systems, the exciter is at the extreme end of the shaft system. Therefore, we do not have 
access to generator speed in most cases. In those cases we can try to extract the modes using 
measured instantaneous signal, p(t), but only to a certain extent. The reason for this is if a mode 
is present with a very small magnitude in the generator-section velocity state, then this small 
magnitude signal of the mode after passing through the electrical stator and the rotor circuits and 
the air-gap flux, which tend to attenuate it, becomes even smaller. It is important to consider 
these constraints as we proceed with the analysis. 
 
   Because certain terms and phrases appear very frequently in all the following sections, they are 
abbreviated to help increase clarity. These terms are as follows: 
HPRSV : High pressure rotor section velocity, 
GNRSV : Generator section velocity, 
HFM : Highest frequency mode, 
76 
 
LFM : Lowest frequency mode, 
STGS : Steam turbine generator shaft, and 
REV : Right eigen vectors. 
 
7.2 System 1: 555 MVA STG System from Kundur pp 1036 
 The schematic for this STG shaft system is provided in Figure 7-1. This shaft system has 
five rotor sections and thus has four natural frequencies: {16.3 Hz, 24.1 Hz, 30.3 Hz and 44.0 
Hz} that are also provided in [16]. The STG shaft system parameters are presented in tables 7-1 
and 7-2. These are first converted to p.u. using proper base quantities as discussed in the 
previous chapter and they are used to compute the system?s state matrix, A. All the relevant 
eigen properties such as left and right eigen vectors and participation factors are then computed 
from this state matrix using Matlab. 
 
 
 
 
 
Figure 7-1 Schematic Layout of Rotor and Shaft Sections of Kundur pp. 1036 system 
 
Table 7-1 Rotor Inertias 
Rotor H (MW-s/MVA) 
HP 0.124 
IP 0.232 
LPB 1.155 
LPA 1.192 
GENR 0.855 
HP IP 
 
LPB LPA 
   
G 
77 
 
Table 7-2 Shaft Spring Constants 
Shaft K (pu. torque/rad) 
HP-IP 21.8 
IP-LPB 48.4 
LPB-LPA 75.6 
LPA-GENR 62.3 
 
7.2.1 Right Eigen Vector Analysis  
The right eigen vectors (REVs) are shown in table 7-3. Note that these are normalized such that 
the maximum absolute value of any given vector element is 1. Because we are currently 
concerned about rotor section velocity states, only the corresponding entries in each eigen vector 
are shown. These eigen vectors are plotted in Figure 7-2. The y-axis is logarithmic because the 
values are of different orders of magnitude. Only the absolute values of right eigen vectors are 
displayed because the logarithmic axis does not allow plotting negative values and, more 
importantly, we are currently interested in the actual magnitudes (not the signs). 
 
Table 7-3  Right Eigen Vectors for Kundur-pp-1036 System 
Rotor 
43.99 Hz 
Mode 
30.30 Hz 
Mode 
24.08 Hz 
Mode 
16.22 Hz 
Mode 
HP -0.76569 1 1 1 
IP 1 -0.09383 0.309226 0.686315 
LPB -0.148 -0.50002 -0.18192 0.363603 
LPA 0.033732 0.708974 -0.15883 -0.14935 
GENR -0.00739 -0.43252 0.23824 -0.61422 
 
78 
 
 
Figure 7-2 Kundur 1036 Right Eigen Vector (mode shape) Absolute Magnitudes Across Rotor 
Sections 
 
 
Shown in Figure 7-2 are four curves with each curve illustrating the relationship between a given 
mode and the five rotor section velocity states in the system. For instance, the first curve shows, 
in a free response case, the relative magnitudes of the first mode with respect to each other as 
they appear in different rotor signals. The observed magnitude of the 44 Hz mode in the 
generator rotor section?s speed signal is about 100 times smaller than the observed magnitude in 
the H.P. rotor section speed signal. It is important to note that these are not indicative of the 
absolute magnitudes because the right eigen vectors are scaled such that the absolute maximum 
value is 1. This means that these charts can be used to infer information about relative 
magnitudes of a given mode in different rotor section speeds, and not regarding magnitudes of 
various modes in a given rotor section velocity signal. For instance, the relative magnitudes of all 
modes in the generator rotor section velocity (GRSV) measurements cannot be predicted. 
However, if the highest frequency modal (HFM) magnitude in one velocity state is known, it can 
be used to predict the HFM?s magnitude in another rotor section velocity state. For instance, in 
0.001
0.01
0.1
1
HP IP LPB LPA GENR
Magn
itu
de
 
Rotor Section 
43.99 Hz
30.3 Hz
24.08 Hz
16.2 Hz
79 
 
the next section, in table 7.3 on impulse response tests, it is given that the magnitude of the 44 
Hz mode in GNRSV is approximately 2e-7. Looking at REV of the 44 Hz mode in Figure 6-5, it 
can be predicted that the magnitude of the 44 Hz mode in HPRSV would be 100 times 2e-7 i.e., 
2e-5. This prediction is confirmed as indicated in table 7-3 where the magnitude of the 44 Hz is 
shown to be approximately equal to 2.2e-5. As mentioned earlier, using the REVs, the relative 
magnitudes of all modes in a given rotor section velocity measurement cannot be predicted. To 
do that, it is recommended to use participation factors, which will be discussed in the next 
section.   
    It can be concluded from Figure 7-2 that the measured magnitude of the highest frequency, 44 
Hz, varies significantly across rotor sections. Especially, as mentioned earlier, it appears very 
small in the generator rotor section velocity signal compared to its measured magnitude in the 
HP rotor section velocity signal. As mentioned in section 7.1, this small magnitude mode in the 
generator rotor section velocity signal would appear even smaller in the instantaneous power, 
p(t) signal.  It can also be inferred that there is only modest variation in measured magnitudes of 
lower frequencies across rotor sections when compared to the highest frequency. Thus, it can be 
deduced that the HP rotor section velocity signal is a much better candidate than p(t) for 
detection of STG shaft torsional modes, especially for higher frequency modal identification. 
 
7.2.2 Participation Factors 
    Right eigen vectors provide valuable information about the manner in which the modes show 
up in different velocity states. But, they do not represent the complete picture of how the modes 
and states are connected. To do that, it is also necessary to take into account the left eigen 
vectors, which depict how the initial conditions on states impact the magnitudes of the modes of 
80 
 
the system in a free response behavior. Participation factors combine both these vectors to 
illustrate completely how the states and modes are associated. Participation factors for the 
Kundur 1036 STGS system computed using (6.18) are given in Table 7-4 and illustrated in 
Figure 7-3. Just like in right eigen vectors charts, the y-axis in this chart is logarithmic because 
the values are of different orders of magnitude.  
 
Table 7-4. Participation Factors of Kundur-pp-1036 System 
Rotor 43.9 Hz  30.3 Hz  24.0 Hz  16.2 Hz  
HP 1.097E-01 5.281E-02 2.357E-01 8.434E-02 
IP 3.500E-01 8.700E-04 4.217E-02 7.433E-02 
LPB 3.817E-02 1.230E-01 7.267E-02 1.039E-01 
LPA 2.046E-03 2.552E-01 5.717E-02 1.808E-02 
GENR 7.047E-05 6.813E-02 9.226E-02 2.194E-01 
  
 
 
 
 
Figure 7-3 Kundur 1036 System?s Participation Factors Across All Rotor Sections 
 
 
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
43.9 HZ 30.3 HZ 24.0 HZ 16.2 HZ
Mag
nitu
de
 
Mode 
HP
IP
LPB
LPA
GENR
81 
 
Each curve in this chart illustrates how strongly a given rotor velocity state is associated with the 
different modes in the system. The participation factor for the 44 Hz mode across various RSV 
states varies tremendously. The lower frequency modes participate approximately the same 
amount across all RSV states. Mode 2 (30 Hz) participates approximately the same amount 
across all RSV states except in the IP RSV. These Figures indicate there will be difficulties 
measuring mode 4 in the generator rotor section velocity signal. To compensate for this, the HP 
rotor section velocity signal could be used as a measurement signal, because in this state all the 
modes participate equally and sufficiently well.  
 
7.2.3 Impulse Response Analysis 
    The third analysis that will be performed on the Kundur 1036 STG shaft system is the 
impulse response test. As discussed in the previous chapter, instead of using frequency domain 
methods to obtain the transfer function, the impulse response includes a time domain simulation 
of the shaft system to obtain the modal magnitudes in different rotor velocity states. This is in 
contrast to the first two analyses, which used right eigen vectors and participation factors, i.e., 
eigen properties of the shaft system. This allows us to compare not only the magnitude of a 
particular mode across rotor sections but also different modal magnitudes in given rotor sections. 
This allows us the opportunity to test the theoretical inferences made using right eigen vectors 
and participation factors analyses.  
   The isolated shaft system is simulated in ATP using an equivalent LC tank circuit where 
the rotor section mass inertias are represented by capacitors, and shaft spring constants are 
represented by inverse inductance values. Just as in full simulations, where the main disturbance 
torque enters the shaft system in the generator-section through an electrical torque signal, here 
82 
 
too the disturbance enters the generator rotor section in the form of a short impulse that 
represents a disturbance torque. Prior to the disturbance, all the states, capacitor voltages and 
inductor currents, are equal to zero. The effect of the disturbance, a short current pulse for about 
20 milliseconds, is that capacitor voltages and inductor currents states display free-response type 
behavior decaying down to zero. These speed signals are measured in all rotor sections and 
modal decomposition is performed using Prony?s method to determine the resulting modal 
frequencies and magnitudes. These results are shown in table 7-5 and Figure 7-4. As predicted 
earlier from eigen analysis, the 44 Hz mode magnitude is very small in the generator rotor 
section velocity signal, as compared to the HP rotor section velocity signal. In addition, as 
mentioned in section 7.1, this small magnitude mode in the generator rotor section velocity 
signal would appear even smaller in the instantaneous power, p(t) signal. Also, all the modes 
appear well in the HP rotor section velocity signal.  
 
Table 7-5. Impulse Response Modal Magnitudes of Kundur-pp-1036 system 
  HP IP  LPB LPA GENR 
Frequency (Hz) Magnitude Magnitude Magnitude Magnitude Magnitude 
43.9975 2.26E-05 2.96E-05 4.4E-06 0.000001 2E-07 
30.31358 0.001817 0.000172 0.000909 0.001292 0.000787 
24.08246 0.005947 0.001838 0.001083 0.000944 0.001415 
16.22699 0.006954 0.004773 0.002529 0.001039 0.004274 
 
 
 
 
83 
 
 
 
Figure 7-4 Kundur 1036 System Impulse Response Results 
 
 
 
7.2.4 Conclusions for Kundur 1036 System 
Looking at the chart for participation factors in Figure 7-3, it can be inferred that the highest 
frequency mode (HFM) does not participate sufficiently well in the generator rotor section 
velocity state. This fact was confirmed by impulse response test results shown in Figure 7-4 
where the measured magnitude of the 44Hz mode (HFM) was 10,000 times smaller than the 16 
Hz mode. In real-life situations, this would mean that it would be extremely hard to detect the 44 
Hz mode in the generator rotor section velocity signal. Moreover, it would be even harder to 
detect in the instantaneous power, p(t) signal. However, there is a solution to this problem. The 
REV or mode shapes, shown in Figure 7-2, indicate that the 44 Hz mode would be 100 times 
stronger in the HP rotor section velocity signal. Moreover, the rest of the modes participate 
sufficiently well in the HPRSV. This indicates that the HP rotor section velocity signal is a much 
1.0E-07
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
43.9 Hz 30.3 Hz 24.0 Hz 16.2 Hz
Mag
nitu
de
 
Mode 
HP
IP
LPB
LPA
GENR
84 
 
better candidate for detecting these torsional modes when compared to the generator rotor 
section velocity signal and a better candidate than the associated instantaneous power signal. 
 
 
7.3 System 2: 960 MVA STG shaft system from Kundur 1038 
 
      This 960 MVA STG shaft system is identical to the previous STG shaft system with respect 
to the layout. Consequently, it also has four natural torsional frequencies. They are 8.4 Hz, 15.2 
Hz, 20.3 Hz and 23.6 Hz. Just as in the previous section, all the eigen properties were computed 
using the system?s state matrix built using the converted p.u. data. 
 
 
 
 
 
Figure 7-5 Schematic Layout of STG Shaft System of Kundur pp. 1038 Problem 
 
Table 7-6 Rotor Section Inertias 
Rotor H (MW-s/MVA) 
HP 0.176 
LPC 1.427 
LPB 1.428 
LPA 1.428 
GENR 0.869 
 
 
 
HP LPC 
 
LPB LPA 
   
G 
85 
 
Table 7-7 Shaft Spring Constants 
Shaft K (pu. torque/rad) 
HP-LPC 17.78 
LPC-LPB 27.66 
LPB-LPA 31.31 
LPA-GENR 37.25 
 
7.3.1 Eigen vectors 
The right eigen vectors for this system are shown in Table 7-8. A plot of these vectors is shown 
in Figure 7-6.  
 Table 7-8  Right Eigen vectors for Kundur-pp-1038 system 
Rotor 
23.62 Hz 
Mode 
20.22 Hz 
Mode 
8.23 Hz 
Mode 
15.18 Hz 
Mode 
HP 1 -0.49422 1 -0.9212 
LPC -0.15693 -0.07506 0.859409 -0.48095 
LPB 0.045615 0.526171 0.139315 1 
LPA -0.01861 -0.99879 -0.58708 0.106318 
GENR 0.010776 1 -0.87799 -0.84163 
 
 
 
Figure 7-6 Kundur 1038 Right Eigen Vector Absolute Magnitudes Across Rotor Sections 
 
 
1.0E-02
1.0E-01
1.0E+00
HP LPC LPB LPA GENR
Mag
nitu
de
 
Rotor Section 
23.62 Hz
20.22 Hz
8.23 Hz
15.18 Hz
86 
 
The highest frequency mode (HFM) is 23.6 Hz. There is extensive variation in the right eigen 
vector (REV) of this mode across rotor section velocity states. Just like in the Kundur 1036 
system, the smallest eigen vector value is about 100 times smaller than the largest value. Also, 
this smallest eigen vector value occurs in the generator rotor section velocity state. Although the 
variation in the rest of the mode shapes is considerable, they are smaller compared to that of the 
23 Hz mode. 
 
7.3.2 Participation Factors 
The participation factors for this system are given in Table 7-9 and a plot of these factors is 
shown in Figure 7-7. It can be inferred from this figure that although the generator rotor section 
velocity state participates in lower frequency modes approximately uniformly, it has very little 
participation in the highest frequency mode (HFM) as indicated by the blue curve with the cross 
marker. The HFM participates very well in the HP rotor section velocity state, as indicated by the 
corresponding participation factor (4.09e-1), which is more than 1,000 times larger than the 
participation factor between the generator rotor section velocity and the 23 Hz mode (2.3e-4). 
Note that the participation factor between the 20 Hz mode and the HP rotor section velocity state 
is 100 times smaller than that of the participation factor between 20 Hz mode and the generator 
rotor section velocity state. However, it is not as small as the participation factor between 23 Hz 
mode and the generator rotor section velocity state, which is 1,000 times smaller than the 
participation factor between the 23 Hz mode and the HP rotor section velocity state. This 
establishes the HP rotor section velocity state as a better candidate to detect all the modes as 
compared to the participation factor between the 20 Hz mode and the generator rotor section 
velocity state. Recall that an identical conclusion was reached for the 1036 system earlier.  
87 
 
 
Table 7-9. Participation factors of Kundur-pp-1038 system 
Rotor 
23.62 Hz 
Mode 
20.22 Hz 
Mode 
8.23 Hz 
Mode 
15.18 Hz 
Mode 
HP 4.0986E-01 7.8448E-03 3.6368E-02 2.9411E-02 
LPC 8.1835E-02 1.4673E-03 2.1778E-01 6.4998E-02 
LPB 6.9192E-03 7.2146E-02 5.7270E-03 2.8120E-01 
LPA 1.1511E-03 2.5996E-01 1.0170E-01 3.1785E-03 
GENR 2.3498E-04 1.5858E-01 1.3842E-01 1.2121E-01 
 
 
 
 
 
 
 Figure 7-7 Kundur 1038 System?s Participation Factors Across All Rotor Sections 
 
 
 
7.3.3 Analysis of Impulse Response Modal Magnitudes  
The impulse response modal magnitudes for this system are given in Table 7-10 and a plot of 
these magnitudes is shown in Figure 7-8. The curve for the generator rotor section velocity 
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
23.62 Hz 20.22 Hz 8.23 Hz 15.18 Hz
Mag
nitu
de
 
Mode 
HP
LPC
LPB
LPA
GENR
88 
 
indicates that all the modes show up equally well in the generator rotor section velocity signal, 
except for the 23 Hz mode whose magnitude is about 1,000 times smaller than the rest. However, 
if the HP rotor section velocity signal is chosen as the measurement signal, then all the modes 
appear approximately equally well. Moreover, these are relatively equal or higher than those 
detected in the generator rotor section velocity signal. Again, the observations in section 7.1 
would imply that the 23 Hz mode would be even weaker in the instantaneous power, p(t), signal 
leading to the conclusion that the HP rotor section velocity signal is again a better candidate for 
modal extraction than the p(t) signal. 
Table 7-10. Impulse Response Modal Magnitudes of Kundur-pp-1038 System 
  HP LPC LPB LPA GENR 
Frequency 
(Hz) Magnitude Magnitude Magnitude Magnitude Magnitude 
23.63 0.000333 5.23E-05 1.52E-05 6.2E-06 3.8E-06 
20.21 0.001348 0.000206 0.001436 0.002725 0.00273 
15.18 0.0026 0.001359 0.002816 0.0003 0.002371 
8.235 0.003453 0.002968 0.000484 0.002024 0.003029 
 
 
 
  
Figure 7-8. Kundur 1038 System Impulse Response Results 
 
 
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
23.6 Hz 20.2 Hz 15.1 Hz 8.2 Hz
Mag
nitu
de
 
Mode 
HP
LPC
LPB
LPA
GENR
89 
 
7.3.4 Conclusions 
As observed previously in the Kundur 1036 system, here too the highest frequency mode 
(23HZ mode) participates very weakly in the generator rotor section velocity state, as validated 
by the impulse response measured magnitude in the generator rotor section velocity signal. In 
contrast, the measured magnitude of the 23 Hz mode in the HP rotor section velocity signal is 
much higher than in the generator rotor section velocity signal. Moreover, the magnitudes of all 
the modal components are sufficiently large enough to be detected in the HP rotor section 
velocity signal. This makes it a better candidate than the generator rotor section velocity signal 
and an even better candidate than the p(t) signal for detecting torsional modal components. 
 
7.4 System 3: 191  MVA STG Shaft System from Kundur pp 1039 
The third STG shaft system that will be studied is the Kundur 1039 system, whose data was 
presented in the previous chapter. For convenience the shaft configuration is presented here 
again in Figure 7-9. This system has four rotor sections and thus three torsional natural 
frequencies (22 Hz, 29 Hz and 52 Hz.)  
 
7.4.1 Eigen vectors Analysis 
The right eigen vectors for this system are shown in Table 6-6. A plot of these vectors is shown 
in Figure 7-10. The mode shape for the 52 Hz mode indicates that it shows up very 
insignificantly in the GNRSV signal, as compared to the rest of the rotor states. The mode shapes 
of the rest of the modes do not vary as much, indicating that they show up well in all the rotor 
section velocity states. 
 
90 
 
 
 
 
 
 
 
Figure 7-9  STG Shaft System Schematic for Kundur 1039 System 
 
 
 
  
 
Figure 7-10 Kundur 1039 Right Eigen Vector Absolute Magnitudes Across Rotor Sections 
 
 
 
7.4.2 Participation Factors 
The participation factors for this system are given in Table 6-7, and a plot of these factors are 
shown in Figure 7-11. 
 
 
1.0E-03
1.0E-02
1.0E-01
1.0E+00
HP IP LP GENR
Mag
nitu
de
 
Rotor Section 
52.7 hz
29.6 hz
22.3 hz
HP IP LP G 
   
91 
 
 
 
Figure 7-11 Kundur 1039 System?s Participation Factors Across All Rotor Sections 
 
 
From Figure 7-11, it can be inferred that the three modes participate very unevenly in the 
generator rotor section velocity state with the highest frequency mode (52 Hz) participating 
approximately 10,000 times weaker than the 22 Hz mode. Also, all the modes participate 
relatively equally and sufficiently in the HP rotor section velocity. 
 
 
7.4.3 Analysis Using Impulse Response Modal Magnitudes  
  The impulse response modal magnitudes for this system are given in Table 6-8, and a plot of 
these factors is shown in Figure 7-12. 
 
 
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
52.7 HZ 29.6 HZ 22.3 HZ
Mag
nitu
de
 
Mode 
HP
IP
LP
GENR
92 
 
 
 
Figure 7-12. Kundur 1039 System Impulse Response Results 
 
 
As predicted by the participation factors and the right eigen vectors, the 52 hz mode magnitude 
in the generator rotor section velocity state is very small (100 times smaller) compared to the 
measured magnitude in the HP rotor section velocity signal.  Moreover, all the measured modal 
magnitudes are higher in the HP rotor section velocity signal than in the generator rotor section 
velocity signal. These facts again point to the conclusion that the HP rotor section velocity signal 
is a better candidate for detecting the torsional modes than the generator rotor section velocity 
signal. In addition, it is an even better candidate than the instantaneous power signal, p(t), 
because of the reasons mentioned in section 7.1. 
 
7.5 System 4: 892.4  MVA STG Shaft System From IEEE-FBM 
The final test case is the so-called First Benchmark (FBM) case, a standardized test case, 
developed by the IEEE?s SSR Task Force based on the Navajo Project [4] for the industry to 
analyze different models, compare their results and validate their simulation programs. It uses 
1.0E-08
1.0E-07
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
52.7 Hz 29.6 Hz 22.3 Hz
Mag
nitu
de
 
Mode 
HP
IP
LP
GENR
93 
 
data that belongs to actual generators and turbines, thus making the case more realistic. The shaft 
system configuration and data are provided in Figure 7-13 and Table 7-11. 
 
 
 
Figure 7-13 STG Shaft System for IEEE-FBM 
 
Table 7-11. Rotor Spring Mass Parameters 
Inertia 
 
Spring Constant 
Mass H (sec) 
 
Shaft K (p.u.) p.u. torque / rad 
HP 0.092897 
 
HP-IP 7277 19.303 
IP 0.155589 
 
IP-LPA 13168 34.929 
LPA 0.85867 
 
LPA-LPB 19618 52.038 
LPB 0.884215 
 
LPA-GEN 26713 70.858 
GEN 0.868495 
 
GEN-EXC 1064 2.822 
EXC 0.0342165 
    
 
 
 
7.5.1 Eigen vectors 
The right eigen vectors for this system are shown in Table 7-12. A plot of these vectors is shown 
in Figure 7-14. 
 
Exc Gen LPB LPA IP HP 
TE ? 
94 
 
Table 7-12  Right Eigen Vectors for IEEE-FBM System 
Rotor 47.45 Hz Mode 32.28 Hz Mode 25.54 Hz Mode 15.71 Hz Mode 20.21 Hz Mode 
HP -0.7874 0.863783 1 -0.77714 -0.10989 
IP 1 -0.0437 0.342151 -0.58376 -0.06464 
LPA -0.11328 -0.50271 -0.22972 -0.34244 -0.01499 
LPB 0.021114 1 -0.09543 0.1117 0.039498 
GENR -0.00446 -0.62049 0.165974 0.373148 0.037344 
EXCR 0.000945 0.376842 -0.25255 1 -1 
 
 
 
Figure 7-14 IEEE-FBM Right Eigen Vector Absolute Magnitudes Across Rotor Sections 
 
 
Shown in the chart are mode shapes of the five modes with frequencies 15.7 Hz, 20 Hz, 25 Hz, 
32 Hz and 47 Hz. The mode shape for highest frequency mode, i.e., the 47 Hz mode, indicates 
that this mode appears 200 times smaller in the generator rotor section velocity state than in the 
HP rotor section velocity state. Moreover, as mentioned in section 7.1, because the generator 
speed is not a measurable quantity because of the presence of the exciter on the end, these modes 
would have be detected in the instantaneous power signal, p(t). Again, as mentioned in section 
7.1, it would be even harder to measure these modes as they would be even smaller, thus leading 
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
HP IP LPA LPB GENR EXCR
Mag
nitu
de
 
Rotor Section 
47.45 Hz
32.26 Hz
25.53 Hz
15.70 Hz
20.21 Hz
95 
 
to the conclusion that the HP rotor section velocity signal is the superior candidate of all these 
signals for torsional modal extraction.  
 
7.5.2 Participation Factors Analysis 
The participation factors for this system are given in Table 7-13, and a plot of these factors is 
shown in Figure 7-15. 
Table 7-13. Participation Factors of IEEE-FBM Machine 
Rotor 
47.45 Hz 
Mode 
32.28 Hz 
Mode 
25.54 Hz 
Mode 
15.71 Hz 
Mode 
20.21 Hz 
Mode 
HP 0.128211 0.02295 0.243714 0.074609 0.014467 
IP 0.346346 9.84E-05 0.047785 0.070507 0.008383 
LPA 0.024527 0.071851 0.118883 0.1339 0.00249 
LPB 0.000877 0.292775 0.021125 0.014671 0.017789 
GENR 3.85E-05 0.110717 0.062766 0.160812 0.015619 
EXCR 6.81E-08 0.001609 0.005726 0.045502 0.441252 
 
 
 
Figure 7-15 IEEE-FBM Participation Factors Across All Rotor Sections 
 
 
1.0E-08
1.0E-07
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
47.45 Hz 32.26 Hz 25.53 Hz 15.70 Hz 20.21 Hz
Mag
nitu
de
 
Mode 
HP
EXCR
GENR
LPB
LPA
IP
96 
 
For the sake of convenience, the two significant curves, participation factors of generator rotor 
section velocity and HP rotor section velocity states from Figure 7-15, are plotted in Figure 7-16. 
It can be inferred that although GNRSV participates equally and sufficiently among the 32 Hz, 
25 Hz, 15.7 Hz and 20 Hz modes, its participation factor with respect to 47 Hz mode is 
extremely small, indicating that there is very slight likelihood that the 47 Hz mode will be 
detected in the generator rotor section velocity signal. Because these generator rotor section 
velocity signals are not measureable, the instantaneous power signal, p(t), would have to be used. 
Again, as mentioned in section 7.1, it would be even harder to measure these modes as they 
would be even smaller in p(t). However, observing the HP rotor section velocity participation 
factor curve, it is clear that it participates equally and sufficiently well in all the modes, 
dramatically increasing the likelihood of the 47 Hz mode being detected in its measured signal. 
 
 
 
 
Figure 7-16 IEEE-FBM Participation Factors Across Generator and HP Rotor Sections 
 
 
 
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
47.45 Hz 32.26 Hz 25.53 Hz 15.70 Hz 20.21 Hz
Mag
nitu
de
 
Mode 
GENR
HP
97 
 
 
 
 
 
7.5.3 Impulse Response Modal Magnitudes  
The impulse response modal magnitudes for this system are given in Table 7-14, and a plot of 
these factors is shown in Figure 7-17. As predicted earlier from eigen analysis, the 47 Hz mode 
magnitude is very small in the generator rotor section velocity signal as compared to the HP rotor 
section velocity signal. Also, all the modes appear well in the HP rotor section velocity signal. 
 
Table 7-14. Impulse Response Modal Magnitudes of the IEEE-FBM System 
Frequency 
(Hz) HP   IP   LPA   LPB   GENR    EXCR   
47.45727 8.19E-06 1.04E-05 1.18E-06 2.19E-07 4.63E-08 1.01E-08 
32.26493 0.001562 7.88E-05 0.000909 0.001809 0.001121 0.000682 
25.53284 0.005395 0.001846 0.001237 0.000515 0.000894 0.001364 
15.70381 0.00649 0.004875 0.002861 0.000932 0.003117 0.00834 
20.21083 0.00079 0.000464 0.000108 0.000284 0.000268 0.007197 
 
 
 
  
Figure 7-17 IEEE-FBM Impulse Response Results: Relative Magnitudes of Modes Measured in 
Different Rotor Section Speeds 
1.0E-08
1.0E-07
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
47.45 Hz 32.26 Hz 25.53 Hz 15.70 Hz 20.21 Hz
Mag
nitu
de
 
Mode 
EXCR
GENR
LPB
LPA
IP
HP
98 
 
 
 
7.6 General Conclusions from Right Eigen Vector Analysis of All Four Systems 
 
For the sake of convenience, the mode-shapes (right eigen vectors) of the highest 
frequency modes of all four systems are plotted in a single chart in Figure 7-18. Because there 
are different numbers of rotor sections in these systems, the x-axis is not labeled. However, all 
the systems share a similar configuration in the sense that the HP rotor section is on the left hand 
side of the shaft system configuration and the exciter is on the right. Thus, we can choose to 
identify the x-axis unit labeled 1 as that of HP rotor section. Similarly, the unit labeled 4 
corresponds to generator rotor section, except in the case of 1039 system where the unit labeled 4 
corresponds to generator rotor section. Note that the exciter section of the FBM system is 
omitted for the sake of clarity. From these curves, it can be inferred that there is huge variation in 
all these mode-shapes and all the modal magnitudes tend to get extremely small in the generator 
rotor section velocity state as compared to the HP rotor section velocity state.  
 
 
  
Figure 7-18 Comparison of Mode Shapes (Absolute Values) of Highest Frequency Modes 
1.0E-03
1.0E-02
1.0E-01
1.0E+00
1 2 3 4 5
Mag
nitu
de
 
Rotor Section 
1036 Sys 43.99 Hz Mode
1038 Sys 23.62 Hz Mode
1039 Sys 52.74 hz
FBM 47.45 Hz Mode
99 
 
 
 
For the sake of convenience, the components corresponding to the HP rotor section velocity and 
generator rotor section velocity states of all the mode-shapes are separately plotted in Figure 7-
19. These are normalized such that the right eigen vector component corresponding to the HP 
state is 1. From this figure, it can be inferred that there is a significant amount of attenuation of 
these modal magnitudes from the HP rotor section velocity to the generator rotor section 
velocity. As mentioned in section 7.1, these modes in the generator rotor section velocity signal 
would be even more attenuated in the instantaneous power, p(t), signal. To remedy this situation, 
it is recommended that the HP rotor section velocity signal be used instead of instantaneous 
power, p(t) signal.  
 
  
 
Figure 7-19 Comparison of Normalized (with respect to HP State) Right Eigen Vectors of 
Highest Mode in Generator State 
 
 
9.7E-03 1.1E-02 
3.8E-03 
5.7E-03 
1.0E-03
1.0E-02
1.0E-01
1.0E+00
Ku
nd
ur 
103
6
Ku
nd
ur 
103
8
Ku
nd
ur 
103
9
IEEE
-FBM
Mag
nitu
de
 
System 
100 
 
Shown in Figure 7-20 is the improvement in magnitude of the highest frequency mode obtained 
by using the HP rotor section velocity signal instead of the generator rotor section velocity 
signal. It can be observed that there is a gain of approximately 100 in all the shaft systems.  
 
 
 
 
Figure 7-20 Improvement Attained using HP Rotor Section Speed Signal Instead of Generator 
Speed Signal for Measuring Highest Modal Component 
 
 
Although using the HP rotor section velocity signal instead of the generator rotor section velocity 
yields enormous gains in the highest mode, we have to be careful to not let the other modes 
attenuate. Shown in Figure 7-21 is the plot of the lowest frequency mode which has the highest 
change compared to others. However, as observed from the figure, these do not attenuate for 
three of the systems and for the system that attenuates, it does so only modestly ( gain = 0.67 ). 
This attenuation is miniscule compared to the attenuation (gain = 3.83e-3) that occurred in the 
same shaft system (Kundur 1039) for the highest mode when using the GNRSV signal instead of 
104 93 
261 
177 
1
10
100
1000
Kundur 1036 Kundur 1038 Kundur 1039 IEEE_FBM
Mag
nitu
de
 
System 
101 
 
the HPRSV signal. It can be concluded from this data that the HP rotor section velocity signal is 
again a much better candidate for detecting the modal components than the generator rotor 
section velocity signal or the associated p(t) signal because of the reasons mentioned in section 
7.1. 
 
 
 
Figure 7-21 Normalized (With Respect To Generator Speed State) Right Eigen Vectors of 
Lowest Frequency Mode in HP Rotor Section Speed State 
 
 
7.7 Overall Analysis and Conclusions Using Participation Factors 
Show in Figure 7-22 is the comparison of participation factors between the four STGS systems 
for the highest frequency mode. 
1.63 
1.14 
0.67 
2.08 
0.10
1.00
10.00
Kundur1036 Kundur 1038 Kundur 1039 IEEE-FBMMag
nitu
de
 
System 
102 
 
 
Figure 7-22 Comparison of Participation Factors of Highest Frequency Mode Between Generator 
and HP Rotor Section Speed States 
 
The first set, labeled Kundur 1036, shows the comparison between participation factors of the 
highest frequency mode and the HP rotor section velocity state, and participation factors of the 
highest frequency mode and the generator rotor section velocity states. The maroon bar indicates 
that the highest frequency mode participates very weakly in the generator rotor section velocity 
state (7e-5). However, it participates relatively strongly in the HP rotor section velocity state i.e., 
approximately 1,000 times better. The same trend continues in the rest of the three shaft systems. 
Overall, it can be inferred that, in all the cases, the highest frequency mode participates very 
strongly in the HP rotor section velocity state compared to the generator rotor section velocity 
state. This again leads to the conclusion that the HP rotor section velocity signal is a better 
candidate for detecting modes than the generator rotor section velocity signal, or its associated 
instantaneous power signal, p(t). 
 
 
1.1E-01 
4.1E-01 2.6E-01 
1.3E-01 
7.0E-05 
2.3E-04 
3.6E-05 3.9E-05 1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
Kundur1036 Kundur1038 Kundur1039 IEEE-FBM
Mag
nitu
de
 
System 
HP
Genr
103 
 
7.8 Overall analysis using Impulse Response Analysis 
The impulse response for all four STGS systems is shown in Figure 7-23. The first set, labeled 
Kundur 1036, indicates the measured magnitudes of the highest frequency modal component (44 
Hz) as detected in the generator rotor section velocity and HP rotor section velocity signals in the 
Kundur 1036 system. This indicates that this mode appears extremely small in the generator 
rotor section velocity signal and the measured magnitude in the HP rotor section velocity signal 
is larger (about 100 times). The remaining sets of pairs show the same relationship for the rest of 
the three STGS systems.  
 
  
Figure 7-23 Comparison of Measured Values of Highest Modal Component Between Generator 
and HP Rotor Section Speed Signal 
 
 
For the sake of convenience, this relationship is plotted in Figure 7-24. As shown in this 
figure there is a significant improvement, a gain of two orders of magnitude, in this modal 
magnitude when it is measured in the HP rotor section velocity signal instead of the generator 
rotor section velocity signal. It can also be inferred that the improvement occurs across all the 
shaft systems and in approximately similar amounts. This gain in modal magnitudes is going to 
1.0E-08
1.0E-07
1.0E-06
1.0E-05
1.0E-04
1.0E-03
1.0E-02
1.0E-01
1.0E+00
KUNDU
R 1036
KUNDU
R 1038
KUNDU
R 1039 IEEE-
FB
M
Mag
nitu
de
 
System 
HP
GENR
104 
 
be even higher when using HP rotor section velocity signal, instead of the instantaneous power 
signal, p(t), because of the reasons mentioned in section 7.1. 
 
 
 
Figure 7-24 Improvement of Using HP Signal Instead of Generator Signal for Measuring Highest 
Modal Component 
 
 
 
7.9 Conclusions 
  In this chapter, four typical STG shaft systems were analyzed using three types of 
analysis: right eigen vector analysis, participation factor analysis and the impulse response test. 
All these shaft systems have a common configuration of rotor masses aligned linearly and 
connected from one end (usually having the high-pressure turbine) to another end (usually 
having the generator/exciter). The generator and the adjacently connected low-pressure turbine 
usually tend to have the highest moments-of-inertia when compared to other rotors in the system. 
In all four cases that were tested, it has been observed that the higher frequency modal 
components tend to show up with much higher order magnitudes (approximately 100 times 
stronger) in the HP rotor section?s velocity signal than that of the generator. This could be 
113.0 87.6 
213.6 176.9 
1.0E+00
1.0E+01
1.0E+02
1.0E+03
KUNDUR 1036 KUNDUR 1038 KUNDUR 1039 IEEE-FBM
Mag
nitu
de
 
System 
105 
 
because the generator and adjacent rotor sections are considerably heavier than the HP and 
connected rotor sections. Moreover, all the modes? magnitudes were of similar order in the HP 
rotor section velocity state. This allows the option of choosing the HP rotor section velocity 
signal when the generator velocity signal or its associated instantaneous power signal, p(t) is 
insufficient to detect the required modal magnitudes. Thus, it can be concluded that the HP rotor 
section velocity signal is the more superior candidate for detecting STG shaft torsional modes 
than the generator rotor section velocity signal (which is not available for measurement in most 
cases) or its associated instantaneous power signal, p(t). 
 
 
 
 
 
 
 
 
CHAPTER 8 
IEEE-FBM FULL TEST SYSTEM APPLICATION AND ANALYSIS 
 
 In this chapter, two test cases, whose data have been provided in literature, will be 
utilized to illustrate and validate the proposed methodology. The first test case, the so called First 
Benchmark (FBM) case, is a standardized test case developed by the IEEE?s SSR Task Force 
based on the Navajo Project [4] for the industry to analyze different models, compare their 
results and validate their simulation programs. All the data used in this case belongs to actual 
generators and turbines, thus making the case more realistic. Moreover, the authors have 
provided enough data to help validate the proposed method sufficiently. The second test case 
consists of a similar STG shaft system, whose torsional data and generator electrical parameters 
have been provided in [16]. However, because the electrical network data was not provided, the 
FBM system?s electrical network configuration and parameters have been used. Both the systems 
have been built and simulated in ATPDraw software. The main concern is to apply the proposed 
method on this simulated model and validate it by verifying SSR frequencies and the shaft spring 
constant values. 
 
 
 
 
107 
 
8.1 Test Case 1: IEEE-FBM 
This test case consists of a six-stage steam turbine-generator connected to a series 
capacitor compensated transmission line, which is connected to an infinite bus on the other end. 
The network and the transmission line parameters are shown in Figure 8-1.  
 
 
 
 
 
 
 
 
Figure 8-1 FBM Test Case?s Network 
The electrical parameters of the 892.4 MVA generator and the rotor masses and shaft 
spring constants of the six turbine sections are given in Tables 8-1, 8-2 and 8-3.  
 
 
 
 
 
 
 
 
X1=j0.06 
X0=j0.06 
 
Infinite Bus 
Ground Fault 
Transmission Line 
Generator 
Filter 
Xl=j0.14 
X1=j0.50 
X0=j0.56 
 
R1=0.02 
R0=0.50 
 
XC= -j0.371 
 
 
Fault Bus 
108 
 
 
 
Table 8-1. Generator Impedances, Time Constants and Operating Conditions 
Generator Impedances (p.u.) 
 
Time Constants (sec) 
 
   
1.79 
 
 
     
4.3 
 
    
0.169 
 
 
      
0.032 
     0.135 
 
 
     
0.85 
 
   
0.13 
 
 
      
0.05 
 
   
1.71 
 
 
   
0 
 
    
0.228 
   
 
     
0.2 
 
Generator Operating Condition 
   
power output 0.9 p.u. 
   
power factor 0.9 lagging 
 
 
 
109 
 
 
Table 8-2. Rotor Circuit Parameters (p.u.) 
Parameter D-axis Q-axis 
Rf 0.00141 0.01406 
Xf 0.062 0.326 
Rk 0.0041 0.0082 
Xk 0.0055 0.095 
Xa 1.66 1.58 
XL 0.13 0.13 
 
 
 
Table 8-3. Rotor Spring Mass Parameters 
Inertia 
 
Spring Constant 
Mass H (sec) 
 
Shaft K (p.u.) p.u. torque / rad 
HP 0.092897 
 
HP-IP 7277 19.303 
IP 0.155589 
 
IP-LPA 13168 34.929 
LPA 0.85867 
 
LPA-LPB 19618 52.038 
LPB 0.884215 
 
LPA-GEN 26713 70.858 
GEN 0.868495 
 
GEN-EXC 1064 2.822 
EXC 0.0342165 
    
 
 
110 
 
This linear system consists of six rotor masses connected by five shafts, as shown in 
Figure 8-2. Although this system would have five resonant frequencies, only four frequencies 
have been provided in [4]. These are shown in the first four rows of Table 8-4. However, because 
the rotor masses and shaft spring parameters are specified in [4], all five frequencies have been 
determined using eigen analysis as described in Chapter 6. The fifth frequency, shown in Table 
8-4, is 47.45 Hz. 
 
 
 
Figure 8-2 Lumped Spring Mass Model of the Steam Turbine 
Table 8-4. List of SSR Frequencies 
Frequency Hz 
1 15.71 
2 20.21 
3 25.55 
4 32.28 
5 47.45 
 
 
 
8.1.1Simulation in ATP/ATPDraw Software 
The test case has been built and simulated in ATPDraw and ATP software, as shown in 
Figure 8-3.  
 
Exc Gen LPB LPA IP HP 
TE ? 
111 
 
 
Figure 8-3 ATPDraw Schematic Showing the IEEE FBM System 
 
The generator is simulated by the UM-1 model, which accommodates all electrical 
parameters of the generator. However, the mechanical rotor components have to be simulated 
external to the electrical model. The software simulates this mechanical system, which consists 
of rotor masses, shafts and torque inputs, with an electrical analog consisting of capacitors, 
inductors and constant current sources. A three-phase ground fault is applied for a duration of 
0.02 seconds at the fault bus as shown in Figure 8-1. The instantaneous terminal power, p(t), is 
shown in Figure 8-4.  
 
 
 
 
 
 
 
 
112 
 
8.1.2 Analysis of IEEE-FBM Simulation Results 
 
Figure 8-4. IEEE-FBM Systems? Instantaneous Power p(t) 
 
Prony?s method was used to estimate the constituent modal components in this signal. 
The resulting modes, sorted using the reordering algorithm discussed in Chapter 4, are shown in 
Table 8-5a.  The input parameters for Prony?s method are: tmin =  3.50, tmax =  4.0, numModes =    
120 and nSkipSample =  12. A comparison of actual vs. estimated signal waveforms is shown in 
Figure 8-4b.   
 
 
 
 
113 
 
Table 8-5a. List of 10 Significant Constituent Modal Components in the IEEE-FBM System?s 
Power Signal, p(t), Found Using Prony?s Method 
 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.276017667 0.0260722 15.8447405 2.1305827 
0.180805621 0.0860599 20.2470613 0.5575678 
0.094140055 0.0472314 11.9562734 0.8357333 
0.086593703 -0.0594629 25.5213181 -0.035871 
0.061468645 -7.4230014 3.2819871 1.4518369 
0.007693741 -0.7214964 31.5658185 -0.4270126 
0.021467517 -918.1780083 0 3.1415927 
0.006047653 -4.0827177 35.8677556 -1.686961 
0.009736733 -17.2561707 41.0963734 1.7500572 
0.008595033 -145.1546038 416.6666783 3.1415927 
 
 
 
Figure 8-4b Comparison of Actual Vs. Estimated Instantaneous Power Signals of IEEE-
FBM System  
 
114 
 
It is to be noted that not all the modes in the list shown in Table 8-5a are shaft torsional 
modes. Most are noise modes that are generated by Prony?s method because of the non-linear 
nature of p(t) signal as discussed in Chapter 3. However, it was observed that 11.9 Hz mode 
appeared consistently across multiple candidate data sets (defined by tmin and tmax). The existence 
of this mode could be attributed to the interaction of the series capacitance in the transmission 
line and other system components. For the sake of convenience, the list of shaft torsional modes 
are shown separately in Table 8-5b. 
 
Table 8-5b. Estimated Shaft Torsional Frequencies Values in the IEEE-FBM System?s Power 
Signal, p(t), Found Using Prony?s Method 
 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.276017667 0.0260722 15.8447405 2.1305827 
0.180805621 0.0860599 20.2470613 0.5575678 
0.086593703 -0.0594629 25.5213181 -0.035871 
0.007693741 -0.7214964 31.5658185 -0.4270126 
 
A comparison between the estimated frequencies and the FBM?s frequencies is given in 
Table 8-6. It is clear from the small percentage errors that the estimated frequency values match 
the FBM?s frequencies very closely. 
 
Table 8-6. Comparison Between Estimated vs. IEEE-FBM System?s Actual SSR frequencies  
SSR frequencies  (Hz)given 
in IEEE FBM paper 
Estimated SSR frequency (Hz) values extracted 
from the  power signal using Prony?s method 
Percentage 
Error 
15.71 15.8447405 0.86 
20.21 20.2470613 0.18 
25.55 25.5213181 0.11 
32.28 31.5658185 2.21 
 
115 
 
But, as shown in Table 8-5b, only four out of the five SSR frequencies could be 
extracted. The highest torsional frequency, 47 Hz, is missing. This is consistent with the results 
of our analysis provided in Chapter 6 using right eigen vectors, participation factors and impulse 
response tests on the IEEE-FBM shaft system. These analyses indicated that the highest 
frequency participated very weakly in the generator rotor velocity state. This causes it to show up 
even more weakly in the instantaneous power signal, p(t) resulting in it not being detected by the 
modal decomposition routine. Another reason might, among others, be the insufficient 
magnitude of the fault disturbance. To compensate for this, the process has been repeated, using 
the HP rotor section speed signal, as shown in Figure 8-5, in lieu of the power signal. It is to be 
noted from Chapter 7 analyses (Table 7-13) that the HP rotor section velocity state has a higher 
participation factor (0.128) with the 47 HZ mode compared to the generator section velocity state 
(3.85E-5).  
 
Figure 8-5. Speed Signal of the IEEE-FBM System?s HP Rotor Section 
116 
 
 
Prony?s method was used to estimate the constituent modal components in this signal. 
These modes, which also include the noise modes along with the shaft torsional frequencies, are 
shown in Table 8-7a.  The input parameters for Prony?s method are: tmin =  3.5, tmax =  4.0, 
numModes = 120 and nSkipSample =  10. For the sake of convenience, the list of shaft torsional 
modes are shown separately in Table 8-7b. 
 
Table 8-7a. List of 10 Significant Constituent Modal Components in the HP Rotor Section 
Velocity Signal Found Using Prony?s Method 
 
Magnitude ExpDecay Frequency (Hz) Phase (Rad) 
0.028999 0.033513 25.49858 -1.57742 
0.017561 0.048661 15.86501 0.541311 
0.00898 0.045958 20.24075 -0.77257 
0.0021 -0.03208 32.25314 0.850057 
0.000353 -10.3972 2.218337 -2.57716 
0.000238 -28.6691 0 0 
5.11E-05 0.004753 47.45543 0.662151 
6.08E-06 -3.96825 36.05672 2.404539 
4.92E-05 -652.728 173.248 1.602849 
5.72E-06 -42.2946 58.6955 1.433849 
 
 
Table 8-7b. Estimated Shaft Torsional Modal Values in the Speed Signal Found Using Prony?s 
Method 
 
Magnitude ExpDecay Frequency (Hz) Phase (Rad) 
0.028999 0.033513 25.49858 -1.57742 
0.017561 0.048661 15.86501 0.541311 
0.00898 0.045958 20.24075 -0.77257 
0.0021 -0.03208 32.25314 0.850057 
5.11E-05 0.004753 47.45543 0.662151 
 
117 
 
A comparison of actual vs. estimated signal waveforms is shown in Figure 8-6. As seen 
in this figure, Prony?s method was able to extract the underlying modal components with 
negligible error.  
 
Figure 8-6. Waveform Comparison of IEEE-FBM System?s Actual vs. Estimated HP signal 
 
Moreover, all five SSR frequencies have been extracted successfully. These SSR 
frequency values and the rotor mass values were used in an iterative least squares method to 
estimate the shaft spring constants. The results of the least squares method are shown in Table 8-
8. It can be seen that this method has been able to estimate the spring constants to within a 
maximum error of 5.6%.  
 
 
 
118 
 
Table 8-8. Comparison Between Estimated vs. IEEE-FBM System?s Actual Shaft Spring 
Constant Values (p.u.) 
Estimated Shaft Spring 
Constants ( found 
using Newton Raphson 
least squares method) 
Actual FBM Shaft Spring 
Constants (values provided 
in paper) 
Percentage 
Error 
7193.975986 7276.96 1.14 
13363.35297 13167.943 1.48 
20711.43256 19617.974 5.57 
25641.41227 26712.949 4.01 
1063.86061 1063.991 0.012 
 
On another note, shown in Table 8-7b, are the exponential decay values of the constituent 
SSR modal components along with the frequency values. Although these might indicate that the 
responses will increase exponentially to infinity, it may not necessarily be true. These damping 
values are simply the estimates determined by Prony's method for the signal for a particular 
chosen time interval. The aim of Prony's method is to minimize the error between the estimated 
and actual signal. These signals match closely with minimum error as seen in Figure 8-6. 
However, because these damping values are so small compared to the frequency values, they 
were not accurately determined.  
 
8.2 Test Case 2: Kundur 1036 System 
The second test case that will be used to validate the presented methodology is the 
Kundur 1036 system. The shaft system torsional and generator data for this system has been 
provided in Kundur [16]. This system consists of a steam turbine generator connected to a 
network similar to the one provided in FBM test case. The only exception is that the 
transmission line?s series capacitances have been reduced such that instead of providing 70 
percent compensation they provide 7 percent compensation.  The shaft system consists of five 
119 
 
rotor masses as shown in Figure 8-7. The torsional data of the shaft system is provided in Table 
8-9. The electrical parameters of the 2-pole, 555 MVA generator are given in Table 8-10. 
 
 
 
 
 
Figure 8-7 Kundur 1036 Shaft System Configuration 
 
Table 8-9 Specifications for Kundur 1036 Shaft System Parameters 
Rotor H ( MW.s/MVA) Shaft K (pu torque/rad) 
HP 0.124 HP-IP 21.8 
IP 0.232 IP-LPB 48.4 
LPB 1.155 LPB-LPA 75.6 
LPA 1.192 LPA-GENR 62.3 
GENR 0.855     
 
Table 8-10 Kundur 1036 System?s Generator Electrical Data 
Paramter p.u. value 
Ll 0.15 
Lad 1.66 
Laq 1.61 
Ld 1.81 
Lq 1.76 
Lafd 1.66 
Lffd 1.825 
Lfd 0.165 
Ra 0.003 
Rfd 0.0006 
 
JHP JIP JLPB JLPA 
KHP-IP KIP-LPB KLPB-LPA KLPA-G 
JG 
120 
 
This shaft system has the following resonant frequencies as provided in the textbook: 
16.3 Hz, 24.1 Hz, 30.3 Hz and 44.0 Hz. The schematic of the entire system consisting of the 
torsional system, the generator model and the connected electrical network, built in ATPDraw, 
is shown in Figure 8-8.   
 
8.2.1 Simulation in ATP/ATPDraw Software 
 
Figure 8-8 ATPDraw Schematic of the Kundur 1036 System 
 
Similar to the FBM test case, a three-phase ground fault is applied at the other end of the 
transmission line for 20 milliseconds and removed. The resulting ring-down signals, generator 
instantaneous power p(t) and HP rotor velocity ?HP(t), are recorded.  
 
8.2.2 Analysis of Kundur 1036 System Simulation Results 
Modal decomposition, based on Prony?s method, is used to identify the underlying modal 
components in both these signals. Results of modal decomposition on the power signal p(t) are 
shown in Table 8-11a, which includes noise modes along with a subset of shaft torsional modes. 
The input parameters for Prony?s method are: tmin =  4.5, tmax =  5.0, numModes =    180 and 
121 
 
nSkipSample =  10. For the sake of convenience, the list of shaft torsional modes are shown 
separately in Table 8-11b. A comparison of actual vs. estimated signal waveforms is shown in 
Figure 8-8b. 
 
Table 8-11a. List of 10 Significant Constituent Modal Components in the Instantaneous  
Power ,p(t), Signal Found Using Prony?s Method 
 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.100335294 -0.214798 1.0278007 1.4407871 
0.026822201 -0.0268801 16.2915718 1.4259349 
0.015395111 0.0635373 20.7025312 2.4360554 
0.007576683 0.0101545 24.1168428 3.1185639 
0.004403952 0.0200221 30.319575 -1.1807898 
0.001614289 -12.5397531 18.5551256 3.1060812 
0.000835548 -19.4412246 6.6322957 1.4899333 
0.000137558 -0.6841653 37.1079194 -2.6494426 
7.82391E-05 -3.971513 45.1327325 -1.1104041 
3.55971E-05 -2.9137386 50.9898275 1.0545063 
 
 
Table 8-11b. Estimated Shaft Torsional Modal Values in the Instantaneous  
Power ,p(t), Signal Found Using Prony?s Method 
 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.026822201 -0.0268801 16.2915718 1.4259349 
0.007576683 0.0101545 24.1168428 3.1185639 
0.004403952 0.0200221 30.319575 -1.1807898 
 
 
 
122 
 
 
Figure 8-8b. Waveform Comparison of Kundur 1036 System Actual vs. Estimated p(t) signal 
 
The 20.7 Hz mode shown in Table 8-11a was observed consistently across multiple 
candidate data sets (defined by tmin and tmax). Just as in the case of FBM system, the existence of 
this mode could also be attributed to the interaction of the series capacitance in the transmission 
line and other system components. As seen in Table 8-11b, the procedure was able to detect 
only three of the four resonant frequencies, just like in FBM case. The highest torsional 
frequency, 44 Hz, is missing. This is again consistent with the results of our analysis provided 
in Chapter 6 using right eigen vectors, participation factors and impulse response tests on the 
Kundur 1036 shaft system. These analyses indicated that the highest frequency participated very 
weakly (participation factor 7.05E-5) in the generator rotor section velocity state. This causes it 
123 
 
to show up even more weakly in the instantaneous power signal, p(t) resulting in it not being 
detected by the modal decomposition routine.  
 However, these same analyses-right eigen vectors, participation factors and impulse 
response tests-on the Kundur 1036 shaft system also indicated that all the modes participated 
strongly in the HP rotor section velocity state (participation factor = 1.1E-1). Therefore, modal 
decomposition has been performed on the HP rotor section velocity signal. The results are 
shown in Table 8-12a and Table 8-12b. The list of modal components in Table 8-12a includes 
both shaft torsional and noise modes. For the sake of convenience, the list of shaft torsional 
frequencies are shown separately in Table 8-12b. A comparison of actual vs. estimated signal 
waveforms is shown in Figure 8-8c. 
 
Table 8-12a Results of Modal Decomposition on the HP Rotor Section Velocity Signal 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.004919104 0.0004694 24.1010872 -1.7483499 
0.004633141 -0.0040039 16.2908815 -0.2803303 
0.001879906 0.0012081 30.324916 -2.9628608 
0.000735043 -0.1783341 1.0251133 2.9835261 
6.60606E-05 0.0004634 43.9975184 0.0776786 
1.1878E-06 -16.5067852 17.756971 1.2672097 
4.963E-07 -2679.291063 0 0 
2.137E-07 -34.1110417 381.9555055 0.9443255 
3.966E-07 -61.5138752 33.0911332 -3.008665 
9.29E-08 -38.2153551 99.2265085 -1.957044 
 
Table 8-12b Results of Modal Decomposition on the Kundur 1036 System?s HP Rotor Section 
Velocity Signal 
 
Magnitude Exp-Decay Frequency (Hz) Phase (Rad) 
0.004919104 0.0004694 24.1010872 -1.7483499 
0.004633141 -0.0040039 16.2908815 -0.2803303 
0.001879906 0.0012081 30.324916 -2.9628608 
6.60606E-05 0.0004634 43.9975184 0.0776786 
124 
 
 
Figure 8-8c. Waveform Comparison of Kundur 1036 System?s Actual vs. Estimated in HP Rotor 
Section Velocity Signal 
 
As seen in Table 8-12b, all four frequencies have been extracted precisely. Next, these 
frequencies were fed into a non-linear Newton-Raphson based algorithm to estimate the shaft 
spring constants. The results are shown in Table 8-13. The units of the spring constants are 
Nm/mech-rad. 
 
Table 8-13 Comparison Between Estimated and Kundur 1036 System?s Actual Spring Constant 
Values 
 
Estimated Actual Error  % 
30760262.08 32092838.2 4.152 
74268733.87 71251989.39 4.234 
112080995.1 111294429.7 0.707 
92176056.7 91714854.11 0.503 
    
125 
 
These results show that two spring constants are estimated to within 1 percent. The maximum 
error in the estimated spring constant values is about 4.2 percent. 
 
8.3 Conclusions 
 The results of time-domain simulations on two test cases were reviewed in this chapter. 
These test cases were built and simulated in ATPDraw software. The ring-down signals 
resulting from a three-phase ground fault on a transmission line were recorded. Modal 
decomposition on the instantaneous power signal revealed that most modes, except the highest 
one, were detected. This is to be expected as per the results of eigen analyses performed in 
Chapter 6. Consequently, alternative ring-down signals from other rotor sections were chosen as 
candidates for modal decomposition. These results of a modal decomposition on IEEE FBM?s 
HP rotor section velocity and Kundur 1036 HP rotor section velocity were sufficient to detect 
all the torsional modes. These modal frequency values were used in a non-linear least squares 
algorithm to estimate, within reasonable (6 percent) error, all the shaft spring constants of the 
respective rotor shaft systems. In the next chapter, real data obtained from a steam turbine 
generator unit will be used for modal decomposition using Prony?s method. 
 
 
 
 
 
 
 
CHAPTER 9 
REAL MEASUREMENT ANALYSIS 
 
Up until now, Prony's method and its refinements, eigen analysis/tools and the behavior 
of typical steam turbine generator shaft (STGS) systems have been examined. In addition, in the 
previous chapter, the proposed methodology was tested using test systems that were simulated in 
ATP software and shaft spring constants were estimated using HP rotor section velocity (RSV) 
signals and not instantaneous power, p(t), signals. However, these data sets were simulated. In 
this chapter, real measurement data will be utilized to test the methodologies presented in 
previous chapters. The effects of the behavior of typical STG shafts, studied in Chapter 7, on the 
proposed methodology will also be examined. 
This real data came from a steam turbine generator unit whose information will be kept 
confidential as requested by the plant owner. However, my advisory committee Chair, Dr. Mark 
Halpin, has full knowledge of the unit?s STGS torsional frequencies and parameters. Note that in 
this dissertation work, the only type of field data available was three-phase instantaneous power 
at the generator terminals. Rotor section speed signals were not available. 
 
9.1 Data Collection  
  A measurement device was placed for a few weeks on the premises to collect and record 
data whenever a disturbance occurred. The measurement device consists of an A/D convertor 
capable of reading six channels, three for instantaneous currents and three for instantaneous 
127 
 
voltages. Later these recorded quantities were converted to an instantaneous power signal, p(t) 
using (9.1).   
p(t) = va(t)*ia(t)+ vb(t)*ib(t)+ vc(t)*ic(t)                                (9.1) 
Because the reading of data from six channels occurred sequentially, there was a 180 Hz 
modulation frequency embedded in the p(t) signal which had to be removed using a filter. 
Whenever a disturbance occurred, a detection algorithm would trigger a procedure to record 30 
seconds of instantaneous voltage and current quantities, which included 2 seconds of data prior 
to the occurrence of disturbance. The data was recorded at 1,920 samples per second 
corresponding to 32 samples per cycle.  The collected data was written to an ascii file. All the 
data were in p.u. 
 
9.2 Analysis of Data Set 1 
A plot of the measured instantaneous power signal, p(t), from first data set is shown in 
Figure 9-1. The reason it looks like a strip of color instead of a single continuous curve is the 
aforementioned 180 Hz modulation frequency embedded in the p(t) signal. This is mainly 
because of the sequential sampling of three voltage quantities followed by three current 
quantities.  
 
128 
 
 
Figure 9-1 Plot of Actual p(t) Measurement Data 
 
To remove this modulation, this data was filtered using a seventh order low-pass 
Butterworth filter with a cutoff-frequency of 59 Hz. A plot of the filtered data is shown in Figure 
9-2. It depicts the ring-down signal as a result of a major disturbance event, which seems to 
begin at 1.7 seconds. To obtain a closer look at the data, the plot is zoomed as shown in Figure 9-
3 
 
129 
 
 
Figure 9-2 Plot of Filtered p(t) Signal from Data Set 1 
 
 
Figure 9-3 Zoomed-in Plot of the Disturbance Event in the p(t) Signal 
130 
 
The next task is to visually inspect the data and choose a proper candidate data set. The 
purpose of this visual inspection is to preclude any non-linear events such as sudden 
discontinuities. This would disqualify the data before 1.9 seconds (approximately) and after 4.5 
seconds (approximately). A data window between 2 seconds and 3 seconds has been chosen and 
Prony?s method has been applied to extract the underlying modal components. A plot of the data 
between 2 and 3 seconds is shown in Figure 9-4  
 
Figure 9-4 Plot of Proper Candidate Data Set for Prony?s Method (Data Set 1) 
 
To extract the correct underlying modal components, it necessary to select suitable 
algorithm parameters of Prony?s method as discussed in Chapter 5. After gradually modifying 
the algorithm parameters and checking the results of Prony?s method, it was found that suitable 
values for numModes and nSkipSample parameters were 90 and 5 respectively. The output of 
131 
 
Prony?s method using these algorithm parameters is shown in table 9-1. Note that only the first 
20 significant modes are shown. In addition, this list is ordered according to the LSE-based 
reordering algorithm presented in Chapter 4. 
 
Table 9-1 Prony?s Method?s Output List of Constituent Modal Components for Data Set 1 
Mode #  LSE Magnitude ExpDecay Frequency(Hz) Phase(rad) 
1 0.145430301 0.012577665 -0.3309921 1.5488809 1.0758312 
2 0.027445814 0.002169144 -0.2515075 12.0571053 -0.3224135 
3 0.007599871 0.00032536 -0.0641375 16.335785 2.9754293 
4 0.007306082 0.000376833 -0.7515305 20.8418966 2.1304798 
5 0.006588591 0.000202878 -0.6714491 59.8825575 2.2886434 
6 0.006550145 0.000184078 -0.5808247 24.952691 -3.086796 
7 0.006499929 0.000460505 -6.4398204 38.9980768 2.0071578 
8 0.006364739 0.000349181 -8.4033961 27.8239136 -0.5632267 
9 0.006357859 0.00020614 -2.9386237 56.6811368 -0.2796263 
10 0.006328438 0.000181467 -3.0516521 52.5733957 -2.0632532 
11 0.006313769 0.000164204 -3.0921935 44.8612443 -0.6163595 
12 0.006301417 0.000273349 -10.2966543 36.4077529 0.6941915 
13 0.00629955 0.00024585 -8.5860998 32.7887102 0.7029586 
14 0.006295845 0.000748658 -110.4579566 60.1552674 -1.030162 
15 0.00629521 0.000112727 -1.9587407 48.0411351 -2.931303 
16 0.006283306 0.000121469 -3.1095516 7.5468043 1.9888544 
17 0.006279016 0.003337735 -592.1246285 152.0977954 1.7286927 
18 0.006259379 5.57457E-05 -1.8085119 63.2491252 0.3082231 
19 0.006254104 0.000176545 -73.2794058 0 0 
20 0.006253283 4.86265E-05 -2.5630077 67.2675652 2.8706005 
 
A comparison of the actual signal and reconstructed signal using all the constituent 90 
modes is shown in Figure 9-5. Also shown in this figure is the difference between these two 
signals. 
 
 
 
132 
 
 
 
Figure 9-5 A Plot of Comparison Between the Actual and Reconstructed Signal Using the Output 
Modes of Prony?s Method for Data Set 1 
 
The LSE-based reordering algorithm, presented in Chapter 4, had been successfully 
utilized in Chapter 7 for identifying the constituent modes in impulse response signals and in 
Chapter 8 for identifying torsional modes from output signals of systems simulated in the 
ATPDraw software. However, because of the poor quality of the actual signals in the present 
cases, it was found that LSE-based reordering algorithm was not effective in identifying these 
modes. However, in the absence of that, an alternative simpler method of identifying these 
torsional modes was chosen. It utilizes a selection criteria of { -4 < exponential-decay < 4 } and 
133 
 
{ 5 < frequency < 59 } to prune away the non-STGS torsional modes and noise modes. Any 
constituent modes outside these bounds were excluded. In addition, because this method is not as 
effective in identifying these torsional modes (because of appreciable measurement error and the 
large number of cases involved), it was necessary to use a scatter plot to identify these modes. 
This method is detailed in the following sections. Applying the selection criteria on the list of 
modes in Table 9-1 gives a new list of modes as shown in Table 9-2. 
 
Table 9-2. Selected List of Constituent Modes in the Signal Detected by Prony?s Method  
magnitude Exp-Decay Frequency(Hz) Phase(rad) 
0.000121 -3.10955 7.546804 1.988854 
0.000113 -1.95874 48.04114 -2.9313 
0.000164 -3.09219 44.86124 -0.61636 
0.000181 -3.05165 52.5734 -2.06325 
0.000206 -2.93862 56.68114 -0.27963 
0.000184 -0.58082 24.95269 -3.0868 
0.000377 -0.75153 20.8419 2.13048 
0.000325 -0.06414 16.33579 2.975429 
0.002169 -0.25151 12.05711 -0.32241 
 
 
Note that as a result of using the selection criteria, Table 9-2 excludes the most dominant 
mode in the signal, which has a frequency of about 1.54 Hz, as seen in row 1 of Table 9-1. It is 
observed that this low-frequency mode, also known as the system mode [16], is present in most 
STGS signals. In the system mode, the entire STGS system oscillates in unison against the power 
system?s electrical network. In most STGS signals, this mode lies in the frequency range of 1 to 
2 Hz [16]. It was found that this mode is also present in most of the present measurement signals. 
Because the main concern is extracting individual torsional frequencies, the selection criteria?s 
lower limit was set as 5 Hz so that this mode would be omitted in further results.  
134 
 
In this procedure, a sliding window analysis was used to detect constituent modes in the 
signal because it is necessary to use it to eliminate any noise modes present in the output list of 
modes. The first window began at 2.0 seconds and ends at 3.0 seconds. The next window was 
chosen to begin at 2.1 seconds and end at 3.1 seconds. The same algorithm parameters chosen 
earlier were chosen for this data window and the rest of the data windows in this analysis. The 
output modal list is shown in table 9-3. This procedure was repeated 10 times to obtain a total of 
10 modal lists. The final list, starting at 3.0 seconds, of constituent modes is shown in Table 9-4.  
 
Table 9-3. List of Constituent Modes in the Signal Detected by Prony?s Method with Algorithm 
Parameters: tmin =  2.10 tmax =  3.10, numModes =   90 and nSkipSample =   5. 
magnitude Exp-Decay Frequency(Hz) Phase(rad) 
0.000118 -2.66588 52.56122 -0.32567 
0.000114 -3.59321 39.85057 0.6168 
0.000239 -2.80993 7.612871 0.734371 
0.000104 -1.79302 48.44691 1.229605 
8.11E-05 -1.77031 44.98968 2.002963 
0.000151 -2.99357 56.51658 -2.17675 
0.000212 -0.67029 25.0694 -0.02513 
0.000261 0.027198 20.92769 2.17706 
0.000408 -0.44486 16.37496 0.524382 
0.002154 -0.2917 12.06021 0.952535 
 
Table 9-4. List of Constituent Modes in the Signal Detected by Prony?s Method with Algorithm 
Parameters: tmin =  3.00 tmax =  4.00, numModes =   90 and nSkipSample =   5. 
magnitude Exp-Decay Frequency(Hz) Phase(rad) 
0.000106 -3.62959 32.86421 -3.01986 
0.000141 -1.94671 55.64619 -1.23221 
0.000232 -3.7667 7.118301 2.090702 
0.000325 -1.93688 25.08226 -2.77231 
0.000264 -0.60996 16.37752 -2.04934 
0.000235 -0.22158 20.80594 1.113286 
0.001679 -0.19653 12.04989 0.031397 
 
135 
 
  In the next step, the frequency data and exponential decay data from all 10 tables were 
chosen and combined into one list. The list is then sorted in ascending order according to 
frequency. The resulting list of modes is shown in Table 9-5. 
Table 9-5 Combined list Modal Components showing only the Exponential Decay and 
Frequency parameters 
Exp-Decay Frequency(Hz) 
 
Exp-Decay Frequency(Hz) 
 
Exp-Decay Frequency(Hz) 
-4.4149138 7.0914403 
 
-1.2673978 20.7462942 
 
-2.0158872 40.055229 
-3.7667008 7.1183011 
 
-0.2215773 20.8059387 
 
-1.4376032 40.1428554 
-4.7144464 7.1467739 
 
-0.1764114 20.8198008 
 
-2.2263055 44.8205265 
-4.4628948 7.2519218 
 
-0.7060205 20.8311284 
 
-3.0921935 44.8612443 
-1.0985038 7.3187341 
 
-0.7515305 20.8418966 
 
-4.2391012 44.9771431 
-4.4055675 7.5183069 
 
-0.8879877 20.8465696 
 
-1.7703147 44.9896791 
-1.714265 7.5407015 
 
-0.563378 20.8829791 
 
-3.6889312 45.0481237 
-3.1095516 7.5468043 
 
0.0271976 20.9276871 
 
-0.2411489 45.1271992 
-3.2565161 7.5586456 
 
-0.891101 20.9488363 
 
-3.8468226 45.1336808 
-2.8099317 7.6128714 
 
-1.2257292 20.9652567 
 
-1.6546833 45.221385 
-0.3113589 12.0355795 
 
-0.5808247 24.952691 
 
-1.4836504 45.2778428 
-0.216627 12.0457647 
 
-0.3921951 24.9758553 
 
-1.9587407 48.0411351 
-0.2209732 12.0474776 
 
-1.0349055 25.0010154 
 
-4.8350115 48.431882 
-0.2157098 12.0491428 
 
-1.4008373 25.0097467 
 
-1.7930235 48.4469058 
-0.1965312 12.0498854 
 
-1.2242982 25.0582747 
 
-4.1797655 48.6226286 
-0.1737293 12.055951 
 
-0.6702933 25.0694005 
 
-1.6994487 48.625585 
-0.2400201 12.0564182 
 
-1.9368757 25.0822632 
 
-0.2213687 52.105882 
-0.2515075 12.0571053 
 
-0.4763268 25.0870394 
 
0.3564415 52.4073332 
-0.2917031 12.0602059 
 
-0.4124591 25.1296092 
 
-3.7171471 52.4537301 
-0.2497538 12.0619462 
 
-0.8554051 25.1898511 
 
-2.6658769 52.5612168 
-0.235234 12.0674471 
 
-0.875131 25.2777751 
 
-3.0516521 52.5733957 
-0.6734198 16.2825171 
 
-4.8668484 32.062936 
 
-4.6596633 52.5916038 
-1.2011537 16.2919621 
 
-3.8070836 32.2774069 
 
-3.5227331 52.7169151 
-1.4117868 16.2957412 
 
-3.7041273 32.558519 
 
-4.7201381 52.7376111 
-0.7071557 16.3104169 
 
-4.6223211 32.7482586 
 
-3.211522 52.7600858 
-1.2077944 16.3278117 
 
-4.9038894 32.8445151 
 
-1.9467082 55.6461888 
-0.5884168 16.3301704 
 
-3.6295861 32.8642094 
 
-2.3711348 55.8850648 
-0.0641375 16.335785 
 
-4.4796474 36.197765 
 
-2.0572016 55.9123034 
-0.6497332 16.34294 
 
-4.0538193 37.5258886 
 
-3.1837891 55.9925392 
-0.4901767 16.3680153 
 
-4.727584 37.6079022 
 
-3.2786228 56.0908796 
-0.4448611 16.3749555 
 
-2.8012967 37.8163192 
 
2.6760975 56.0990834 
-0.6099618 16.3775186 
 
-3.5932084 39.8505682 
 
-2.9935713 56.5165834 
-0.5316564 20.732887 
 
-2.2440638 39.8671388 
 
-2.9386237 56.6811368 
136 
 
 
 
This list was then used to create a scatter plot as shown in Figure 9-6. From this figure it 
is obvious that there is are four clusters with large numbers of modal components, as illustrated 
in Figure 9-7, where  these four clusters were marked and labeled. It appears that there are four 
dominant modal components in the signal with frequencies: 12 Hz, 16.3 Hz, 20.7 Hz and 25 Hz.  
To illustrate the number of occurrences in these clusters a histogram based on frequency 
components is shown in Figure 9-8.  
 
 
Figure 9-6 Scatter Plot of Constituent Modes of Data Set 1 
0
10
20
30
40
50
60
-6 -5 -4 -3 -2 -1 0 1 2 3 4
Fre
qu
en
cy 
( H
z) 
Exp-Decay 
Scatter Plot of Modal components For Data Set 1  
137 
 
 
 
Figure 9-7 Scatter Plot of Constituent Modal Components with Major Clusters Demarcated to 
Show the Four Dominant Modes 
 
 
Figure 9-8. Histogram of Modal Components for Data Set 1 with Bin Width = 1 Hz 
0
2
4
6
8
10
12
1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55 57 59 61
Nu
mb
er 
of Occ
ur
ren
ce
s 
Bin Frequency  (Hz) 
Histogram of Data Set 1 Modal Components 
138 
 
 
9.4 Analysis of Data Set 2  
The above process was repeated on another measured instantaneous power signal, p(t), 
which consists of 30 seconds of measurement data. Just as in previous case, this data was filtered 
to remove the 180 Hz and other high frequency components. Shown in Figure 9-9 is a plot of the 
resulting waveform.  
 
Figure 9-9 Plot of Filter p(t) signal from Data Set 2 
After visually inspecting the data, a data window between 5.5 seconds and 6.5 seconds 
was chosen, shown in Figure 9-10, as the proper candidate for Prony analysis. Prony?s method 
was then applied on this window, using the same algorithm parameters as used for data set 1 ( 
numModes = 90 and nSkipSample = 5.) Shown in Figure 9-11 is a comparison of the specified 
signal and the reconstructed signal using the output modal list of Prony?s method.     
139 
 
 
Figure 9-10 Plot of Proper Candidate Data Set for Prony?s Method (Data Set 2) 
 
Figure 9-11 A Plot of Comparison Between the Actual and Reconstructed Signal Using the 
Output Modes of Prony?s Method for Data Set 2 
140 
 
 
Shown in Table 9-6 are the first 20 constituent modes obtained from Prony?s method. Note these 
are sorted in descending order of the LSE norms.  
Table 9-6 Prony?s Method?s Output List of Constituent Modal Components for Data Set 2 
Mode #  LSE  Magnitude Exp-Decay Frequency(Hz) Phase(Rad) 
1 0.142834448 0.014368914 -0.8583511 1.5861547 -0.6056632 
2 0.015434505 0.000837627 -0.4167396 12.0482712 0.4570315 
3 0.014450656 0.001738426 -60.3963416 0 3.1415927 
4 0.012976919 0.000968104 -15.3660889 50.958838 2.4137976 
5 0.012738226 0.000274668 -1.1273558 25.1106227 1.2089022 
6 0.012531351 0.000472737 -9.2393103 56.9428588 -2.827555 
7 0.01251901 0.000348595 -6.6352855 43.1166366 0.9396941 
8 0.012511692 0.000628562 -226.8675358 192 3.1415927 
9 0.012493123 0.000226302 -0.4553859 21.1035697 -0.5255472 
10 0.012452034 0.000123901 -4.7570831 47.119961 1.8736162 
11 0.012438872 0.000271993 -3.3331461 52.8516535 -2.15545 
12 0.012427474 0.000160298 -39.8492623 95.0812441 2.767976 
13 0.012423201 0.000100095 -9.3272313 63.55788 1.7068752 
14 0.01241987 0.000211142 -0.6739873 16.7001269 -0.261739 
15 0.012418476 7.72883E-05 -5.2187667 67.0044016 2.4191549 
16 0.012415291 6.00817E-05 -3.55913 29.0431211 1.086974 
17 0.012399607 3.99114E-05 -20.1855681 115.7552143 2.8246925 
18 0.012398451 7.39638E-05 -35.0501387 149.7066291 1.6086517 
19 0.012395666 2.02879E-05 -6.8153595 79.0506098 -3.0683173 
20 0.01239508 1.48799E-05 -11.7363145 111.638931 1.9637448 
 
As discussed in the previous section, a selection criteria was applied to the entire list of 
constituent modes obtained from Prony?s method to obtain the list of modes shown in Table 9-7. 
Table 9-7 Selected List of Constituent Modes in the Signal Detected by Prony?s Method 
Magnitude ExpDecay Frequency (Hz) Phase (rad) 
6.00817E-05 -3.55913 29.0431211 1.086974 
0.000211142 -0.6739873 16.7001269 -0.261739 
0.000271993 -3.3331461 52.8516535 -2.15545 
0.000123901 -4.7570831 47.119961 1.8736162 
0.000226302 -0.4553859 21.1035697 -0.5255472 
0.000274668 -1.1273558 25.1106227 1.2089022 
0.000837627 -0.4167396 12.0482712 0.4570315 
 
141 
 
 
 
This process was repeated for a set of 14 sliding windows. The starting and ending 
position of each window was incremented by 0.05 seconds from one window to the next, with 
the final window beginning at 6.15 seconds and ending at 7.15 seconds. Each time Prony?s 
method was applied to extract the constituent modes, the resulting output list was pruned 
according to parameters discussed in the previous section and stored in a separate list. Finally, 
the output list of 14 data windows were all combined together. As before, the exponential decay 
and frequency parameters were chosen and inserted into a scatter plot as shown in Figure 9-12 
 
 
Figure 9-12 Scatter plot of Modal Components with Four Major Clusters Demarcated 
 
142 
 
Again, a histogram was plotted to identify the most frequent modes. Just as before, it was 
observed that there are four major clusters of modes as identified in Figure 9-13. It appears that 
there are four dominant modes with frequencies: 12 Hz, 16.7 Hz, 20.8 Hz and 25 Hz. in these 
clusters.  
 
 
 
Figure 9-13 Histogram of Modal Components in Data Set 2 with Bin Width = 1 Hz 
 
9.4 Comparison of Results and Analysis 
Shown in the first two columns of table 9-8 are the lists of the four dominant modes 
extracted from the two given measurements. It appears that the results of both cases match very 
closely. Shown in the third column are the frequencies provided by the plant owner. As shown in 
0
2
4
6
8
10
12
14
16
18
1.5 3.5 5.5 7.5 9.5 11.5 13.5 15.5 17.5 19.5 21.5 23.5 25.5 27.5 29.5 31.5 33.5 35.5 37.5 39.5 41.5 43.5 45.5 47.5 49.5 51.5 53.5 55.5 57.5 59.5 61.5
Nu
mb
er 
Of
 Occ
ur
ren
ce
s 
Frequency (Hz) 
Histogram of Data Set 2 Modal Components 
143 
 
this table, the maximum percentage difference is less than 8 percent and the rest of the 
difference-percentages are under 4 percent. However, it is to be noted that these are differences, 
and not necessarily errors. There could be many reasons for these differences, including  
different operating conditions, i.e., the plant owner may have derived these figures based on 
design data.  
 
Table 9-8 Results of Modal Decomposition using Prony?s Analysis on the two data sets 
Frequencies of Major 
Modal Constituents of 
Data Set 1 Signal 
Frequencies of Major 
Modal Constituents of 
Data Set 2 Signal 
Frequencies of the 
Shaft System 
Provided by the Plant 
Owner 
Percentage 
Difference Between 
the Estimated 
(average) and Owner 
Specified 
Frequencies 
12 Hz 12 Hz 13 Hz 7.69% 
16.3 Hz 16.7 Hz 17 Hz 2.94% 
20.7 Hz 20.8 Hz 21 Hz 1.19% 
25 Hz 25 Hz 26 Hz 3.84% 
 
This measurement data was obtained from a unit that had the STGS configuration shown 
in Figure 9-14, consisting of six rotor sections. This implies that this shaft system has five natural 
frequencies. However, we have found only 4 of the five frequencies: 12 Hz, 16.7 Hz, 20.8 Hz 
and 25 Hz. Prony?s method was unable to detect and extract the fifth frequency. One of the 
reasons could be that there were too many non-linearities in the signal, in addition to noise in 
measurements, precluding Prony?s method from detecting the fifth frequency. Another reason 
could be that the magnitude of disturbance was insufficient to excite the torsional modes 
significantly. However, a more probable reason could be the eigen behavior of the STGS system. 
Recall that all the analyses performed in Chapter 7 using right eigen vectors, participation factors 
and impulse response test, indicated that for typical STG shaft systems, the highest frequency 
144 
 
mode participates extremely weakly in the generator-section velocity state signal. This small 
magnitude signal, after passing through the generator stator & rotor circuits and the air-gap flux, 
gets even more attenuated, and thus appears indistinguishably small in the instantaneous power 
signal, p(t). The results obtained in this chapter appear to be consistent with the results of the 
analyses performed in Chapter 7 and lead to the hypothesis that HP rotor section velocity signal 
is a more superior candidate, compared to instantaneous power signal, for extracting the torsional 
modes of the steam turbine generator shaft system.       
 
 
 
 
Figure 9-14 Steam Turbine Generator Shaft System Configuration of the Studied System 
 
9.5 Conclusions 
 In this chapter, Prony?s method was applied to extract the underlying modal constituents 
of actual instantaneous power measurements, p(t), of an actual STG unit. In the analyses of two 
separate p(t) signals of the same unit, it was found that there were four dominant modal 
components in the signal related to the STGS system. These frequencies matched closely 
(difference-percentage < 8%) those that were provided by the plant owner. It was known that the 
STGS system had six rotor sections, implying it had five natural frequencies. This meant that this 
procedure was not able to identify the fifth natural frequency of the STGS system. However, 
according to the analysis in Chapter 7 regarding eigen properties of typical STGS systems, this 
was to be expected. In that chapter it was found that the highest frequency mode always 
Exc Gen LP2 LP1 IP HP 
TE ? 
145 
 
participated very insignificantly in the generator velocity state, precluding it from being detected 
by Prony?s method. This again emphasizes the previous finding that to estimate the shaft spring 
constants, it is highly recommended to use the HP rotor section velocity signal for modal 
decomposition.       
 
 
 
 
 
 
 
  
 
CHAPTER 10 
CONCLUSIONS 
Over the past half century, power systems have grown tremendously which brought forth 
new problems in regard to shaft systems of steam turbine generator (STG) power plants ( which 
contribute greater than 60 percent of generating capacity in the U.S.A.) interacting with long 
transmission lines consisting of series capacitors, causing potentially irreplaceable damage to 
these shaft systems. Current work presented here adds to a growing body of solutions and 
preventive measures developed by power industry professionals to address this problem by 
proposing a new methodology to estimate the shaft spring constants non-intrusively using 
measurement signals from the shaft systems.  
This work uses Prony?s method to identify the modal components in a given signal as it 
has many advantages over the FFT. The effect of various parameters of Prony?s method - 
number of modes, number of samples and window position - have been thoroughly investigated. 
Also, an algorithm has been developed to reorder the resulting output modes of Prony?s method 
to identify the power system modes.   
In Chapters 6 and 7, a thorough analysis of four typical STG shaft systems presented in 
the literature using eigen analysis concepts and metrics, such a right eigen vectors and 
participation factors, was performed. Time domain simulations of impulse responses of the shaft 
systems have also been performed to further understand these LTI shaft systems? behavior. The 
purpose of this has been to bring forth the unique properties of STG shaft systems, with regards 
to the manner in which the modal components show up in different rotor velocity states. It has 
147 
 
been observed that high frequency modal components show up with lesser magnitudes in the 
generator and adjacent rotors section velocities relative to the high pressure (HP) turbine and 
adjacent rotor section velocities. It could be because the generator and adjacent rotors are 
considerably heavier than the HP and connected rotor sections. It has also been observed that all 
the modes show up well in the HP and connected rotor sections. This allows us options of 
choosing a proper rotor velocity signal to identify the modes. So far, the most suitable candidate 
for extracting modal components has been the HP rotor section velocity signal. 
In Chapter 8 we have tested, using ATP simulation software, a fully connected IEEE-
FBM model consisting of (1) an STG shaft system, (2) a generator d-q-o model and (3) an 
electrical network consisting of series capacitor, transmission line and infinite voltage source at 
the other end. After simulating a 20 millisecond ground fault disturbance on the line at the 
generator terminals, the resultant ring-down generator power output signal, p(t), and HP rotor 
velocity signal were analyzed using Prony?s method. The results of the modal decomposition 
revealed that while four modes showed in the generator signal, all five modes show up in the HP 
velocity signal. Using the frequency values and specified rotor moment-of-inertia values, rotor 
shaft spring values were estimated using a non-linear Newton-Raphson based least squares 
algorithm, thus validating the developed algorithm. 
An analysis of instantaneous power measurements, p(t), obtained from an STG unit, was 
presented in Chapter 9. The shaft system configuration consisted of six rotor sections, implying 
that it would have five natural frequencies. Prony's Method was applied, using the methodologies 
described in Chapters 3 and 5, to extract the underlying modal components. Four frequencies 
were identified from the list of the Prony?s output modes by applying a simple selection criteria.  
Repeating this process on two separate data sets yielded a set of four frequencies that matched 
148 
 
closely with the manufacturer specified frequencies. However, it is be noted that the shaft system 
had six rotor sections, implying that it had five natural frequencies. Because the identification 
methodology was able to identify only four frequencies, it meant that it was unsuccessful in 
identifying the fifth frequency. But, this is not particularly surprising according to the analysis 
performed in Chapters 6 & 7. Recall that the mode corresponding to the highest natural 
frequency always participated very weakly in the generator and exciter rotor sections. It can be 
hypothesized that, in this case as well, the highest mode of this shaft system is participating very 
weakly in the generator section, thereby precluding it from being identified strongly in the 
instantaneous power, p(t), measurement signals.  
In conclusion, an algorithm to estimate shaft spring constants non-intrusively using STG 
shaft measurements and an algorithm to reorder Prony?s method?s output modes have been 
presented in this work. Another significant contribution of this work has been in determining, 
after a thorough investigation of eigen properties of several steam turbine generator shaft 
systems, which of the non-intrusively available measurement signals is the best candidate to be 
used for extracting the shaft torsional modes. Taken together, these advances aim to assist the 
power engineer in designing better preventive measures and solutions to the complex and 
multifaceted subsynchronous resonance problem in power systems.  
 
 
 
 
 
 
149 
 
 
 
 
BIBLIOGRAPHY 
 
[1] IEEE Subsynchronous Resonance Working Group, ?Terms, Definitions And 
Symbols For Subsynchronous Oscillations,? IEEE Transactions on Power 
Apparatus and Systems, Vol. PAS-104, No. 6, June 1985. 
 
[2] M. C. Hall and D. A. Hodges, ?Experience with 500 kV Subsynchronous 
Resonance and Resulting Turbine Generator Shaft Damage at Mohave Generating 
Station,? IEEE Publication 76 CH1066-PWR. New York IEEE Press, 1976. pp 
22-29. 
 
[3] IEEE Subsynchronous Resonance Working Group, ?Proposed Terms And 
Definitions For Subsynchronous Oscillations,? IEEE Transactions on Power 
Apparatus and Systems, Vol. PAS-99, No.2 March/April 1980. 
 
[4] IEEE Subsynchronous Resonance Task Force of the Dynamic System 
Performance Working Group, Power System Engineering Committee,  ?First 
benchmark model for computer simulation of subsynchronous resonance,? IEEE 
Transactions on Power Apparatus  and Systems, Vol. 96, No. 5, pp. 1565-1572, 
September 1977. 
 
[5] IEEE Subsynchronous Resonance Task Force of the Dynamic System 
Performance Working Group, Power System Engineering Committee,  ?Second 
benchmark model for computer simulation of subsynchronous resonance,? IEEE 
Transactions on Power Apparatus  and Systems, Vol. PAS-104, No. 5, pp. 1057-
1066, May 1985.  
 
[6] Butler, J.W.; Concordia, C., "Analysis of Series Capacitors Application Problem", 
AIEE Transactions, Vol. 56, pp. 975-988, 1937. 
 
[7] Carter, G.K.; Concordia, C.;, "Negative Damping of Electrical Machinery", AIEE 
Transactions, Vol. 60, pp. 116-119, March 1941. 
 
 [8] Ramey, D.G.; Sismour, A.C.; Kung, G.C.; , ?Important Parameters in Considering 
Transient Torques on Turbine-Generator Shaft Systems,?  IEEE Transactions on 
Power Apparatus and Systems, vol.PAS-99, no.1, pp.311-317, Jan. 1980. 
 
150 
 
 [9] ?Introduction to numerical analysis,? by Hildebrand, F. B.  New York, McGraw-
Hill, pp. 457-462, 1956.  
 
 [10] Hauer, J.F.; Demeure, C.J.; Scharf, L.L.; , "Initial results in Prony analysis of 
power system response signals ," IEEE Transactions on Power Systems, vol.5, 
no.1, pp.80-89, Feb 1990. 
 
[11] Kay, S.M.; Marple, S.L., Jr.; , ?Spectrum analysis?A modern perspective,? 
Proceedings of the IEEE , vol.69, no.11, pp. 1380- 1419, Nov. 1981. 
 
 [12] G.R.B. Prony, ?Essai experimental et analytic ...,? J. l'Ecole Polytech. (Paris), vol. 
1, pp. 24-76, 1795.   
 
[13] Kumaresan, R.; Feng, Y., ?FIR Prefiltering improves Prony?s method,? IEEE 
Transactions on Signal Processing, Vol. 39, No. 3, pp. 736 ? 741, March 1991. 
 
[14] Pitstick, G.M.; Cruz, J.R.; Mulholland, R.J., ?A novel interpretation of Prony's 
method,? Proceedings of the IEEE, Vol. 76, No. 8, pp. 1052 - 1053, August 1988. 
 
[15]  ?Subsynchronous Resonance in Power Systems,? by Anderson, P.; Agrawal, B.; 
Ness, J.;  IEEE Press, 1990. 
 
[16] ?Power System Stability and Control,? by Kundur, P.  McGraw-Hill, 1993.  
 
[17] ?Modern Control Theory,? by Brogan, W.L.; Prentice Hall, 1990. 
 
[18] Kumaresan, R.; Tufts, D.W.; Scharf, L.L.; , "A Prony method for noisy data: 
Choosing the signal components and selecting the order in exponential signal 
models," Proceedings of the IEEE , vol.72, no.2, pp. 230- 233, Feb. 1984. 
 
[19] Zivanovic, R.; Schegner, P.; , "Pre-filtering improves Prony analysis of 
disturbance records," Developments in Power System Protection, 2004. Eighth 
IEE International Conference on , vol.2, no., pp. 780- 783 Vol.2, 5-8 April 2004. 
 
[20] J. V. Milanovic, ?The influence of shaft spring constant uncertainty on torsional 
modes of turbogenerator,? IEEE Transactions on Energy Conversion, Vol. 13, 
No. 2, June 1998. 
 
[21]  Subsynchronous Resonance Working Group of the System Dynamic Performance 
Subcommittee, ?Reader's Guide To Subsynchronous Resonance,? IEEE 
Transactions on Power Systems. Vol. 7, No. 1, February 1992. 
 
[22] Katz, E.; Tang, J.; Bowler, C.E.J.; Agrawal, B.L.; Farmer, R.G.; Demcko, J.A.; 
Selin, D.A.; Cruz, C.; Sims, J.; , "Comparison of SSR calculations and test 
results," IEEE Transactions on Power Systems, vol.4, no.1, pp.336-344, Feb 
1989. 
151 
 
 
[23] Members of the Bibliography Group of the Subsynchronous Resonance Task 
Force of the Dynamic Performance Working Group of the IEEE Power System 
Engineering Committee, ?A bibliography for the study of subsynchronous 
resonance between rotating machines and power systems,? Power Apparatus and 
Systems, IEEE Transactions on , vol.95, no.1, pp. 216- 218, Jan. 1976. 
 
[24] Bibliography Group of the Subsynchronous Resonance Task Force of the 
Dynamic Performance Working Group of the IEEE Power Systems Engineering 
Committee.; , ?First Supplement to a Bibliography for the Study of 
Subsynchronous Resonance Between Rotating Machines and Power Systems,? 
IEEE Transactions on Power Apparatus and Systems, vol.PAS-98, no.6, pp.1872-
1875, Nov. 1979. 
 
[25] Bibliography Task Force of the Subsynchronous Resonance Working Group of 
the Dynamic System Performance Subcommittee of the IEEE Power Systems 
Engineering Committee, ?Second Supplement to A Bibliography For the Study 
Of Subsynchronous Resonance Between Rotating Machines And Power 
Systems,? IEEE Transactions on Power Apparatus and Systems, vol.PAS-104, 
no.2, pp.321-327, Feb. 1985. 
 
[26] IEEE Subsynchronous Resonance Working Group, ?Third Supplement to a 
Bibliography for the Study of Subsynchronous Resonance Between Rotating 
Machines and Power Systems,? IEEE Transactions on Power Systems, vol.6, 
no.2, pp.830-834, May 1991. 
 
[27] Iravani, M.R.; Agrawal, B.L.; Baker, D.H.; Bowler, C.E.J.; Farmer, R.G.; Hedin, 
R.A.; Larsen, E.H.V.; Tang, J.F. , ?Fourth Supplement to a Bibliography for the 
Study of Subsynchronous Resonance Between Rotating Machines and Power 
Systems,? IEEE Transactions on Power Systems, vol.12, no.3, pp.1276-1282, 
Aug 1997. 
 
[28] G. C. Verghese, I.J. Perez-Arriaga, and F.C. Schweppe, , ?Selective Modal 
Analysis with Applications to Electric Power Systems, PART I: Heuristic 
Introduction, Part II: The Dynamic Stability Problem,? IEEE Transactions on 
Power Apparatus and Systems, vol.PAS-101, No. 9, pp.3117-3134,  September 
1982.