Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Modeling alternating current rotating electrical machines using constant-parameter RL-branch interfacing… Chapariha, Mehrdad 2013

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
24-ubc_2014_spring_chapariha_mehrdad.pdf [ 6.4MB ]
Metadata
JSON: 24-1.0165681.json
JSON-LD: 24-1.0165681-ld.json
RDF/XML (Pretty): 24-1.0165681-rdf.xml
RDF/JSON: 24-1.0165681-rdf.json
Turtle: 24-1.0165681-turtle.txt
N-Triples: 24-1.0165681-rdf-ntriples.txt
Original Record: 24-1.0165681-source.json
Full Text
24-1.0165681-fulltext.txt
Citation
24-1.0165681.ris

Full Text

  MODELING ALTERNATING CURRENT ROTATING ELECTRICAL MACHINES USING CONSTANT-PARAMETER RL-BRANCH INTERFACING CIRCUITS   by  Mehrdad Chapariha  B.Sc. The Isfahan University of Technology, 2006 M.Sc. The Isfahan University of Technology, 2009    A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF  DOCTOR OF PHILOSOPHY  in  The Faculty of Graduate and Postdoctoral Studies  (Electrical and Computer Engineering)    THE UNIVERSITY OF BRITISH COLUMBIA  (Vancouver)  November 2013   ? Mehrdad Chapariha, 2013  ii  Abstract Transient simulation programs are used extensively for modeling and simulation of various electrical power and energy systems that include rotating alternating current machines as generators and motors. In simulation programs, traditionally, the machine models are expressed in qd-coordinates (rotational reference frame) and transformed variables, and the power networks are modeled in abc-phase coordinates (physical variables), which represents an interfacing problem. It has been shown in the literature that the method of interfacing machine models and the electric network models plays an important role in numerical accuracy and computational performance of the overall simulation.  This research considers the state-variable-based simulation programs and proposes a unified constant-parameter decoupled RL-branch circuit in abc-phase coordinates (with optional zero-sequence). The proposed circuits are based on voltage-behind-reactance (VBR) formulation and can be used for interfacing both induction and synchronous machine models. The new models achieve a direct and explicit interface with arbitrary external electrical networks, which results in many computational advantages. Extensive computer studies are presented to verify the proposed models and to demonstrate their implementation in several commonly-used simulation programs. The new models are shown to offer significant improvements in accuracy and numerical efficiency over the existing state-of-the-art models due to their direct interface. It is further envisioned that the proposed models will receive a wide acceptance in research community and simulation software industry, and may enable the next generation of power systems simulation tools. iii  Preface Some of the research results presented in this thesis have been already published as journal articles, conference proceedings and/or submitted for peer review. In all publications, I was responsible for developing the mathematical formulations, implementing the models, conducting the simulations, compiling the results and conclusions, as well as preparing majority of the manuscripts. My research supervisor, Dr. Juri Jatskevich, has provided supervisory comments and corrections during the process of conducting the studies and writing the manuscripts. My research co-supervisor, Dr. Hermann W. Dommel, also has provided supervisory feedback, comments, and corrections, for the research and written manuscripts. The contributions of other co-authors of published and submitted papers are explained below.  A version of Chapter 2 and 3 has been published. M. Chapariha, L. Wang, J. Jatskevich, H. W. Dommel, and S. D. Pekarek, ?Constant-Parameter RL-Branch Equivalent Circuit for Interfacing AC Machine Models in State-Variable-Based Simulation Packages," IEEE Transactions on Energy Conversion, vol. 27, no. 3, pp. 634?645, September 2012. Dr. L. Wang and Dr. S. D. Pekarek have provided comments and constructive feedback about the proposed interfacing circuit which is a generalization and improvement of their previous research work.  A part of Chapter 2 was presented at a conference. M. Chapariha, F. Therrien, J. Jatskevich, and H. W. Dommel, ?Implementation of Induction Machine VBR Model with Optional Zero-Sequence in SimPowerSystems, ASMG, and PLECS Toolboxes?, In proceedings of iv  International Conference on Power Systems Transients (IPST 2013), Vancouver, BC, Canada, 18?20 July 2013. F. Therrien has provided comments, suggestions, and constructive feedback. A part of Chapter 3 was presented at a conference. M. Chapariha, F. Therrien, J. Jatskevich, and H. W. Dommel, ?Implementation of Constant-Parameter Directly-Interfaced VBR Synchronous Machine Models in SimPowerSystems, ASMG, and PLECS Toolboxes?, In proceedings of Power and Energy Society General Meeting, 2013 IEEE, Vancouver, BC, Canada, 21?25 July 2013. F. Therrien has provided comments, suggestions, and constructive feedback and revised the manuscript. A version of Chapter 4 has been published. M. Chapariha, F. Therrien, J. Jatskevich, and H. W. Dommel, "Explicit Formulations for Constant-Parameter Voltage-Behind-Reactance Interfacing of Synchronous Machine Models," IEEE Transactions on Energy Conversion, vol. 28, no. 4, pp. 1053-1063, December 2013. F. Therrien has provided comments, suggestions, and constructive feedback and revised the manuscript. He has also helped me in the selection and design of the case-study system and the numerical analysis of the results. A version of Chapter 5 is has been submitted for peer review. M. Chapariha, F. Therrien, J. Jatskevich, and H. W. Dommel, "Constant-Parameter Circuit-Based Models of Synchronous Machine". F. Therrien has provided comments, suggestions, and constructive feedback and revised the manuscript. He also helped me in selection and design of the case studies. v  Table of Contents Abstract .................................................................................................................................................... ii Preface ..................................................................................................................................................... iii Table of Contents .................................................................................................................................... v List of Tables ............................................................................................................................................ x List of Figures ........................................................................................................................................ xi List of Abbreviations ........................................................................................................................ xvii Nomenclature ................................................................................................................................... xviii Acknowledgments............................................................................................................................. xxii Dedication ........................................................................................................................................... xxiv CHAPTER 1: INTRODUCTION ............................................................................................................ 1 1.1 Motivation .................................................................................................................................................... 1 1.2 Background ................................................................................................................................................. 3 1.2.1 Power systems transient simulation tools ................................................................................. 3 1.2.2 Mathematical models of ac rotating machines....................................................................... 4 1.2.3 Voltage-behind-reactance models ................................................................................................ 7 1.2.3.1 Induction machine models ....................................................................................... 7 1.2.3.2 Synchronous machine models ................................................................................ 8 1.3 State-of-the-Art Research ...................................................................................................................... 9 1.4 Research Objectives and Anticipated Impacts ........................................................................... 10 CHAPTER 2: CIRCUIT INTERFACING OF INDUCTION MACHINE MODELS ........................ 15 2.1 Induction Machine Modeling ............................................................................................................. 15 vi  2.1.1 Coupled-circuit phase-domain model ...................................................................................... 15 2.1.2 The qd0 model in arbitrary reference-frame ........................................................................ 19 2.2 Voltage-Behind-Reactance Formulations ..................................................................................... 21 2.2.1 Explicit formulation with zero-sequence branch in the interfacing circuit ............ 27 2.3 Interfacing of Induction Machine Models in State-Variable-Based Simulation Programs ................................................................................................................................................. 29 2.3.1 Examples of implementation in SimPowerSystems, ASMG, and PLECS toolboxes 30 2.4 Computer Studies .................................................................................................................................. 32 CHAPTER 3: CIRCUIT INTERFACING OF SYNCHRONOUS MACHINE MODELS ................. 40 3.1 Synchronous Machine Modeling ...................................................................................................... 41 3.1.1 Coupled-circuit phase-domain model ...................................................................................... 41 3.1.2 The qd0 model in rotor reference-frame ................................................................................ 44 3.2 Voltage-Behind-Reactance Formulations ..................................................................................... 46 3.2.1 Variable-parameter formulation ............................................................................................... 47 3.2.2 Standard state-space form for the rotor subsystem .......................................................... 50 3.3 Approximation of Dynamic Saliency in VBR Models ................................................................ 53 3.3.1 Additional winding method .......................................................................................................... 54 3.3.2 Singular perturbation method .................................................................................................... 56 3.4 Generalized Constant-Parameter VBR Formulation ................................................................ 60 3.4.1 Possible explicit and implicit formulations ............................................................................ 62 3.5 Interfacing of Synchronous Machine Models in State-Variable-Based Simulation Programs ................................................................................................................................................. 63 vii  3.6 Implementation of Synchronous Machine VBR Formulations in Simulation Programs ..................................................................................................................................................................... 66 3.6.1 Examples of implementation in SimPowerSystems, ASMG, and PLECS toolboxes 67 3.7 Computer Studies .................................................................................................................................. 69 3.7.1 Verification of the implicit VBR formulation ........................................................................ 69 3.7.1.1 Constant time-step simulation study ................................................................ 70 3.7.1.2 Variable time-step simulation study ................................................................. 71 3.7.2 Accuracy of the approximated constant-parameter VBR formulation ..................... 72 3.7.3 Small machine case-study ............................................................................................................. 75 3.7.4 Large machine case-study ............................................................................................................. 80 CHAPTER 4: NUMERICAL METHODS TO ACHIEVE CONSTANT-PARAMETER VBR FORMULATIONS ......................................................................................................... 86 4.1 Method of Using Current Derivatives............................................................................................. 86 4.2 Method of Using Algebraic Feed-through..................................................................................... 88 4.3 Algebraic Equivalence of Implicit VBR Formulations .............................................................. 90 4.4 Numerically Efficient Explicit Implementation .......................................................................... 93 4.4.1 Approximation of current derivative ....................................................................................... 93 4.4.2 Relaxation of algebraic loop ........................................................................................................ 95 4.5 Summary of Approximation Techniques ...................................................................................... 98 4.6 Computer Studies .................................................................................................................................. 99 4.6.1 Continuous-time approximation techniques ...................................................................... 102 4.6.2 Discrete-time approximation techniques ............................................................................ 106 viii  4.6.3 Comparison of models and approximation techniques ................................................. 108 CHAPTER 5: DIRECT INTERFACING OF SYNCHRONOUS MACHINE MODELS FROM STATOR AND ROTOR TERMINALS .................................................................... 112 5.1 Rotor and Stator Interfacing of Synchronous Machine Models ......................................... 113 5.2 All-Circuit Formulation ...................................................................................................................... 116 5.2.1 Stator voltage equations ............................................................................................................ 116 5.2.2 Rotor voltage equations .............................................................................................................. 119 5.3 Stator-and-Field-Circuit Formulation .......................................................................................... 122 5.3.1 Stator voltage equations ............................................................................................................ 122 5.3.2 Rotor voltage equations .............................................................................................................. 126 5.4 Numerically Efficient Explicit Implementation ........................................................................ 129 5.5 Computer Studies ................................................................................................................................ 130 5.5.1 Single-phase-to-ground fault in the network .................................................................... 134 5.5.2 Diode failure in the exciter system ......................................................................................... 139 CHAPTER 6: SUMMARY OF CONTRIBUTIONS AND FUTURE WORK ............................... 144 6.1 Conclusions and Contributions ...................................................................................................... 144 6.1.1 Objective 1 ......................................................................................................................................... 144 6.1.2 Objective 2 ......................................................................................................................................... 145 6.1.3 Objective 3 ......................................................................................................................................... 145 6.1.4 Objective 4 ......................................................................................................................................... 146 6.2 Potential Impacts of Contributions ............................................................................................... 147 6.3 Future Work ........................................................................................................................................... 148 ix  6.3.1 Doubly-fed induction machine model with direct interface ........................................ 148 6.3.2 Inclusion of magnetic saturation ............................................................................................ 148 6.3.3 Application to dynamic phasor solution and shifted frequency analysis .............. 149 References .......................................................................................................................................... 150 Appendices ......................................................................................................................................... 157 Appendix A: Parameters of the Induction Machine Model in Section ?2.4 ............................. 157 Appendix B: MATLAB Ordinary Differential Equation (ODE) Solvers ................................... 157 Appendix C: Parameters of the Synchronous Machine in Subsections ?3.7.1 to ?3.7.3 ....... 159 Appendix D: Parameters of the Power System in Subsection ?3.7.4 ......................................... 159 Appendix E: Parameters of the Power System in Section ?4.6 .................................................... 159 Appendix F: Parameters of the Power System in Section ?5.5 .................................................... 160  x  List of Tables Table ?2?1  Comparison of Numerical Efficiency of VBR Formulations for Single-Phase Fault Study ....................................................................................................................................................................... 39 Table ?3?1  Simulation Efficiency for Single-Phase Fault Study for a Synchronous Machine Models .................................................................................................................................................................... 80 Table ?3?2  Numerical Efficiency of CP-VBR Models Versus qd0 with Snubber for the Large Machine Study ..................................................................................................................................................... 85 Table ?4?1  Summary of Approximation Techniques to Achieve the Interfacing Circuit of Figure ?2-2 .............................................................................................................................................................. 99 Table ?4?2  Comparison of Continuous-Time Approximation Methods (With Added Eigenvalue of 1,000) ....................................................................................................................................... 103 Table ?4?3  Comparison of Discrete-Time Approximation Methods ............................................. 107 Table ?4?4  Eigenvalues of the System when using Different Approximation Techniques .. 110 Table ?5?1  Classification of Synchronous Machine Formulations based on Interfacing with External Inductive Networks....................................................................................................................... 114 Table ?5?2  Numerical Performance of the Models for the Single-Phase-to-Ground Fault Study ..................................................................................................................................................................... 138 Table ?5?3  Numerical Performance of the Models for the Exciter Diode Failure Study ....... 143  xi  List of Figures Figure ?2-1  Magnetically coupled circuit model of induction machine. ........................................ 16 Figure ?2-2  General interfacing circuit for induction machine VBR formulation. ..................... 28 Figure ?2-3  Interfacing of induction machine qd0 models with external electrical network. ................................................................................................................................................................................... 30 Figure ?2-4  Interfacing of induction machine VBR models with external electrical network. ................................................................................................................................................................................... 30 Figure ?2-5  Example of implementation of the proposed induction machine VBR model in Simulink using the SimPowerSystems toolbox. ..................................................................................... 31 Figure ?2-6  Example of implementation of the proposed induction machine VBR model equivalent interfacing circuit with external electrical network inside the PLECS Circuit block. ....................................................................................................................................................................... 32 Figure ?2-7  Induction generator connected to the Th?venin equivalent circuit of a network for the single-phase-to-ground fault transient study. ......................................................................... 33 Figure ?2-8  Simulation results for the single-phase-to-ground fault study. From top: source voltage, source current, machine neutral current, and electromechanical torque. ................. 36 Figure ?2-9  Detailed view of source current ics in steady-state excerpt from Figure ?2-8. ...... 37 Figure ?2-10  Detailed view of source current ics during transient excerpt from Figure ?2-8. 37 Figure ?2-11  Detailed view of electromagnetic torque Te during transient excerpt from Figure ?2-8. ............................................................................................................................................................. 38 Figure ?3-1  Magnetically coupled circuit model of synchronous machine. ................................. 42 xii  Figure ?3-2  Interfacing synchronous machine classical VBR models with external electrical network-circuit with algebraic feed-through of current. ................................................................... 65 Figure ?3-3  Interfacing synchronous machine constant-parameter VBR models with external electrical network-circuit with algebraic feed-through of current and possible algebraic loop for voltage. .............................................................................................................................. 65 Figure ?3-4  Example of implementation of the synchronous machine CP-VBR model in SimPowerSystems toolbox. ............................................................................................................................ 68 Figure ?3-5  Example of implementation of the synchronous machine CP-VBR model in PLECS toolbox. .................................................................................................................................................... 68 Figure ?3-6  Transient response to the three phase fault predicted by the variable impedance and the constant-parameters implicit VBR models. ..................................................... 71 Figure ?3-7  Detailed view of current iqs from the three phase fault study shown in Figure ?3-6; (a) constant time-step solution; and (b) variable time-step solution. ................................. 72 Figure ?3-8  Current iqs from the three phase fault study predicted by the approximate constant-parameters VBR formulation: (a) overall transient; and (b) magnified view of the window in part (a). ............................................................................................................................................ 73 Figure ?3-9  The effect of additional winding resistance on: (a) approximate model error in iqs; and (b) the system largest eigenvalue magnitude. ......................................................................... 75 Figure ?3-10  Synchronous generator connected to an inductive network. The snubbers are required only for the classical qd0 model. ............................................................................................... 76 xiii  Figure ?3-11  Transient response to a single phase to ground fault predicted by qd0 model with snubber versus approximated constant-parameters VBR and exact implicit VBR models. ................................................................................................................................................................... 78 Figure ?3-12  Detailed view of current ibs in steady-state from Figure ?3-11. ............................... 79 Figure ?3-13  Detailed view of current ics during the single phase fault from Figure ?3-11. .... 79 Figure ?3-14  A large synchronous generator connected to an inductive network. The snubbers are required only for the classical qd0 model. .................................................................... 81 Figure ?3-15  Simulation results for large machine single-phase-to-ground fault study as predicted by CP-VBR models and the conventional qd0 model. ...................................................... 83 Figure ?3-16  Detailed view of phase current in steady-state shown in Figure ?3-15. ............... 84 Figure ?3-17  Detailed view of phase current during transient shown in Figure ?3-15. ........... 84 Figure ?3-18  Detailed view of electromagnetic torque shown in Figure ?3-15. .......................... 84 Figure ?4-1  Implementation of explicit constant-parameter VBR model using filter Hi in continuous-time (high-pass) or discrete-time (backward difference) to approximate the current derivative. ............................................................................................................................................. 95 Figure ?4-2  Implementation of explicit constant-parameter VBR model using filter Hv in continuous time (low-pass) or discrete time (zero- or first-order-hold) to break the algebraic loop. ..................................................................................................................................................... 96 Figure ?4-3  Approximation of voltage vqs using zero-order-hold (delay) or first-order-hold (linear prediction). ............................................................................................................................................ 97 xiv  Figure ?4-4  Test system consisting of a grounded steam turbine generator connected to the network via a unit transformer. A single phase-to-ground fault is applied at the machine terminals. ............................................................................................................................................................ 101 Figure ?4-5  Transient response to a single phase-to-ground fault predicted by continuous-time approximated models. From top to bottom: terminal voltages vabcs, phase a stator current ias, neutral grounding current ing, field current ifd, and electromagnetic torque Te.103 Figure ?4-6  Detailed view of current ias shown in Figure ?4-5 for continuous-time approximated models..................................................................................................................................... 104 Figure ?4-7  Detailed view of current ifd shown in Figure ?4-5 for continuous-time approximated models..................................................................................................................................... 104 Figure ?4-8  Detailed view of electromagnetic torque Te shown in Figure ?4-5 for continuous-time approximated models. ......................................................................................................................... 105 Figure ?4-9  Cumulative error in predicting currents iabcs versus the magnitude of the added eigenvalue for continuous-time approximated models. ................................................................... 105 Figure ?4-10  Detailed view of current ias shown in Figure ?4-5 for discrete-time approximated models..................................................................................................................................... 106 Figure ?4-11  Detailed view of current ifd shown in Figure ?4-5 for discrete-time approximated models..................................................................................................................................... 107 Figure ?4-12  Detailed view of electromagnetic torque Te shown in Figure ?4-5 for discrete-time approximated models. ......................................................................................................................... 107 Figure ?4-13  Cumulative error in predicting currents iabcs versus maximum time-step size for discrete-time approximated models. ................................................................................................ 108 xv  Figure ?5-1  All-circuit constant-parameter VBR synchronous machine model (AC-CP-VBR) with direct interfacing to arbitrary external ac and dc networks. ................................................ 121 Figure ?5-2  Rotor subsystem for the stator-and-field-circuit constant-parameter voltage-behind-reactance (SFC-CP-VBR) model wherein the damper windings are represented as a state model and the field winding is made available as an interfacing circuit (the stator interfacing circuit is the same as in Figure ?5-1). .................................................................................. 128 Figure ?5-3  Algebraic loops in AC-CP-VBR resulting in an implicit formulation. The H blocks indicate where the low-pass filters may be inserted to break the algebraic loops. ............... 130 Figure ?5-4  A wye-grounded steam turbine generator with a static 12-pulse rectifier-based exciter system. The stator snubbers are required only for the classical qd0 model and the field snubber is required only for the CP-VBR model from Chapter 4. ....................................... 132 Figure ?5-5  Single-phase-to-ground fault case-study transient responses as predicted by the considered models. From top to bottom: bus voltages vabc and currents iabc, machine neutral current ing, three-phase ac exciter current iabcex, rectifier current idc, and electromagnetic torque Te. ............................................................................................................................................................. 135 Figure ?5-6  Magnified fragment from Figure ?5-5: phase current ic (in steady-state). ........... 136 Figure ?5-7  Magnified fragment from Figure ?5-5: phase current ic (during transient). ........ 136 Figure ?5-8  Magnified fragment from Figure ?5-5: ac exciter phase current icex. ...................... 137 Figure ?5-9  Magnified fragment from Figure ?5-5: rectified current idc. ....................................... 137 Figure ?5-10  Left: the diode bridge connected to the delta transformer in Figure ?5-4. Right: simulation results of the diode currents in pu: 1) after the initial failure (short-circuit) of xvi  D3; and 2) after D1 and D3 (or the corresponding protective diodes) become open due to the ensuing excessive currents. .................................................................................................................. 141 Figure ?5-11  Transient responses as predicted by the considered models for the static exciter diode failure case-study. From top to bottom: exciter delta transformer line voltages vT?, delta transformer currents iT?, wye transformer line voltages vTY, wye transformer currents iTY, rectifier output voltage vfd, and rectifier current idc. ................................................. 142 xvii  List of Abbreviations ac alternating current AC-CP-VBR all-circuit CP-VBR ASMG Automated State Model Generator (toolbox) CC-PD coupled-circuit phase-domain CPU central processing unit CP-VBR constant-parameter VBR DAE differential-algebraic equation d-axis direct axis  dc direct current  IVBR implicit VBR IM induction machine ODE ordinary differential equation pu per-unit q-axis quadrature axis  s second(s) SI the International System (of Units)  SFC-CP-VBR stator-and-field-circuit CP-VBR SM synchronous machine SPS SimPowerSystems (toolbox) UBC the University of British Columbia VBR voltage-behind-reactance xviii  Nomenclature Throughout this thesis, matrix and vector quantities are boldfaced (e.g. abcsv ), and scalar quantities are italic non-boldfaced (e.g. si0 ). All machine variables are referred to the stator side using the appropriate turns-ratio.  cos  Cosine function )diag(x  A diagonal matrix containing vector x  on the main diagonal  abce ??  Sub-transient voltage vector abce ???  Modified (approximated) sub-transient voltage vector xfde Scaled excitation voltage equal to ? ? fdfdmd vrX /  (synchronous machines) H  Inertia constant si0  Stator zero-sequence current abcsi  Stator actual current vector ngi Current of the VBR interfacing circuit zero-sequence branch  J  Moment of inertia 0L  Inductance of the zero-sequence branch of the VBR interfacing circuit DL  Inductance of the three-phase branch of the VBR interfacing circuit  lfdL Field winding leakage inductance (synchronous machines) lkdjL, Nj ?1?  The d-axis rotor damper leakage inductances (synchronous machines) xix  lkqjL, Mj ?1?  The q-axis rotor damper leakage inductances (synchronous machines) lrL  Rotor winding leakage inductances (induction machines) lsL  Stator winding leakage inductances (induction and synchronous machines) mL  Magnetizing inductance (induction machines) mdL  Magnetizing inductance of the d-axis (synchronous machines) mqL Magnetizing inductance of the q-axis (synchronous machines) p  Heaviside?s operator (differentiation with respect to time) P  Number of poles BP  Base power 0r  Resistance of the zero-sequence branch of the VBR interfacing circuit Dr  Resistance of the three-phase branch of the VBR interfacing circuit  fdr Field winding resistance (synchronous machines) gr External grounding resistance  kdjr, Nj ?1?  The d-axis rotor damper resistances (synchronous machines) kqjr, Mj ?1?  The q-axis rotor damper resistances (synchronous machines) rr  Rotor winding resistance (induction machines) sr  Stator winding resistance (induction and synchronous machines) s  Laplace variable  sin  Sine function xx  BT  Base torque  eT  Electromagnetic torque  lT  Load torque (motor)  MT  Mechanical input torque (generator) abcnv  Voltage vector from stator terminals to the point n in VBR interfacing   circuit abcsv  Stator actual voltage vector drv  Rotor voltage transformed to d-axis (induction machines) dsv  Stator voltage transformed to d-axis (induction and synchronous   machines) fdv Excitation voltage (synchronous machines) ngv Voltage of the VBR interfacing circuit zero-sequence branch  qrv Rotor voltage transformed to q-axis (induction machines) qsv Stator voltage transformed to q-axis (induction and synchronous   machines) gX External grounding reactance z  A discrete random variable  )(x?  2-norm cumulative relative error (calculated for all the points of x ) )(x?  Average 2-norm cumulative relative error of all x  elements (three  phases) xxi  ?  Flux linkage d? ??  Sub-transient flux linkages of the d-axis q? ?? Sub-transient flux linkages of the q-axis ?  Angular velocity of the rotational reference frame b?  Base angular velocity r?  Angular velocity of the rotor (electrical)    xxii  Acknowledgments I would like to express my sincere gratitude to my research supervisor, Dr. Juri Jatskevich, for his research direction and excellent vision, great inspiration, patient guidance, caring, and generous support throughout my PhD studies. I also would like to express my very great appreciation to my research co-supervisor Dr. Hermann W. Dommel for his exceptional expertise, great help and guidance, admirable manners, and continuous support. The financial support for this research was made possible through the Natural Science and Engineering Research Council (NSERC) of Canada under the Discovery Grant entitled ?Modelling and Analysis of Power Electronic and Energy Conversion Systems? and the Discovery Accelerator Supplement Grant entitled ?Enabling Next Generation of Transient Simulation Programs? lead by Dr. Juri Jatskevich as a sole principal investigator.  I would like to offer my special thanks to the members of my examining committees of the departmental and university exams, Drs. Ludovic Van Waerbeke, Jose Marti, Farrokh Sassani, Edmond Cretu, Narayan Kar, William G. Dunford, and Martin Ordonez for dedication of their valuable time and expertise and assessment of my dissertation. I would like to thank all the faculty members, instructors, and professors at the Electrical and Computer Engineering Department whom I attended their classes and lectures. Specifically I would like to thank Dr. Luis Linares, Dr. William Dunford, Mr. Nathan Ozog, and my supervisor Dr. Juri Jatskevich for the great teaching experience I gained assisting them. xxiii  I wish to acknowledge the help and support of Dr. Liwei Wang and Dr. Steven D. Pekarek. This thesis would not have been possible without their prior achievements in the area. I am also particularly grateful for the help of Francis Therrien, for perfecting this research by dedicating his time, providing valuable feedback, and careful investigations into expanding the proposed techniques and inclusion of magnetic saturation. My special thanks are extended to my colleagues and fellow graduate students in the Power and Energy Systems research group at the University of British Columbia. Particularly, I would like to thank Jaishankar Iyer, Amir Rasuli, and Sina Chiniforoosh, for all their help and support in my academic life and teaching and research activities. I also would like to thank Milad Gougani, Mohammed Talat Khouj, Mehmet Sucu, Hamid Atighechi, Kamran Tabarraee, Shahrzad Rostamirad, and all the other former and current members of our research group for their help and good memories we have together during the years at the University of British Columbia. I also thank my friends, many of whom are all around the world, for their moral support and some of the greatest moment of my life that I spend with them. The last but not the least, I would like to thank my family for their continuous support and encouragement. My greatest gratitude goes to my parents Mohammad and Ata, for their unconditional love, to my brother Alireza, for his encouragements, and to my sister Mojgan for her caring while I was thousands of miles away.   xxiv  Dedication  To My Parents Mohammad and Ata      1  CHAPTER 1:  INTRODUCTION 1.1     Motivation Synchronous and induction machines are responsible for almost all of electromechanical energy conversion in today?s world. Substantial portion of electrical energy is produced by synchronous generators run by hydro turbines, steam turbines, or diesel engines. Induction machines are the workhorse of electric power industries ?[1]. They are commonly used as motors with squirrel cage rotors in industrial and commercial applications. More recently, however, induction machines have been also used as smaller-scale generators in distributed energy systems and wind turbines in particular ?[2], ?[3]. In addition to the energy conversion, the synchronous machines have a major influence on the power systems stability, and are also used for the voltage and frequency control at the system level ?[4] ? ?[6]. Detailed three-phase transient simulation of small-to-large scale power systems is needed for various design, testing, and research purposes. Almost all types of large-and-small scale power systems today include ac rotating machines: utility grids, small community micro-grids, wind farms, vessel power systems, off-grid telecommunication sites power systems, electrical drives, etc.  2  Mathematical modeling of rotating electrical machines has been an active area of research over the century. The complex physical structures of machines and windings in relative rotation make this modeling quite involved. The mathematical formulation describing the machine electromechanical phenomena are therefore non-linear and time variant, making the machine model a bottleneck for many transient simulation programs. A change of variables to rotational qd0-reference frame is typically used to simplify the modeling and simulation of electric machines ?[7]. However, this technique has difficulties when it comes to interfacing of the machine models with external power networks which are typically represented in physical variable and abc-phase coordinates ?[8]. Improving the electric machine models and their interface to external power networks was shown to have considerable impact on the overall simulation time and accuracy ?[8]. General purpose machine models are found as built-in components in almost all simulation software packages that are used today in industry and academia. These models are usually represented in qd0-rotational reference frames, and need to be interfaced with the rest of the network represented in abc-phase coordinates, which represents an interfacing problem ?[8]. The traditional approaches used to solve the interfacing problem include using snubbers or time-delays, which lead to inaccuracies and may cause numerical stiffness and problems with convergence of solutions. The state-of-the-art solutions to solve the interfacing problem (without using three-phase snubbers or time-delays) result in models with variable parameters and/or models that are implicit and require iterative solutions. These 3  solutions require computationally expensive methods and algorithms; therefore they are not implemented in most software packages today. The preferred machine model of choice should be computationally efficient, explicit, accurate, and use basic circuit components for interfacing with remaining power networks without complicating the overall simulation program solver. Proposing and developing such models will have significant benefit for the software developers, who will be able to easily implement the models into their packages, and thousands of users, who would be able to use, implement, or modify the models. The result will be equipping the new generation of transient simulation tools with highly accurate and numerically stable models. These tools can be applied to much larger systems and achieve considerably faster simulation results without overly-constrained integration time-steps. 1.2     Background 1.2.1    Power systems transient simulation tools The methods to simulate electrical power systems can be broadly categorized into two general types:  1) Nodal-analysis-based electromagnetic transient programs (EMTP-type) ?[9], and 2) State-variable based simulation programs ?[8]. 4  In EMTP-type solution, the circuit and dynamic components are first discretized using a particular integration rule (e.g. trapezoidal), and then the resultant system of algebraic equations is solved for every time step ?[9]. EMTP type simulation programs are widely used for simulation of transients in power systems. Examples of such programs include: ATP ?[10], MicroTran ?[11], PSCAD ?[12], EMTP-RV ?[13], and PSIM ?[14]. In state-variable-based simulation programs, the entire system is represented as a state-space model. The resulting system of ordinary differential equations (OEDs) or differential algebraic equations (DAEs) is discretized at the system level and solved using either fixed or variable integration techniques. State-variable-based simulation programs are readily used in simulation of smaller power systems and electrical drive systems. These programs have the flexibility of using different ODE solvers and settings that may be used to optimize the simulation performance and time, while allowing the use of the state-space model for the system-level analysis and control designs either in time or in frequency domains. Examples of such programs include MATLAB/Simulink ?[15], ?[16] and toolboxes such as PLECS ?[17], ASMG ?[18], SimPowerSystems ?[19], as well as other programs including Modelica ?[20], acslX ?[21], EASY5 ?[22], and EUROSTAG ?[23]. In this thesis, only the state-variable-based simulation programs are considered. 1.2.2    Mathematical models of ac rotating machines Modeling of ac machines has been an active area of research for many decades. Depending on the objectives of studies and required accuracy, various models have been proposed in 5  the literature. The general purpose models, which are considered in this thesis, are based on magnetically coupled circuit representation of the machine?s windings, which leads to a relatively small number of equations ?[1]. Such models can be generally represented in either direct abc-phase coordinates using the physical variables ?[1], ?[24] ? ?[27], the so-called coupled-circuit phase-domain (CC-PD) models; or in rotating qd-coordinates using the transformed variables (the classical qd0-models) ?[1], ?[7]; or the hybrid of those known as the voltage-behind-reactance (VBR) models ?[28] ? ?[37]. The commonly used qd0 machine models offer a computationally efficient alternative to the original lumped parameter coupled-circuit models in phase domain. However, as was mentioned in section ?1.1, the qd0 models suffer from interfacing issues in connection to inductive networks. In typical state-variable-based programs, these models are represented as voltage-input current-output subsystems and are interfaced as dependent current-sources to abc phase-domain systems ?[8]. Such interfacing method requires the addition of snubber circuits when the machine model is connected directly to inductive elements. As a result, the numerical advantages of the qd0 model may deteriorate very quickly by adding the snubbers. The snubber circuits add error and can make the system numerically stiff, decreasing the numerical efficiency and stability of the entire simulation. Since many elements in power systems are inductive (transformers, transmission lines, etc.), the interfacing problem in simulation of power systems is quite common. To solve the interfacing problem, the direct implementation of machine model in physical abc-phase coordinates, the CC-PD models, may be reconsidered ?[24] ? ?[27]. Although the 6  CC-PD models offer direct interfacing with any external circuit, they are numerically inefficient due to their rotor-position-dependent (variable) and coupled inductances and poorly-scaled eigenvalues ?[28]. In general, when implemented in conventional state-variable-based simulation packages (e.g. MATLAB/Simulink ?[15], ?[16]), the variable inductances require costly reformulation (or update) of state-space equations at every time-step. At the same time, the poorly-scaled eigenvalues may also force the integration solver to use smaller time-steps ?[28]. These factors generally lead to longer simulation times. The inclusion of magnetic saturation generally increases the accuracy of the machine modes and their range of application to various transient studies and scenarios, but also comes at the cost of increased complexity of the models.  In many power systems studies, the saturation curves of many machines in the system are not available; therefore such studies are commonly carried out using the magnetically-linear machine models with sufficient level of accuracy. For these reasons, many simulation programs have the option of enabling or disabling the saturation in machine models, wherein the user can decide to use the magnetically-linear machine models if the saturation curves are not available or to achieve a faster simulation. This thesis considers only the magnetically-linear machine models and their interfacing as the primary focus of research. Development of constant-parameter saturable models with direct interface is left for the future research, and is currently undertaken by other members of our research group at the University of British Columbia (UBC).  7  1.2.3     Voltage-behind-reactance models The VBR models have been derived to have a direct interface of the stator circuits similar to the CC-PD models. In VBR models ?[28] ? ?[34], the transformation of variables to the rotational reference frame ?[1], ?[7] is used to increase numerical efficiency, while still leaving part of the model as RL circuits (instead of pure dependent current-sources). The changed variables in rotational reference frame are coupled with the circuit part of the model through dependent voltage-sources behind an impedance (equivalent resistance and inductance). Many proposed VBR models possess rotor-position-dependent variable-parameters. However, compared to the original CC-PD models, the VBR models generally have better scaled eigenvalues and are more numerically efficient ?[28] ? ?[34].  1.2.3.1    Induction machine models  For induction machine modeling, a simple VBR formulation was presented in ?[32], which has a constant and decoupled interfacing circuit. This circuit is algebraically equivalent to the original CC-PD model and the qd0 model, and it has the combined advantage of constant-parameter model with direct interface. These properties make the VBR model ?[32] superior to the qd0 model when connected to arbitrary networks in abc-physical coordinate.  In a general case, the induction machines (especially generators) may be grounded for protection and monitoring purposes ?[35]. In the VBR model ?[32] an algebraic loop is created when the machine is grounded and the zero-sequence is included. The algebraic 8  loop makes this model formulation implicit and complicates its implementation and numerical solution.  1.2.3.2    Synchronous machine models  Due to asymmetries in the rotor, the development of a constant-parameter VBR model is especially challenging for the synchronous machines with either round or salient rotor. For this reason, classical VBR models and many of the most recent ones for synchronous machines have time-variant interfacing circuits ?[28], ?[29], ?[33], ?[34]. The pioneering approach to achieve a constant parameter VBR model consists of adding a fictitious high-frequency damper winding to the model ?[30]. To minimize error, the additional winding parameters are chosen specifically as not to affect the typically low frequencies of interest in the power system under study. However, the added pole in this method may make the system numerically stiff. The high-frequency winding and the consequential error can be removed by pushing the effective frequency of the extra winding to infinity while keeping the constant-parameter property of the formulation ?[31]. It will be shown in this thesis that the final formulation is algebraically equivalent to the qd0 model; however, this process creates an algebraically implicit formulation that requires an iterative solution. It will also be shown that an iterative solution is very computationally demanding. Additionally, convergence is not guaranteed for implicit formulations ?[16]. 9  1.3     State-of-the-Art Research  Recently, there has been a considerable amount of research from various groups in the world focused on the modeling of ac machines for transient simulation programs. A number of ac machine models have been proposed by the research group at Purdue University.  S. D. Pekarek originally proposed the exact VBR model for synchronous machines (1998) ?[28]. This research group has also proposed the additional winding method for approximation of dynamic saliency (1999) ?[30] and removal of the error by singular perturbation (2002) ?[31]. D. C. Aliprantis at Iowa Sate University, USA, has made considerable contributions with introducing the arbitrary network rotor representation and inclusion of saturation in the VBR formulations (2008) ?[33]. A. C. Cramer, who is with the University of Kentucky, USA, has continued Aliprantis? work with introduction of rotor and stator VBR formulation for the first time (2012) ?[34].  ?cole Polytechnique de Montr?al in Canada has also been active in the transient modeling of electric machines leaded by U. Karaagac. In recent publications the CC-PD, qd0, and VBR models are used to increase the efficiency of EMTP-type simulations (2011) ? (2013) ?[38] ? ?[40].  The University of Glasgow and the University of Manchester are two of the center for research on modeling of electric machines in the United Kingdom. Their research on VBR models are focused on the EMTP-type solutions. They have worked on modeling the 10  internal fault in the stator of synchronous machines using the VBR representation (2009) ?[41]. In a more recent publication (2012) ?[42], they claim a new efficient implementation of machine models for EMTP-type programs. The Electric Power and Energy Systems research group at the University of British Columbia is taking an active and leading role in this direction (since 2004). A significant part of this research effort has been on advance models of electrical machines for the EMTP-type programs (2006) ? (2011) ?[36], ?[37], ?[43] ? ?[45], also considering of magnetic saturation (2010) ? (2013) ?[44] ? ?[46], and constant conductance matrix implementations for induction and synchronous machines ?[43], ?[47],  respectively.   The VBR modeling approach is getting recognition and wide acceptance. Fore example, recently, a VBR model for five-phase induction machines was presented in ?[48]. This model, which is derived for state-variable-based simulation programs, includes the effect of magnetic saturation with some approximations; however, the presented model is numerically stiff. 1.4     Research Objectives and Anticipated Impacts Although the interfacing problem is well-known and various VBR models are already developed, the qd0 models are still very commonly used for induction and synchronous machines modeling. The VBR models provide direct interface, but they may have variable parameters, contain algebraic loops, or be numerically stiff as will be shown in this thesis. 11  Although some simulation tools are designed for implementation of variable elements and solving algebraic loops ?[18], many program do not have variable impedance elements ?[19], or are not able to solve algebraic loops ?[17]. Some VBR models are recently included in programs such as PLECS ?[17]. However, most programs cannot take advantage of the VBR formulation with variable parameters and cannot solve implicitly formulated VBR models. In summary, to achieve new generation of transient simulation tools with highly accurate and numerically stable models, the objectives of this research are to develop new models and improve the interfacing of ac machine models specifically in the state-variable-based simulation programs. More specifically, the following objectives are considered: Objective 1: Elimination of algebraic loop in induction machine VBR formulation  The state-of-the-art VBR induction machine model has a simple direct interface to arbitrary networks ?[32]. Although not very common in practice, in a general case, induction machines may be grounded for monitoring or protection purposes ?[35]. As was mentioned earlier in section ?1.2.3.1, the VBR model ?[32] in such cases imposes an algebraic constraint to the formulation, and therefore may not be suitable as a general purpose universal component model for many simulation programs. Removal of this algebraic constraint is the first objective of this thesis.    12  Objective 2: Unified interfacing circuit for models of ac machines The pioneering methods to achieve constant-parameter VBR model for the synchronous machines ?[30] and ?[31] do not achieve a simple interfacing equivalent circuit, and therefore are also difficult to integrate into many simulation programs. The second objective of this research is to achieve a simple and general constant-parameter interfacing circuit for both induction and synchronous machine models. In fact, the constant-parameter equivalent circuit developed in Objective 1 can be used as a reference to develop new formulations for the synchronous machine model that indeed can have the same equivalent circuit for interfacing of the machine's stator terminals to the rest of the power network.  Objective 3: Numerically efficient approximations for constant-parameter models The state-of-the-art existing methods to achieve constant-parameter interface for synchronous machines ?[30], ?[31] are based on modification of the original coupled-circuit model or result in implicit formulation. Proper selection of the parameters of the modified circuits and/or dealing with implicit algebraic constraints represents additional undesirable challenges. Exploring and proposing alternative methods for achieving the constant-parameter interface with better numerical properties and/or other modeling advantages is the third objective of this thesis.   13  Objective 4: Direct interfacing of rotor and stator machine terminals Direct interfacing of machine models with the external ac power networks has been the concern of most VBR models ?[28]??[31], ?[33] . However, if the rotor terminals are also required for interfacing, for example to study and design the thyristor rectifiers of the exciter systems ?[49], ?[50], the conventional VBR models do not provide a direct interface of the rotor circuit and therefore are not suitable. A recently published model ?[34] incorporates a direct rotor interface into the VBR formulation. However, this model has rotor-position-dependent (variable) inductances with a 4-by-4 full inductance matrix (3 stator windings and 1 field winding). Therefore, the final objective of this thesis is to develop a decoupled constant-parameter models that have capability of direct interfacing of stator and rotor circuits to arbitrary ac side grids and dc exciter systems, respectively.  The outcome of this research will be constant-parameter general purpose ac machine models especially for state-variable-based simulation programs. The equivalent interfacing circuit for either induction or synchronous machine will be composed of conventional circuit elements such as constant and decoupled RL branches, and controllable voltage-sources. The remaining differential and algebraic equations can also be readily implemented using the conventional blocks such as integrators, gains, summers, etc. which are available in most simulation programs.  14  The numerical advantages of the new VBR models with direct interface will offer more accurate and less computationally demanding solutions than the existing classical and built-in models. Availability of fast and accurate models will save many hours of time of uncountable number of power system engineers and researchers. More efficient models will also enables the users to model larger systems with more details and/or to run multiple simulations in a fraction of time.  The favourable properties of the new models will facilitate their adoption in various simulation programs. Software developers and even novice users should be able to implement such models in their programs and software applications. Moreover, having the same general structure and the interfacing circuit for both synchronous and induction machines will make it easier to develop electrical machines components with similar user interfaces and parameter entry methods.  15  CHAPTER 2:  CIRCUIT INTERFACING OF INDUCTION MACHINE MODELS Induction machines are usually used as motors with squirrel cage rotors. The wound-rotor machines are also popular especially in renewable energy systems. For the purpose of modeling, both types of induction machines are represented similarly. Having a symmetrical stator and rotor structures facilitates mathematical modeling of induction machines. To provide the reader with the background material helpful in understanding the proposed advanced models, the classical modeling of induction machines for power systems transient studies is briefly reviewed. The state-of-the-art VBR model is presented, and the new explicit general interfacing circuit including zero-sequence is proposed in this chapter. The numerical efficiency and accuracy of the proposed model is verified by extensive simulation studies included later in the chapter. 2.1     Induction Machine Modeling 2.1.1    Coupled-circuit phase-domain model  To make the reader familiar with the fundamentals of induction machine modeling, the classical magnetically-coupled circuit-model is shown in Figure ?2-1. This model has a 16  relatively small number of equation but is considered to offer sufficient accuracy for general purpose system-level transient simulations studies ?[1]. The VBR model is derived from the coupled-circuit model shown in Figure ?2-1.  Figure ?2-1 considers a basic three-phase wye-connected induction machine which is modeled as magnetically coupled circuits. The rotor may be squirrel cage or wound type; for the squirrel cage, the rotor terminals are short-circuited. Throughout this thesis, all electrical variables are assumed to be referred to the stator side using appropriate turns-ratio. The windings are also assumed to be sinusoidally distributed and the effects of slotting and magnetic saturation are neglected. The corresponding voltage equations in matrix form are ?[1] ? ? ?????????????????????abcrabcsrrTsrsrssabcrabcsppppiiLrLLLrvv  (?2?1)   Figure ?2-1  Magnetically coupled circuit model of induction machine. 17  where p is the Heaviside?s operator (differentiation with respect to time) and the resistance and inductance matrices are ? ?ssss rrrdiag?r  (?2?2) ? ?rrrr rrrdiag?r  (?2?3) ??????????????????????????mslsmsmsmsmslsmsmsmsmslssLLLLLLLLLLLL212121212121L (?2?4) ??????????????????????????mrlrmrmrmrmrlrmrmrmrmrlrrLLLLLLLLLLLL212121212121L (?2?5) ???????????????????????rrrrrrrrrsrsr L???????????????cos)32cos()32cos()32cos(cos)32cos()32cos()32cos(cosL. (?2?6) In the above equations, msL  and mrL  are the magnetizing inductance of the stator and rotor winding respectively; lsL  and lrL  are the leakage inductance of the stator and rotor winding respectively; and diag[X] is an nn?  matrix having vector X in its diagonal.  Finally, the electromagnetic torque equation in machine variable is given by 18  ? ? abcrsrrTabcsePT iLi ???? 2. (?2?7) The extended form of (2?7) is given in [1, see p. 147].  Throughout this thesis, a simple single rigid body mechanical system is assumed that is a simplified model of all the masses connected to the shaft including rotor, mechanical loads, and/or prime movers. If the external load torque is LT , the mechanical angular velocity is calculated by ? ?Ler TTJPp ?? 2? (?2?8) where J is the inertia of the rotor and connected load and r?  is the electrical angular velocity of the rotor. If the system is modeled in per-unit, (?2?8) becomes ?[1] ? ?Lebr TTHp ??????????21??  (?2?9) where b?  is the base angular velocity and H is the inertia constant (in seconds).  H is defined as the ratio of kinetic energy stored in the rotor at rated speed (in Joules) over the nominal power (in Watts) as BbBbPJPTJPH22221221 ?? ????????. (?2?10) 19  In (?2?10), TB and PB are the base torque and power respectively. The CC-PD model defined by (2?1) ? (2?7) can be implemented in simulation programs such as PLECS ?[17] and ASMG ?[18], that have ability to implement variable inductances. This model provides direct interface to external circuits on the stator and rotor sides. However, direct implementation of this model is computationally very expensive and the model has poorly scaled eigenvalues ?[32]. 2.1.2    The qd0 model in arbitrary reference-frame A change of variable is commonly used in analysis of ac machines to replace the physical variables (voltages, currents, and flux linkages) with transformed variables expressed in rotating reference frame. The matrix that transforms a set of three-phase variables to the arbitrary reference frame is expressed as ?[1] ?????????????????????212121)32sin()32sin(sin)32cos()32cos(cos32 ??????????sK. (?2?11) Transforming the stator and rotor variables of induction machine described in subsection ?2.1.1 to the arbitrary reference frame simplifies the voltage equations to the following ?[1] qsdsqssqs pirv ??? ??? (?2?12) dsqsdssds pirv ??? ??? (?2?13) 20  ssss pirv 000 ???  (?2?14) qrdrrqrrqr pirv ???? ???? )( (?2?15) drqrrdrrdr pirv ???? ???? )( (?2?16) qrrrr pirv 000 ??? (?2?17) where mqqslsqs iL ?? ?? (?2?18) mddslsds iL ?? ??  (?2?19) slss iL 00 ??  (?2?20) mqqrlrqr iL ?? ?? (?2?21) mddrlrdr iL ?? ??  (?2?22) rlrr iL 00 ??  (?2?23) The q- and d-axis magnetizing flux linkages are defined as  )( qrqsmmq iiL ??? (?2?24) )( drdsmmd iiL ???  (?2?25) where msm LL 23? (?2?26) 21  The circuit representation for the above equations is given in [1, see p. 151]. The electromagnetic torque can be calculated by ?[1] )(223dsqsqsdse iiPT ?? ?? (?2?27) which in per-unit is expressed as )( dsqsqsdsbe iiT ??? ??. (?2?28) This qd0-model defined by (2?12) ? (2?26) in state-variable form is a voltage-input and current-output system. Therefore, this model is typically interfaced with the external network represented in abc-phase coordinates using dependent current-sources.  2.2     Voltage-Behind-Reactance Formulations The voltage-behind-reactance formulation takes advantage of the reference frame transformation but leaves enough elements in abc phase coordinates to obtain direct interfacing to arbitrary network. Due to the symmetry of induction machines structure, the VBR formulation excluding magnetic saturation has constant interfacing circuit ?[32]. Inclusion of zero-sequence in ?[32] makes the formulation implicit. In this thesis, the interfacing circuit of ?[32] is changed by adding a zero-sequence branch to avoid forming an algebraic loop when zero-sequence current may exist. 22  Starting from the full-order induction machine model given by (?2?1) ? (2?7), after algebraic manipulation the following VBR form is derived ?[32] abcabcsabcsabcsabcsabcs p eiLirv ?????????  (?2?29) where the inductance matrix is ???????????SMMMSMMMSabcsLLLLLLLLL''L  (?2?30) and the entries are defined as mlsS LLL ???? 32  (?2?31) mM LL ???? 31 . (?2?32) The sub-transient magnetizing inductance is defined as 111 ????????? ????lrmm LLL. (?2?33) The resistance matrix is given by ???????????SMMMSMMMSabcsrrrrrrrrr''r  (?2?34) 23  and the entries are defined as rlrmsS rLLrr 22''32?? (?2?35) rlrmM rLLr 22''31??. (?2?36) The sub-transient back EMF voltages in (?2?29) in the stationary reference-frame are defined as  ? ? Tdqsabc ee ]0[1'' ?????? ?Ke  (?2?37) where ? ? qrlrmqrqlrrmdrq vLLLrLe''''2'''''' ???? ???? (?2?38) ? ? drlrmdrdlrrmqrd vLLLrLe''''2'''''' ????? ????. (?2?39) The sub-transient flux linkages of q- and d-axis are defined by lrqrmq LL?? ????? (?2?40) lrdrmd LL?? ?????. (?2?41) 24  The rotor subsystem state model is defined in terms of the rotor flux linkages as the state variables, as following qrdrrmqqrlrrqr vLrp ?????? ?????? )()( (?2?42) drqrrmddrlrrdr vLrp ?????? ?????? )()( (?2?43) where the magnetizing fluxes are  qqsmmq iL ?? ?????? (?2?44) ddsmmd iL ?? ?????? . (?2?45) There are a few possible equations to calculate the developed electromagnetic torque ?[1]. For this formulation, the torque can be calculated by using the magnetizing fluxes as )(223dsmqqsmde iiPT ?? ??. (?2?46) which in per-unit is expressed as )( dsmqqsmdbe iiT ??? ??. (?2?47) Note that in (?2?29), the inductance matrix (?2?30) and resistance matrix (?2?34) are independent of reference frame and constant, which is a consequence of machine symmetry. These are very desirable numerical properties that make the VBR model more efficient than the CC-PD model (see VBR-I in ?[32]). However, using (?2?29) directly to 25  interface the machine to the network requires implementation of mutual inductances and resistances which are not available in the majority of simulation programs.  A further simplification is achieved by diagonalizing the inductance and resistance matrices (?2?30) and (?2?34) using the stator current and zero-sequence current as ?[29] scsbsas iiii 03??? . (?2?48) Taking derivative of (?2?34) gives scsbsas pipipipi 03??? . (?2?49) The off-diagonal terms in (?2?30) and (?2?34) may be eliminated by expressing the stator currents associated with the off-diagonal entries in terms of the zero-sequence current and the remaining phase (diagonal) currents. After algebraic manipulations, (?2?29) can be written as ? ?ssabcabcsDabcsDabcs pLrpLr 0000'' 33 iieiiv ????? . (?2?50) where  Tssss iii ][ 0000 ?i  (?2?51) rlrmsD rLLrr 22???? (?2?52) mlsD LLL ????  (?2?53) 26  rlrmM rLLrr 220 31 ????? (?2?54) mM LLL ????? 310. (?2?55) If the stator winding neutral is floating, then the zero-sequence will not be present and (?2?50) will be simplified. However, for the general purpose interfacing circuit, the zero-sequence has to be included. The corresponding voltage equation is   ssss pirv 000 ???  (?2?56) where slss iL 00 ??  (?2?57) that gives the additional state equation )(1 000 ssslss irvLpi ??. (?2?58) Substitution of (2?58) into (?2?50) makes the following VBR formulation (See VBR-III in ?[32]) slsslssabcabcsDabcsDabcs LLLrLrpLr 00000''333 vieiiv ??????? ????? (?2?59) where 27  Tssss vvv ][ 0000 ?v  (?2?60) Implementation of (?2?59) requires only constant parameter decoupled RL-branches and basic built-in circuit elements, which is a major advantage. However, if the zero-sequence is included, the formulation (?2?59) becomes implicit with respect to the stator branch voltages abcsv  through the corresponding zero-sequence voltage sv0 , which are algebraically related. Therefore, a direct implementation of (?2?59) with the zero-sequence (as may be required in a general case) will result in an algebraic loop if the terminal voltage abcsv  is unknown, e.g. when the machine is in series with inductive elements. Presence of algebraic loops requires iterative solutions that extensively increase the computational burden. However, the algebraic loop can be avoided by adding a zero-sequence branch to the interfacing circuit. 2.2.1    Explicit formulation with zero-sequence branch in the interfacing circuit To avoid an implicit formulation, the stator voltage equation (?2?50) can be separated into two parts ngabcnabcs vvv ??. (?2?61) In the above equation, abcnv  is the phase-to-neutral voltages for the RL interfacing branches, which is 28  ''abcabcsDabcsDabcn pLr eiiv ??? . (?2?62) It is important to distinguish abcnv  from the stator phase voltages of the original induction machine, abcsv . Moreover, the voltage Tngngngng vvv ][?v  defines a new branch for representing the optional zero-sequence of the interfacing circuit, that is ngngssng piLiripLirv 000000 )3()3( ????. (?2?63)  Here, the variable ngi is the zero-sequence branch current. Equations (?2?61), (?2?62), and (?2?63) , altogether, define the four constant RL branches depicted in Figure ?2-2 that are needed for direct and explicit interface of machine models with external electrical circuits. The dependent voltage-sources abce ??  are functions of the branch current and determined by the rotor subsystem (?2?37), (?2?38), and (?2?39). If the zero-sequence is not needed, the branch ng  is simply removed leaving the neutral point n  floating.   Figure ?2-2  General interfacing circuit for induction machine VBR formulation. 29  2.3     Interfacing of Induction Machine Models in State-Variable-Based Simulation Programs  In the state-variable-based environment (e.g. MATLAB/Simulink ?[15], ?[16]), the qd0 model of induction machines (and synchronous machines as will be explained in section ?3.5) are interfaced to the external network (typically represented in abc phase coordinates) using three-phase voltage-controlled current-sources ?[8]. These models have voltage-input, current-output formulation as depicted in Figure ?2-3. Whenever it is possible to feed the machine model from either capacitors or the voltage-sources, for example, this interface becomes input-output compatible and explicit. However, the direct interfacing of machines models as current-source in series to inductive branches (which is in fact very common) represents a challenge. The state model generation algorithm is unable to formulate a proper state-space equations when the external circuit is inductive ?[8]. In this case, the interfacing is typically resolved using a resistive (e.g. see in Figure ?2-3) or capacitive snubber circuit. The snubber leads to numerical disadvantages, namely lose of accuracy and increased numerical stiffness.  The main goal of the VBR formulation is to achieve a direct interface of the machine model with the external (inductive) networks as depicted in Figure ?2-4. Moreover, to make the interface very effective and easy to use in most simulation packages, the interfacing circuit must be composed of very simple and constant-parameter branch elements, which is achieved for induction machine as presented in section ?2.2.  30   Figure ?2-3  Interfacing of induction machine qd0 models with external electrical network.   Figure ?2-4  Interfacing of induction machine VBR models with external electrical network. 2.3.1    Examples of implementation in SimPowerSystems, ASMG, and PLECS toolboxes The formulation proposed in subsection ?2.2.1 is easy to implement in almost any available simulation program. The implementation in the SimPowerSystems (SPS) ?[19] toolbox is 31  shown in Figure ?2-5, wherein the interfacing circuit is shown inside the box. For compactness, the machine is connected to the Th?venin equivalent circuit of a network; however, the network could also be represented in detail. In this figure, the machine is grounded through the resistance gr. If the machine is not grounded, the neutral point is left unconnected without any modification to the model. To focus on the electrical model, the mechanical system and the additional inputs and outputs are omitted in the figure. For the other two toolboxes, PLECS ?[17] and ASMG ?[18], the electrical network is simply replaced with one instance of ASMG-System or PLECS Circuit, respectively. The interfacing circuit   Figure ?2-5  Example of implementation of the proposed induction machine VBR model in Simulink using the SimPowerSystems toolbox. 32   Figure ?2-6  Example of implementation of the proposed induction machine VBR model equivalent interfacing circuit with external electrical network inside the PLECS Circuit block. implementation in PLECS ?[17] is shown in Figure ?2-6. It is also possible to implement the control blocks inside the PLECS Circuit instead of using Simulink blocks. The ASMG version 2 used in this thesis does not have a graphical user interface. It has been verified that all three implementations give identical results. 2.4     Computer Studies  To assess the numerical efficiency of the proposed explicit VBR model, a simple distributed generation system is considered here. It is assumed that a 4-pole 50-hp 60Hz induction machine is connected to a prime mover (e.g. a wind turbine) that in the course of the study maintains a constant speed of 1.027 pu. The system is illustrated in Figure ?2-7, wherein the 33  generator is connected to a Th?venin equivalent circuit that represents the rest of the ac grid. The machine parameters are taken from ?[51] and are also summarized in Appendix A. To emulate a severe unbalanced transient condition, a single-phase-to-ground fault is assumed in the system close to the generator feeder and the grounding resistance is set to zero (rg = 0).  Initially, the machine is in steady-state, and the fault occurs after one electrical cycle. The fault is emulated by decreasing voltage of phase a of the equivalent source to zero (va = 0). The VBR model is directly connected to the RL network similar to Figure ?2-5. However, the classical qd0 model requires a three-phase snubber for interfacing as shown in Figure ?2-7. To minimize the error, a 10 pu resistive snubber is used here. The snubber adds error and makes the system numerically stiff.   Figure ?2-7  Induction generator connected to the Th?venin equivalent circuit of a network for the single-phase-to-ground fault transient study. 34  As a reference, the system (including the machine) is represented in the synchronous reference frame and the whole model is solved using the ode45 solver with the maximum time-step limited to 1 ?s. A list of MATLAB ordinary differential equation (ODE) solvers is given in Appendix B. The reference model is implemented using basic Simulink blocks. Having a very simple and linear system, the conversion to the rotational reference frame would be straightforward in this example.   The implicit VBR-III ?[32] [(?2?59)] model is implemented in the SimPowerSystems (SPS) toolbox. The SimPowerSystems toolbox allows for solving algebraic loops in Simulink. The explicit VBR (subsection ?2.2.1) with the interface shown in Figure ?2-2 is implemented in the SimPowerSystems, ASMG, and PLECS toolboxes. The three implementations are shown to give identical results. Since the built-in qd models in the library of SimPowerSystems and PLECS do not include the zero-sequence, the qd0 machine model is implemented using basic Simulink blocks. This model is then interfaced to an instance of PLECS Circuit by means of controlled current-sources.  To show the consistency of the explicit VBR and the implicit VBR-III models, both models are run with the MATLAB solver ode45 using the same settings. The maximum and minimum time-steps are set to 1 ms and 0.1 ?s, respectively. The models are implemented in pu and the relative and absolute tolerances of the solver are set to 10?4. For the qd0 35  model with snubbers, the stiffly-stable MATLAB solver ode15s is used with settings that are identical to those used for the other VBR models.  The simulation results shown in Figure ?2-8 demonstrate the consistency between the VBR models and the reference solution. To show more details, the three fragment windows highlighted in Figure ?2-8 are enlarged and shown in Figure ?2-9 to Figure ?2-11, respectively. Figure ?2-9 shows phase c stator current in steady-state, where it can be seen that both VBR models yield exactly the same results. Moreover, the models have chosen identical time- steps. For the qd0 model, the snubbers sink part of the machine output current and therefore produce some error. Comparatively, the qd0 model with snubbers has chosen several times more time-steps and has a visible steady-state error. Figure ?2-10 and Figure ?2-11, which show phase c stator current and electromagnetic torque during the fault, confirm that the VBR models have no visible error in transient as well. For the qd0 model with snubbers, the transient response has some error, although it is less visible due to the comparatively large fault current. The stiffly-stable solver ode15s has chosen even smaller time-steps during the transient period (after the fault) than in steady-state. A quantitative evaluation of the considered transient study is summarized in Table ?2?1. The time-steps and calculation data are obtained from Simulink Profiler ?[16]. The relative error is calculated by comparing the predicted trajectory with the reference solution using the 2-norm error and normalizing the difference ?[52]. For example, the error for asi  trajectory (including all the solution points) is given by 36    Figure ?2-8  Simulation results for the single-phase-to-ground fault study. From top: source voltage, source current, machine neutral current, and electromechanical torque.  37   Figure ?2-9  Detailed view of source current ics in steady-state excerpt from Figure ?2-8.   Figure ?2-10  Detailed view of source current ics during transient excerpt from Figure ?2-8.  %100~~)(22 ???asasasas iiii? (?2?64) where asi~  is the reference solution trajectory. The average three-phase current error )( abcsi? , which is shown in Table ?2?1, is evaluated using the following 38   Figure ?2-11  Detailed view of electromagnetic torque Te during transient excerpt from Figure ?2-8. ? ?)()()(31)( csbsasabcs iii ???? ???i. (?2?65) Table ?2?1 verifies that both VBR formulations are algebraically identical to the reference. It also reveals the difference between the computational costs of the two VBR formulations by comparing the number of sub-transient voltage calculations. Practically, the implicit VBR-III is significantly slower since it requires iterations in each time-step for the algebraic loop solution (3865 calculations compared to 764 for the explicit VBR). The qd0 model is explicit but numerically stiff, and thus it has used several times more time-steps (989 compared to 110 times for the VBR models). The qd0 model has also more number of internal current calculations (7,070 times) than the implicit VBR sub-transient voltage calculations (3,865 times). As shown in Table ?2?1, the largest eigenvalue of the qd0 model with snubber is several orders of magnitude bigger than the largest eigenvalue of the VBR models. The eigenvalues of system are found by linearizing the model around operating point by using MATLAB/Simulink functions. 39  Based on the study presented in this section, as well as the ones in ?[32], the explicit VBR model is suggested as the general purpose model for squirrel-cage induction machines in state-variable-based simulation programs. This model yields identical results to the reference and does not require a snubber when connected to an inductive network/system, while offering high accuracy, numerical stability, and simulation efficiency.    Table ?2?1  Comparison of Numerical Efficiency of VBR Formulations for Single-Phase Fault Study Simulation Parameter Considered Formulations Implicit VBR-III Explicit VBR qd0 Number of Major Outputs/Steps calculations 110 110 989 Number of Internal Minor Calculations* 3,865 764 7,070 Current abcsi Prediction Error 0.000 % 0.000 % 2.861 % Largest eigenvalue ?199 ? j118 ?199 ? j118 ?1.18?105 * This row shows the number of the sub-transient voltage calculations for the VBR models and the number of injected current calculations for the qd0 model.  40  CHAPTER 3:  CIRCUIT INTERFACING OF SYNCHRONOUS MACHINE MODELS Synchronous machines typically have asymmetrical rotor structure and one field winding. The filed winding is essential for the control of the machine voltage and stability of power systems. The field winding and rotor saliency makes the synchronous machine models more complex than induction machine models. Specifically, the rotor asymmetry makes it more challenging to achieve the constant-parameter VBR model for the synchronous machines compared to symmetrical induction machines. In this Chapter, classical modeling of synchronous machines for power systems transient studies and the formulation of the models for efficient numerical implementation is provided. The classical voltage-behind-reactance (VBR) formulation is reviewed for completeness as well. A state-of-the-art VBR model is used as the basis for deriving the new synchronous machine VBR models with the same constant-parameter interface as the induction machine VBR models proposed in Chapter 2. Extensive computer studies presented in this chapter will demonstrate the accuracy and benefits of the proposed models.  41  3.1     Synchronous Machine Modeling 3.1.1    Coupled-circuit phase-domain model To familiarize the reader with basics of synchronous machine modeling for transient studies of power system, the classical model is reviewed. Similarly to the induction machines, synchronous machines are represented by lumped-parameter magnetically-coupled circuit-based models. Such general-purpose models have the advantage of simplicity while yielding a relatively small number of equations ?[1]. A synchronous machine may have a round- or salient-pole rotor and usually has one field winding. A basic three-phase wye-connected synchronous machine is considered here. It is assumed that the machine has M damper windings in the q-axis and N damper windings plus a field winding in the d-axis. The magnetically coupled-circuit representation of the synchronous machine model is shown in Figure ?3-1. For simplicity of equations, in this subsection only, it is assumed that the machine has two damper windings in q-axis, and one damper winding in d-axis. The standard assumption of sinusoidal distribution of winding is also used here. All variables are referred to the stator side. The voltage equations in terms of machine physical variable are ?[1] ? ? ???????????????????????qdrabcsrrTsrsrssqdrabcsppppiiLrLLLrvv32 (?3?1) where the voltage and current vectors are defined as 42   Figure ?3-1  Magnetically coupled circuit model of synchronous machine. ? ?Tcsbsasabcs fff?f  (?3?2) ? ?Tkdfdkqkqqdr ffff 121?f  (?3?3) where f can represent v or i. The resistance and inductance matrices are ? ?ssss rrrdiag?r  (?3?4) ? ?121diag kdfdkqkqr rrrr?r  (?3?5) 43  ???????????????????????????????????????????)32(2cos)(2cos21)3(2cos21)(2cos21)32(2cos)3(2cos21)3(2cos21)3(2cos212cos?????????????????rBAlsrBArBArBArBAlsrBArBArBArBAlssLLLLLLLLLLLLLLLLLLLLLL (?3?6) ?????????????????mdlkdmdmdmdlfdmqlkqmqmqmqlkqrLLLLLLLLLLLL12100000000L  (?3?7) ?????????????????????????)32sin()32sin()32cos()32cos()32sin()32sin()32cos()32cos(sinsincoscos????????????????????rmdrmdrmqrmqrmdrmdrmqrmqrmdrmdrmqrmqsrLLLLLLLLLLLLL. (?3?8) In the above equations, mqL and mdL  are the magnetizing inductance of the q- and d-axis respectively; lsL  is the leakage inductance of the stator, and ljL  is the leakage inductance of the rotor windings. Here, the index j  denotes the windings, 1kq , 2kq , fd , and 1kd . The electromagnetic torque equation in machine variables is given by ?????????????? qdrsrrTabcsabcslssrTabcse LPT iLiiILi ][)(][)(212 ??. (?3?9) The extended form of the above equation is given in [1. see p. 198]. The mechanical system can be modeled similar to induction machine discussed in Chapter 2 by using (?2?8) or (?2?9). 44  Equations (3?1) ? (3?9) define the classical CC-PD model. For direct interfacing to arbitrary networks, this model may be implemented in simulation programs that have ability to include variable inductances such as PLECS ?[17] and ASMG ?[18]. The CC-PD model is computationally expensive and is less numerically efficient than the VBR models ?[28] (also see section ?3.2). 3.1.2    The qd0 model in rotor reference-frame  The time-varying inductances from the voltage equations are eliminated if the variables are transformed to the rotor reference-frame by the well-known Park?s transformation ?[7]. The transformation matrix is given by (?2?11), where dtd rr /?? ? . Using Park?s transformation, the equations become more compact and the rotor-position-dependency is removed. Here, considering simplicity of equations unlike subsection ?3.1.1, the general case with M damper windings in q-axis and N damper windings in d-axis is considered. The change of variables applied to (?3?1) results in Park?s voltage equations qsdsrqssqs pirv ??? ??? (?3?10) dsqsrdssds pirv ??? ??? (?3?11) ssss pirv 000 ???  (?3?12) Mjpirv kqjkqjkqjkqj ?,1, ??? ? (?3?13) fdfdfdfd pirv ??? (?3?14) Njpirv kdjkdjkdjkdj ?,1, ??? ? (?3?15) 45  where mqqslsqs iL ?? ?? (?3?16) mddslsds iL ?? ??  (?3?17) slss iL 00 ??  (?3?18) MjiL mqkqjlkqjkqj ?,1, ??? ?? (?3?19) fdfdlfdfd iL ?? ?? (?3?20) NjiL mdkdjlkdjkdj ?,1, ??? ??  (?3?21) The magnetizing fluxes are )(1????Mjkqjqsmqmq iiL? (?3?22) )(1?????Njkdjfddsmdmd iiiL? (?3?23) The circuit representation for the above equations is given in [1, see p. 202] and is not included here. The electromagnetic torque can be calculated by )(223qsdsdsqse iiPT ?? ?? (?3?24) which in per-unit becomes )( qsdsdsqsbe iiT ??? ??. (?3?25) 46  The mechanical system is similar to the one Chapter 2 for induction machines, however since synchronous machines are usually used as generators, equation (?2?8) becomes ?[1] ? ?eMr TTJPp ?? 2? (?3?26) where TM is the mechanical input torque from the prime mover. Equation (?3?26) in per-unit becomes ? ?eMbr TTHp ??????????21??  (?3?27) Similar to the induction machine model, the classical qd0 synchronous machine model in state-variable form has voltage-input and current-output. This model is therefore interfaced to abc networks by the means of voltage-controlled current-sources.  3.2     Voltage-Behind-Reactance Formulations The synchronous machine voltage-behind-reactance model ?[28] partially solves the interfacing issue by changing the model structure to branches consisting of voltage-sources behind series RL elements. However, the asymmetry of the rotor of synchronous machine makes their VBR formulation more complex compared to the induction machines presented in Chapter 2. This asymmetry is not only due to the possible uneven paths for the magnetic flux in the q- and d-axis (salient-pole rotor vs. round-pole rotor), but also due 47  to the presence of the rotor field and damper windings, which together result in unequal parameters of the sub-transient equivalent circuits known as dynamic saliency.  3.2.1    Variable-parameter formulation The original algebraically-exact voltage-behind-reactance synchronous machine model was proposed in ?[28]. Unless there is no dynamic saliency, in general, this formulation has variable impedance interfacing circuit. The resistance matrix of the interfacing circuit can be made constant by transferring its time-variant components into the sub-transient voltage-source ?[29]. The final stator voltage equation is  '''' ])([ abcsabcsrabcsabcssabcs pr eiLiv ??? ?. (?3?28)  where p denotes Heaviside?s operator (differentiation with respect to time). The sub-transient inductance matrix is ?????????????????????????????????????????????????????????????????????????????????)32(2cos)(2cos21)3(2cos21)(2cos21)32(2cos)3(2cos21)3(2cos21)3(2cos212cos?????????????????rBAlsrBArBArBArBAlsrBArBArBArBAlsabcsLLLLLLLLLLLLLLLLLLLLLL  (?3?29)  where  3mqmdALLL ???????? (?3?30)  48  3mqmdBLLL ????????. (?3?31)  The sub-transient magnetizing inductances are determined by 1111?? ?????????????? ?Mj lkqjmqmq LLL  (?3?32)  11111?? ????????????? ?Nj lkdjlfdmdmd LLLL . (?3?33)  The sub-transient inductance of the q- and d-axis are mqlsq LLL ?????? (?3?34)  mdlsd LLL ?????? . (?3?35)  The inductance matrix (?3?29) is rotor-position dependent and time variant if 0???BL  or mdmq LL ????? (or dq LL ?????). The inequality of sub-transient inductances is referred to as dynamic saliency. In general, the sub-transient magnetizing inductances mqL ?? and mdL ??  are not equal even if the rotor is round due to influence of field and damper windings [see (?3?32) and (?3?33)].  The sub-transient sources are given by ? ? Tdqsabc ee ]0[1 ???????? ?Ke  (?3?36)  49  where sK  is Park?s transformation matrix from stationary abc to rotor qd0 reference frame given by (?2?11) ?[1]. The sub-transient voltages for this model are qsmqMj lkqjkqjMjkqjqlkqjkqjmqdrq iLLrLrLe 21212 )( ?????????????????????????????? ???????? (?3?37)  dsmdNj lkdjkdjlfdfdfddlfdfdmdfdlfdmdNjkdjdlkdjkdjmdqrd iLLrLrLrLvLLrLe 2122212 )()( ????????????????????????????????????? ??????????. (?3?38)  The sub-transient flux linkages are defined by ????????????? ??Mj lkqjkqjmqq LL1?? (?3?39)  ?????????????? ??Nj lkdjkdjlfdfdmdd LLL1???. (?3?40)  The rotor state equations are  MjLrp mqkqjlkqjkqjkqj ?,1,)( ???? ??? (?3?41)  fdmdfdlfdfdfd vLrp ???? )( ??? (?3?42)  NjLrp mdkdjlkdjkdjkdj ?,1,)( ???? ???. (?3?43)  50  The magnetizing fluxes are calculated by qqsmqmq iL ?? ??????  (?3?44)  ddsmdmd iL ?? ?????? . (?3?45)  Electromagnetic torque in SI units is calculated using the same equation as induction machine given by (?2?46). This VBR formulation has two parts. The first part defines the interfacing circuit given by (?3?28) and the variable coupled inductance matrix (?3?29). The second part includes the sub-transient voltages (?3?37) and (?3?38) that are the link between the rotor subsystem [with state equation given by (?3?41) ? (?3?43)] and the interfacing circuit. While the numerical properties of the classical VBR model are improved compared to the CC-PD model (subsection ?3.1.1), the variable mutual inductances of (?3?29) are still a concern that will be considered in this chapter and the following one.  3.2.2   Standard state-space form for the rotor subsystem For an easier implementation, the rotor subsystem state equations can be written in a standard state-space form: ?????????DUCeBUA??)('' rqdp? (?3?46)  51  where the state vector is ? ?kdNkdfdkqMkqT ????? ?? 11??  (?3?47)  and the inputs are ? ?fddsqsT vii?U . (?3?48)  The matrices A and B are defined as ???????????22)1()1(11A00AAMNNM   (?3?49)  ??????????221)1(211B00BBNM (?3?50)  where the sub-matrices are defined by ??????????????????????????????????????????????????????????????????????11121222221221121111111lkqMmqlkqMkqMlkqlkqMmqkqMlkqlkqMmqkqMlkqMlkqmqkqlkqmqlkqkqlkqlkqmqkqlkqMlkqmqkqlkqlkqmqkqlkqmqlkqkqLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLr???????A (?3?51)  52  ?????????????????????????????????????????????????????????????????????????????????????????????????111121222221222211211111112122lkdNmdlkdNkdNlkdlkdNmdkdNlkdlkdNmdkdNlfdlkdNmdkdNkdNlkdmdkdlkdmdlkdkdkdlkdmdkdlfdlkdmdkdkdNlkdmdkdkdlkdmdkdlkdmdlkdkdlfdlkdmdkdlkdNlfdmdfdlkdlfdmdfdlkdlfdmdfdlfdmdlfdfdLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLrLLLr?????????A (?3?52)  ???????? ???????lkqMmqkqMlkqmqkqlkqmqkqTLLrLLrLLr ?221111B (?3?53)  ?????????? ?????????0001221122??lkdNmdkdNlkdmdkdlkdmdkdlfdmdfdT LLrLLrLLrLLrB . (?3?54)  If G and H are defined as ?????Mj lkqjkqjmqLrLG122  (?3?55)  ????????Nj lkdjkdjmdlfdfdmdLrLLrLH12222 . (?3?56)  then )( rT ?C  and D  will become 53  ???????????????????????????????????????????????????????????????????????????????????????22222221111222222212111)(lkdNmdkdNlkdNlkdNmdrlkdmdkdlkdlkdmdrlkdmdkdlkdlkdmdrlfdmdfdlfdlfdmdrlkqMmqrlkqMmqkqMlkqMlkqmqrlkqmqkqlkqlkqmqrlkqmqkqlkqrTLLrLHLLLLrLHLLLLrLHLLLLrLHLLLLLLrLGLLLLrLGLLLLrLG????????????C (?3?57)  ???????????lfdmdLLHG000D (?3?58)  3.3     Approximation of Dynamic Saliency in VBR Models  The inductance matrix (?3?29) will be constant if BL ??  [(?3?31)] is equal to zero. This is achieved if the sub-transient inductances are equal, that is dq LL ????? (or mdmq LL ?????) which is not true in a general case. The sub-transient inductances are related to the operational impedances ?[30]. In particular, in very high frequencies range, the leakage inductances of the rotor windings (circuits) will dominate the resistances. Consequently, the value of the operational impedances in the q- and d-axis equivalent circuits (normalized with respect to 54  frequency) will become equal to the sub-transient inductances qL ?? and dL ?? , respectively. The traditional approximation based on averaging qL ?? and dL ??  was shown to reduce the model accuracy however other accurate approximation methods have been suggested in the literature ?[28], ?[30].  To better explain the proposed models, and for consistency of derivations, the methods for removing the dynamic saliency are briefly explained here.  3.3.1    Additional winding method A pioneering technique based on addition of an artificial damper winding was proposed in ?[30]. This method is based on adding one extra damper winding to the rotor circuit with the purpose of enforcing numerical equality of qL ?? and dL ?? . In doing so, the added winding should have sufficiently high resistance as not to impact the low-frequency operational impedance of the rotor circuit. But otherwise, the added winding does not have any physical meaning with respect to the original parameters of a given machine. Since the sub-transient inductance in d-axis is typically smaller, the additional winding is normally added to the q-axis equivalent circuit. Therefore, denoting this added winding as (M +1), its leakage inductance is calculated based on (?3?32) and (?3?33) as 1)1(11?? ?????????????? mqmdMlkq LLL. (?3?59)  55  This also adds one more state equation to (?3?41), where 1,1 ?? Mj ? . To distinguish the sub-transient quantities of the approximate model that uses the (M +1)th  winding, triple-prime sign (''' ) instead of double-prime sign ('' ) is used here. Thus, including this additional winding gives the desired result mdmdmq LLL ?????????? (?3?60)  which makes 0???BL  in (?3?29) and results in a constant inductance matrix ?????????????MSSSMSSSMabcsLLLLLLLLLL  (?3?61)  where  mdlsS LLL ???? 32  (?3?62)  mdM LL ???? 31 . (?3?63)  Thus, the voltage equation (?3?28) becomes ''abcsabcsSMMMSMMMSabcssabcs pLLLLLLLLLr eiiv ????????????? . (?3?64)  Here, the subscript M  in ML  should not be confused with the winding indexing.  56  For either round rotor or salient pole synchronous machines, the parameters are such that qd LL ?. However, for the sub-transient inductances (reactances) this relationship is typically reversed, i.e. dq LL ????? (or dq XX ?????). A table with typical parameters for the synchronous generators (turbo and hydro), condensers, and motors is found in [5, see p. 40, Table 2], wherein the condition dq XX ????? holds for the low, average, and high ranges of typical parameters. The same conclusion can be made for the synchronous machine parameters found in ?[6], as well as in ?[1] (where for the hydro generator, dq XX ?????). The condition dq XX ????? follows from the fact that the field winding leakage inductance is typically small, which makes the dX ??  smaller even when qd XX ? . Based on this observation, the additional artificial winding would be normally added to the q-axis equivalent circuit and its leakage inductance calculated using (?3?59) will have positive value (see Appendix C and Appendix D). However, if for some reason one has dq XX ?????, then the winding can be added to the d-axis instead in order to remove the numerical saliency and achieve a constant-parameter voltage equation as given by (?3?64). 3.3.2    Singular perturbation method  The artificial winding adds an additional state variable and a fast dynamic mode. The disadvantageous effect of this winding can be ?undone? at the expense of making this formulation implicit, which is achieved using the singular perturbation ?[31]. 57  It is observed that when the resistance of additional winding is approaching infinity, the largest eigenvalue of the system is also going to infinity. In this method, the aim is to remove the fast state from the model. To adequately separate the fast dynamic due to additional winding from the other dynamics of the model, a new state variable is defined as the voltage across the additional winding resistance as ?[31] )1()1()1( ??? ? MkqMkqMkq irv (?3?65)  which yields to the following state equation after required algebraic manipulation  ???? ????????????????????????????? ?????MjmqkqjlkqjkqjmdqssqdsdrqsdmdMkqlkqMmdMkq LrLireiLvLLvLLpv12)1()1( )()(1 ???? (?3?66)  where )1()1( / ??? MkqMkq rL?. The new sub-transient voltage of the q- and d-axis become qsMj lkqjkqjmdMkqMlkqmdMjkqjqlkqjkqjmddrq iLrLvLLLrLe ??????????????????????????????? ?????? 122)1()1(12 )( ???? (?3?67)  dsNj lkdjkdjlfdfdmdfddlfdfdmdfdlfdmdNjkdjdlkdjkdjmdqrd iLrLrLLrLvLLrLe ??????????????????????????????????????? ???? 1222212 )()( ??????  (?3?68)  where the q-axis sub-transient flux is defined by  ??????????????????????? ?? ??Mj MlkqMkqlkqjkqjmdq LLL1 )1()1(???. (?3?69)  58  If the winding resistance approaches infinity, ?  goes to zero which makes (?3?66) an algebraic equation 0)1()1( ??? ?? MkqMkq pvthenrif ?. (?3?70)  Finding )1( ?Mkqv from (?3?66) and substituting qe ??? from (?3?67) after simplification gives  ???? ?????????????????????????????? )(11)1(2)1()1( drdsdrqsdmdMlkqdmdMlkqmdMkq iLvLLLLLLv ???????????????????????????????????????????????????????? ????MjqssdmddmdlkqjmdkqjMjmqkqjdmdlkqjmdkqj irLLLLLLrLLLLr12212 1)(1 ?? (?3?71)  where the approximate magnetizing flux is defined as qqsmdmq iL ?? ??????????. (?3?72)  To remove the additional winding parameters from the voltage equations, it is considered that mqmqMkq thenrif ?? ??????? )1(. (?3?73)   Therefore, it can be assumed that mqmq ?? ????, which considering (?3?44) and (?3?72) gives qqsmdmqq iLL ?? ???????????? )(. (?3?74)  Combining (?3?39) with the above equation yields 59  ???????????????????? ??Mj lkqjkqjmqqsmdmqq LLiLL1)(?? (?3?75)  which is the flux linkage equation independent of additional winding parameters. In summary, this method uses the same number of state-variables as the variable impedance VBR model, but nonetheless it has a constant-parameter interface as (?3?64). The sub-transient flux linkage of q-axis are given by (?3?75) and the sub-transient voltage equations are defined by (?3?67) and (?3?71), and (?3?68). If the machine is connected to an infinite bus, then the bus voltage qsv is readily available, thus (?3?71) and consequently (?3?67) do not require the solution for the branch voltages and the overall formulation becomes explicit. However, in a general case, voltage qsv is calculated from the stator branch voltages abcsv  that is solved together with the external network which has the sub-transient voltages abce ??  as input. This implicit relationship between abcsv  and abce ??  forms an algebraic loop that has to be solved iteratively by the program solver. It should be pointed out that because the effect of additional (artificial) winding has been removed (up to infinite frequency), this formulation becomes algebraically equivalent to the original VBR model with rotor-position-dependant inductances. Further algebraic derivations and computer studies confirm this conclusion.  60  3.4    Generalized Constant-Parameter VBR Formulation The interfacing circuit introduced in Chapter 2 (see Figure ?2-2) consists of constant parameters decoupled RL branches and includes the zero-sequence. These are much desired properties for simple and direct interface of ac electric machine models to arbitrary networks. In this section, the constant-parameter VBR formulations (achieved by adding an artificial winding or using singular perturbation) are modified to present the same interfacing circuit as it was achieved for the induction machine. Here, again the off-diagonal elements of inductance matrix are removed by incorporation of the zero-sequence current given by (?2?48).  More specifically, considering (?3?64) and defining MSD LLL ??  (?3?76)  MLL ?0  (?3?77)  yields the following voltage equation  ''00 )3( abcsabcsDabcssabcs pLpLr eiiiv ????  (?3?78)  Consideration of zero-sequence for modeling of synchronous machines is more important than it was in case of induction machines. In many instances of synchronous generators, the neutral point of the stator winding is grounded through external impedance for protection and monitoring purposes. Therefore, it is important to develop a general purpose interfacing circuit with zero-sequence that will be fully capable of predicting 61  unbalanced faults and grounding currents for different grounding classes (low/high resistance/impedance, etc.) and protection requirements ?[53].  Similar to the induction machine formulation, a fourth branch for representing zero-sequence is expressed using (?2?61), where  abcabcsDabcssabcn pLr eiiv ?????  (?3?79)  The zero-sequence branch is defined by the voltage equation ngsng piLipLv 000 )3( ??. (?3?80)  Unlike the induction machine circuit, the zero-sequence branch here does not have a resistive part. Instead, the resistance is left in the sub-transient voltage-sources defined by (?3?37) and (?3?38). It can be pointed out that it is possible to move equal parts of resistance from the sub-transient voltage-sources to the voltage equation, which would add a resistive component to the zero-sequence branch and change the value of the interfacing circuit resistances. However, this step does not eliminate the direct feed-through of the stator currents (through sub-transient voltage equation) due to unequal resistances of the rotor q- and d-axis circuits. Such feed-through will require a current-dependent voltage-source which is allowed in most simulation programs.  Eventually, (?2?61), (?3?79), and (?3?80) define the four constant and decoupled RL-branches depicted in Figure ?2-2. This circuit is used to interface the synchronous machine model with the external network. The circuit parameters are 62  sD rr ?  (?3?81)  mdlsD LLL ????  (?3?82)  00 ?r  (?3?83)  mdM LLL ????? 310. (?3?84)  In summary, the sub-transient voltage-sources defined by (?3?37) and (?3?38) with the interfacing circuit given in Figure ?2-2 and state equations (?3?41), (?3?42), and (?3?43) define the constant-parameter synchronous machine VBR formulation.  To implement a constant parameter VBR using the additional winding technique, the standard state-space A, B, C, and D matrices given in subsection ?3.2.2 can be used. However, using the singular perturbation method requires modification of the output equations to incorporate the additional winding voltage )1( ?Mkqv into the sub-transient voltages. 3.4.1    Possible explicit and implicit formulations The approximation of dynamic saliency using additional winding described in subsection ?3.3.1 (and used before in this section to achieve a generalized interfacing circuit) achieves the explicit implementation. The additional winding inductance is given by (?3?59), but the user must choose the winding resistance. In general, choosing this resistance very large will improve the accuracy of this approximation at the expense of making the system 63  numerically stiffer. Otherwise, this method gives zero error in steady-state since the effect of damper windings in this condition diminishes. As mentioned in subsection ?3.3.2, by removing the fast state using singular perturbation method the algebraically exact model is achieved but an algebraic loop is added that makes the formulation implicit. Additional iterations will be required to solve the algebraic loop. Therefore, this approach can make the solution computationally very expensive.  3.5     Interfacing of Synchronous Machine Models in State-Variable-Based Simulation Programs Interfacing of synchronous machine qd0 models in state-variable-based simulation programs is similar to the induction machines case shown in Figure ?2-3 in Chapter 2. However, due to asymmetrical rotor structure and dynamic saliency, interfacing of synchronous machine VBR models is more complex as shown in Figure ?3-2 for the classical VBR formulation with variable parameters and Figure ?3-3 for the constant-parameter VBR formulation presented earlier in this chapter. Specifically, without approximation of dynamic saliency, as shown in Figure ?3-2, the synchronous machine VBR models will contain coupled and variable inductances (and even resistances ?[28]), which precludes forming the linear time invariant (LTI) state-space model with constant matrices A, B, C and D. A constant resistance matrix is obtained by algebraic feed-through of the stator currents to the sub-transient voltage sources ?[29] [see 64  (?3?37) and (?3?38)], i.e. ? ?rqdsqdrqd g ?,i,?e ??? . This algebraic feed-through, however, does not generally represent any issue for simulation programs since the current-controlled voltage-source can be easily connected with inductive branches. A basic inductive branch equation that is solved within the external circuit-system has the form ? ? edtLidriv ??? / . Therefore, an algebraic feed-through is being formed within the circuit solver from the input controlled voltage-sources e to the output branch voltages v. In some constant-parameter VBR formulations as illustrated in Figure ?3-3, to calculate the remaining machine state-variables, in addition to the branch (stator) currents i, the (stator) branch voltages v may also be required for the following reasons: ? To approximate the dynamic saliency ?[31] (also see Chapter 4, section ?4.2); ? To include the zero-sequence ?[32]; ? To include saturation ?[33]; etc. If the sub-transient voltages in machine subsystem are calculated in terms of the rotor state-variables and the stator currents, i.e. ? ?rqdsqdrqd g ?,i,?e ??? , then the algebraic loop is not formed as it is the case in ?[28], ?[33]. However, in these models the interfacing branches have variable impedances in a general case. Conversely, whenever (in addition to the rotor states) the branch voltages v are required for calculating the controlled voltage-sources e, i.e. ? ?rqdsqdsqdrqd g ?,v,i,?e ??? , an algebraic loop is being formed as depicted in Figure ?3-3.  65   Figure ?3-2  Interfacing synchronous machine classical VBR models with external electrical network-circuit with algebraic feed-through of current.   Figure ?3-3  Interfacing synchronous machine constant-parameter VBR models with external electrical network-circuit with algebraic feed-through of current and possible algebraic loop for voltage. The existence of such algebraic loop makes the overall system of equations of the external circuit and the machine subsystem implicit, which generally results in invoking an iterative 66  solution process and leads to increase of computing time ?[16], ?[54]. It is therefore very desirable firstly to have a constant-parameter decoupled interfacing circuit such that the LTI state-space model with constant matrices A, B, C and D can be formed; and secondly, to have an explicit formulation that does not require the branch voltages for the machine subsystem. 3.6     Implementation of Synchronous Machine VBR Formulations in Simulation Programs Constraints specific to individual toolboxes may affect the implementation of the entire system model and simulation accuracy, especially when other components such as switches and transformers are used. For example, in SimPowerSystems (SPS) ?[19], the ideal transformers are only allowed in simple cases when no series inductive elements exists. In PLECS, the ideal transformers are also available, but this toolbox cannot solve algebraic loops and prevents the use of constant time step ODE solvers if switches (such as diodes) are used in the system. There is no ideal switch model in SimPowerSystems, and typically some snubber circuits are required. The PLECS toolbox is the only toolbox that presently has built-in VBR machine models, but these models have variable inductances. The ASMG toolbox can implement the variable impedances for a machine of arbitrary configuration. The SimPowerSystems toolbox does not have variable impedance elements. 67  The rotor model for the explicit formulation could be implemented using the standard state-space A, B, C, D matrices given in subsection ?3.2.2, where an additional winding is added to the system. For the implicit formulation, A and B matrices stay the same as in the original VBR formulation (without the extra winding). However, the output equations are given by (?3?67), (?3?71), (?3?68), and (?3?75). Therefore, C and D matrices should be modified accordingly. 3.6.1    Examples of implementation in SimPowerSystems, ASMG, and PLECS toolboxes An example of implementing the explicit VBR model in SimPowerSystems V5.1 toolbox ?[19] is shown in Figure ?3-4, wherein the rotor subsystem is implemented using conventional Simulink blocks. To clearly show the electrical circuit interfacing, some of the outputs and the mechanical system are not shown. For the ASMG toolbox ?[18], the electrical circuit of Figure ?3-4 is replaced with one instance of the ASMG-System, wherein the circuit may be defined using text (or graphical representation in the newer version of this toolbox). The algebraic feed-through in the rotor subsystem causes false detection of algebraic loops in SimPowerSystems if any voltmeter block is used. For the ASMG toolbox, the model settings are changed to allow algebraic feed-through. A false detection of algebraic loops needlessly calls a routine that may significantly decrease simulation speed. 68   Figure ?3-4  Example of implementation of the synchronous machine CP-VBR model in SimPowerSystems toolbox.   Figure ?3-5  Example of implementation of the synchronous machine CP-VBR model in PLECS toolbox. 69  An example of implementing the CP-VBR model in the PLECS toolbox Version 3.3.1 ?[17] is shown in Figure ?3-5. To prevent false detection of algebraic loops, the rotor subsystem is included inside the PLECS Circuit and it is implemented using its basic components. 3.7     Computer Studies  A 125-kW synchronous machine is considered in the first three following subsections. The parameters of this machine are summarized in Appendix C. In the last subsection (?3.7.4), a large 555-MVA power systems generator is considered to show the accuracy of models for different machine sizes. The parameters of the larger machine are given in Appendix D. All considered models have been implemented using MATLAB/Simulink ?[16] and the PLECS ?[17], ASMG ?[18], or SimPowerSystems ?[19] toolboxes. The simulation studies are conducted on a personal computer with 2.53GHz Intel CPU. 3.7.1    Verification of the implicit VBR formulation The 125-kW machine is originally modeled with one damper winding in the q-axis and one damper winding in the d-axis. Wherever appropriate, one artificial damper winding 2kq  is added to the q-axis to approximate the dynamic saliency.  In the following study, a single machine infinite bus system is assumed in order to verify the consistency of the implicit voltage-behind-reactance (IVBR) model discussed in subsection ?3.3.2. The machine is assumed to initially operate as a generator in a steady-70  state defined by the field excitation voltage ? ? 15.2/ ?? fdfdmdxfd vrXe pu and mechanical torque 5.0?mT pu. A three phase fault is applied at time 05.0?t s to the machine terminals. The fault is subsequently removed after 3 electrical cycles (50 ms). The reference solution is produced by the conventional qd0 model implemented using standard Simulink blocks and run with very small time-step of s?1  to ensure high accuracy. The classical VBR model is implemented in the ASMG toolbox ?[18]. 3.7.1.1    Constant time-step simulation study To verify the consistency among the explicit/exact VBR model with variable inductance (VBR) (subsection ?3.2.1) and the implicit model using singular perturbation (IVBR) (subsection ?3.3.2), these models were run with large time-step of 1 ms using the MATLAB solver ode4. The corresponding transient responses predicted by various models are shown in Figure ?3-6. As can be seen in Figure ?3-6, all models produce responses that are visibly consistent and very close to each other. A magnified view of the current qsi shown in Figure ?3-7 (a) reveals that at such a large time-step there is a slight difference between the subject VBR models and the reference solution. However, as can also be seen in Figure ?3-7 (a), the exact variable-parameter explicit VBR and IVBR produce exactly the same solution point by point. This is expected result since these models are algebraically equivalent despite their different implementation as explained in subsection ?3.3.2.  71   Figure ?3-6  Transient response to the three phase fault predicted by the variable impedance and the constant-parameters implicit VBR models. 3.7.1.2    Variable time-step simulation study To compare the model efficiency, the same transient study was also run using variable integration time-step. To achieve reasonably accurate results, here the ode15s (see Appendix B) was used with relative and absolute error tolerances set to 10?4, and the maximum and minimum time-step set to 10?3 and 10?7 seconds, respectively. The result of this study looks very similar to Figure ?3-6. For comparison, the same magnified view of the current qsi is also shown in Figure ?3-7 (b), which shows that both VBR models take about72   Figure ?3-7  Detailed view of current iqs from the three phase fault study shown in Figure ?3-6; (a) constant time-step solution; and (b) variable time-step solution.  the same number of time-steps to satisfy the required tolerance and the solution is visibly very close to the reference. 3.7.2    Accuracy of the approximated constant-parameter VBR formulation To evaluate the accuracy of the approximate explicit model (AVBR) that has constant parameter interfacing circuit, and specifically the choice of resistance in the additional artificial winding 2kq , the same three phase fault study is considered here. For 73  conciseness, only the current qsi and its magnified view are shown in Figure ?3-8. Moreover, it was found that among other variables, the current qsi has the largest error due to the extra winding. The AVBR model is run with additional winding resistance 2kqr  set to 1 pu and 10 pu As is seen in Figure ?3-8 (a), when 2kqr is equal to 1 pu, there is a noticeable difference between the predicted transient in current qsi. This difference becomes much smaller when the resistance is increased to 10 pu.    Figure ?3-8  Current iqs from the three phase fault study predicted by the approximate constant-parameters VBR formulation: (a) overall transient; and (b) magnified view of the window in part (a). 74  To systematically evaluate the impact of winding resistance 2kqr on the solution accuracy, the resistance value is changed from 10?3 to 103 pu and the solution error is calculated. Here also the normalized 2-norm cumulative error defined by (?2?64) in section ?2.4 is considered. To ensure consistency with the reference solution and since the encountered error is primarily due to the saliency approximation, a time-step of s?1  was used here as well. The result of this evaluation is shown in Figure ?3-9 (a). For example, for  2kqr equal to 10 pu, the cumulative error in qsifor the chosen time span is 1.7%. As can be seen in Figure ?3-9 (a), increasing the additional winding resistance 2kqr improves the solution accuracy. This accuracy improvement comes at the price of increasing the magnitude of the largest eigenvalue of the system which results in numerical stiffness. In Figure ?3-9 (b), the largest eigenvalue calculated as a function of resistance 2kqr is also plotted. For comparison, the largest eigenvalue of the conventional (time invariant) qd0 model and the exact VBR model are shown in that figure as well. Similar to Chapter 2, section ?2.4, the eigenvalues of system are found by using MATLAB/Simulink functions by linearizing the models around the steady-state operating point. The eigenvalues of qd0 model are the smallest and the most favourable from numerical point of view if snubbers are not required. The exact VBR model is time varying and uses the stator branch currents as the independent variables, which results in higher magnitude of the largest eigenvalue. However, using 2kqr as large as 10 pu in the approximate VBR model raises the largest eigenvalue to the 104 range 75   Figure ?3-9  The effect of additional winding resistance on: (a) approximate model error in iqs; and (b) the system largest eigenvalue magnitude.  3.7.3    Small machine case-study To demonstrate the improvement achieved by the proposed constant-parameter interfacing circuit, the synchronous machine is connected to a source with 5% inductive impedance as depicted in Figure ?3-10. The considered source with inductive impedance76   Figure ?3-10  Synchronous generator connected to an inductive network. The snubbers are required only for the classical qd0 model.  represents the Th?venin equivalent of the external power system which is implemented as part of the power network using the transient simulation program and its circuit-interface (and therefore should not be combined with the machine?s stator leakages). To verify the models in a commonly-used commercially available tool, the electrical circuit has been implemented in the SimPowerSystems toolbox ?[19]. As explained in section ?3.5 and also in ?[8], in order to interface the conventional qd0 model with the external network (that is represented as a physical circuits in abc variables/coordinates), a snubber circuit has to be used. Due to lack of zero-sequence in the built-in qd0 model in SimPowerSystems, the classical qd0 model has been implemented using standard Simulink blocks and interfaced with controlled current-sources as shown in Figure ?2-3. Aiming for reasonably high accuracy (on the order of 1%), a 100 pu resistive snubber is used as shown in Figure ?3-10. In order to emulate a severe unbalanced condition, the machine?s neutral was assumed grounded and a single phase fault is applied 77  to phase a as depicted in Figure ?3-10. This study also verifies the applicability of the proposed interfacing circuit for modeling the unbalanced conditions with the zero-sequence, which may be needed particularly for modeling synchronous generators ?[53]. The reference solution has been obtained using the exact VBR model (with variable inductances, implemented in the ASMG toolbox ?[18]) that was run with time-step of s?1  to ensure very high accuracy. Calculating the cumulative error using (?2?64), it was found that the qd0 model with snubber results in about 1.2% error in the phase current csi . It was further determined that the same level of accuracy (cumulative error) can be achieved by the approximate VBR model if 2kqr is chosen equal to 15 pu. The transient response of the system resulting from applying and removing this fault is shown in Figure ?3-11. The magnified view of the phase current trajectory before and during the fault is shown in Figure ?3-12 and Figure ?3-13, respectively. The qd0 model has some steady-state error due to the snubbers, which can be readily seen in Figure ?3-12. However, the approximate VBR model has zero steady-state error since the effect of the damper windings disappears in steady-state. Moreover, since the transient phenomenon is predominantly low frequency, the accuracy of predicting the fault currents also remains very good as can be seen in Figure ?3-13. Also, as shown in Figure ?3-12 and Figure ?3-13, the IVBR model has even smaller error (about 0.1% cumulative) in transient and steady-state. To compare the numerical efficiency of considered models, the MATLAB solver ode15s with relative and absolute error tolerances set to 10?4, and the maximum and minimum78   Figure ?3-11  Transient response to a single phase to ground fault predicted by qd0 model with snubber versus approximated constant-parameters VBR and exact implicit VBR models.  time-step set to 10?3 and 10?7 seconds respectively is used here (as in sub-subsection ?3.7.1.2). The time-steps and CPU time taken by each model together with the magnitude of their respective largest eigenvalues are summarized in Table ?3?1 As can be seen in this table, the 100 pu snubber interface of the qd0 model results in significant numerical stiffness and the largest eigenvalue of 1.2?106. At the same time, the proposed interfacing circuit and the approximate VBR model with 2kqr equal to 15 pu results in the largest eigenvalue of 2.4?104, which is significantly smaller. This in turn explains very large number of steps and the CPU time required by the conventional qd0 model as compared to79   Figure ?3-12  Detailed view of current ibs in steady-state from Figure ?3-11.    Figure ?3-13  Detailed view of current ics during the single phase fault from Figure ?3-11. the proposed approximate VBR model. This is also consistent with Figure ?3-12 and Figure ?3-13, wherein it can be seen that the stiffer model requires smaller time-steps to satisfy the same error tolerances under the variable step integration. Table ?3?1 also verifies that the implicit VBR (IVBR) model is equivalent to the exact VBR model in terms of accuracy. The IVBR model gives accurate results and even requires fewer time-steps as shown in Table ?3?1. However, in this model, an algebraic loop is formed80  Table ?3?1  Simulation Efficiency for Single-Phase Fault Study for a Synchronous Machine Models Performance measure qd0 with 100 pu snubber  (section ?3.1.2) Approximate VBR with 15 pu rkq2 (section ?3.3.1) Implicit VBR  (section ?3.3.2) Number of Time-Steps 1871 557 400 Simulation Time 0.35 s 0.20 s 0.47 s Cumulative Error of ics 1.2% 1.2% 0.1% Largest Eigenvalue ?1.2?106 ?2.4?104 ?1.7?102   because the sub-transient voltages abce ??  (in addition to the rotor states) also depend on the branch voltages abcsv  (e.g. ?[31]). Therefore, the high accuracy of the IVBR model comes at the price of iterative solution and generally results in longer CPU time per step or longer overall simulation (computing) time as seen in the third row of Table ?3?1. At the same time, the AVBR model is explicit, non-iterative, and for the same local error tolerances requires less than half of the CPU time as compared to the IVBR model. Therefore, considering the number of integration time-steps by itself is not an accurate measure for comparing the efficiency of models. 3.7.4    Large machine case-study The studies presented so far were using a relatively small machine (125-kW, Appendix A). To verify the model for larger scale power systems machines, and also to compare 81  implementation results in different simulation programs side-by-side, a larger synchronous generator (555-MVA, Appendix C) connected to a network is implemented using three different toolboxes PLECS ?[17], ASMG ?[18], and SimPowerSystems (SPS) ?[19]. Figure ?3-14 shows the test system. Here, the network including a short transmission line and a unit transformer modeled as a Th?venin equivalent circuit with impedance ZS = (2+j16)%. The mechanical input torque is 0.9 pu, and the excitation voltage xfde is 2.35 pu. The generator is grounded through a 0.3 pu resistance. A single-phase-to-ground fault occurs in the source behind impedance ZS. To obtain a constant-parameter VBR (CP-VBR) model, a winding is added to the q-axis. Its inductance is calculated using (?3?59) and its resistance is chosen to be 2.0 pu to get approximately 1% cumulative error in the stator current. Since the SimPowerSystems and PLECS built-in qd0 models do not include the stator zero-sequence circuit, the qd0 model was implemented using standard Simulink blocks, and then that model was interfaced  Figure ?3-14  A large synchronous generator connected to an inductive network. The snubbers are required only for the classical qd0 model. 82  using controlled current-sources in the SimPowerSystems toolbox. A 15 pu three-phase resistive snubber was used therein. This snubber gives a similar stator current cumulative error (around 1%) as the CP-VBR model mentioned above. As a reference, the entire system implemented in conventional Simulink blocks (without any snubbers) is solved using the ode45 solver (see Appendix B) with a very small time-step of s0.1 ???t . The qd0 model with snubber and CP-VBR models (implemented in different toolboxes) are numerically stiff. Thus, the ode15s solver is selected for the studies. The minimum and maximum time-step settings are 1 ms and 0.1 ?s, respectively. The initial time-step is set to the minimum time-step. All models are represented in pu and the relative and absolute tolerances of the solver are set to 10?4. The simulation results of the considered fault study are shown in Figure ?3-15, where it can be seen that all models produce consistent transient trajectories. To show more details, three fragments of Figure ?3-15 are enlarged and shown in Figure ?3-16 to Figure ?3-18, where it can be seen that the CP-VBR models use much fewer time-steps than the qd0 model with snubber. Figure ?3-16  shows that unlike the qd0 model, the CP-VBR models have no error in steady-state.  A more quantitative evaluation of this simulation study is presented in Table ?3?2. As defined by (?2?64) in section ?2.4, the 2-norm relative error is chosen as the mean to evaluate the accuracy of the variables. For the stator three-phase current iabcs, the average error of the three currents is given in Table ?3?2. This table shows that both formulations have similar stator current cumulative errors. However, the CP-VBR models are faster and require less than one third the numbers of time-steps of the qd0 model. The largest 83  eigenvalue of the qd0 model is about 46 times bigger than it is in the CP-VBR models, which shows that this model is significantly stiffer. There are also considerable differences between the run times of the different toolboxes, highlighting the numerical efficiency of their respective algorithms. The error in electromagnetic torque eT  is larger for the qd0 model due to the snubber that adds error to both q- and d-axis.   Figure ?3-15  Simulation results for large machine single-phase-to-ground fault study as predicted by CP-VBR models and the conventional qd0 model.  84   Figure ?3-16  Detailed view of phase current in steady-state shown in Figure ?3-15.   Figure ?3-17  Detailed view of phase current during transient shown in Figure ?3-15.   Figure ?3-18  Detailed view of electromagnetic torque shown in Figure ?3-15.    85  Table ?3?2  Numerical Efficiency of CP-VBR Models Versus qd0 with Snubber for the Large Machine Study Implementation Performance measures Time-steps Error iabcs Error Te Run Time Max. |?| * qd0 with snubber SPS 797 1.074 % 1.332 % 250 ms 74,111 CP-VBR SPS 251 1.069 % 0.559 % 188 ms 1,594 CP-VBR ASMG 255 1.054 % 0.559 % 73 ms 1,594 CP-VBR PLECS 262 1.052 % 0.558 % 180 ms 1,594 * ? = Eigenvalue        86  CHAPTER 4:  NUMERICAL METHODS TO ACHIEVE CONSTANT-PARAMETER VBR FORMULATIONS In this chapter, alternative formulations to achieve the constant-parameter VBR interfacing circuit of Figure ?2-2 are proposed and verified by simulations in commercial simulation programs. Although these formulations achieve the same goal as the added damper winding described in Chapter 3, section ?3.3, their derivation and application appear very intuitive and easy-to-use, which can offer additional advantages. The method presented in this chapter is based on mathematical equalization of q- and d-axis fictitious circuits putting the unequal parts to the sub-transient voltage-sources. This method is also used in the next chapter to extend direct interfacing to the rotor circuit. The proposed method is validated by extensive computer studies presented at the end of this chapter. The same generic mechanical system as the previous chapters is considered here.  4.1     Method of Using Current Derivatives To obtain an equivalent formulation, considering (?3?16) to (?3?18), (?3?44), (?3?45), (?3?34), and (?3?35), the stator voltage equations given by (?3?10) to (?3?12) are rewritten as ?[28] 87  )()( qqsqddsdrqssqs iLpiLirv ??? ????????????? (?4?1) )()( ddsdqqsqrdssds iLpiLirv ??? ????????????? (?4?2) )( 000 slssss iLpirv ?? . (?4?3) To match the interfacing circuit of Figure ?2-2, (?4?1) and (?4?2) are rewritten as qqsddsdrqssqs eiLpiLirv ?????????? ? (?4?4) ddsdqsdrdssds eiLpiLirv ?????????? ? (?4?5) where  qqsdqdrq piLLpe ??? ?????????????? )( (?4?6) dqsdqrqrd piLLe ???? ??????????????? )(. (?4?7) Following the classical VBR model derivation procedure ?[28], and considering (?3?39) and (?3?40), the terms qp? ??  and dp? ??  in (?4?6) and (?4?7) are replaced by the rotor state equations (?3?41) to (?3?43) to get ?? ???????????????????????Mjkqjmqlkqjkqjmqqsdqdrq LrLiLLpe12 )()( ???? (?4?8) ?????????????????????????? ??fdlfdmdNjkdjmdlkdjkdjmdqsdqrqrd vLLLrLiLLe12 )()( ????? 88  )(2 fdmdlfdfdmdLrL ?? ???. (?4?9) Finally, transforming (?4?4), (?4?5), and (?4?3) into stationary abc-phase coordinates gives the branch voltage equation (?3?78). The interfacing circuit of the proposed model is thus identical to Figure ?2-2, and its parameters are defined by (?3?81) to (?3?84). However, a direct and straightforward implementation of this model in a state-variable-based simulation package requires the current derivative term qsdq iLLp )( ????? in order to calculate qe ?? using (?4?6). This, in turn, leads to a non-proper state model. Approximation techniques will be presented in subsection ?4.4.1 to resolve this issue. The electromagnetic torque equation is not changed and the torque is calculated using the same equations given by (?3?24) or (?3?25) in Chapter 3. Mechanical system also stays unchanged. 4.2    Method of Using Algebraic Feed-through If the machine?s terminal voltages are known (i.e. as external inputs or state variables), the stator current derivative qspi in (?4?6) is calculated using state variables and inputs from the circuit subsystem. This results in an algebraically exact VBR model (without any approximations). In order to express qspi in terms of inputs and states, (?4?4) and (?4?8) are first combined to obtain 89  ???????????????????????????? ??Mjkqjmqlkqjkqjmqqdrdsdrqssqsqqs LrLLiLirvLpi12 )(/)(1 ?????. (?4?10) Substituting (?4?10) into (?4?8) again yields ? ? )()()(12 qssqsqdqMjkqjmqlkqjkqjmqdsdqdrqdq irvLLLLrLiLLLe ???????????????????????????????????????????? ??????. (?4?11) Finally, after replacing (?4?8) by (?4?11) and transforming the sub-transient voltages to abc-phase coordinates (similar to section ?4.1), the second VBR formulation is completed. In a general case, the terminal voltage abcsv  is unknown (e.g. machine connected to an inductive network). In such cases, this VBR formulation (including the circuit model) creates an implicit set of equations that results in an algebraic loop which contains input voltage qsv and sub-transient voltages abce ?? : in this loop, branch voltage abcsv , which is an output of the circuit subsystem, is a function of the rotor subsystem output abce ?? , while abce ??  itself is a function of abcsv  (both having algebraic feed-through). It is further shown in the next section (section ?4.3) that this formulation is algebraically equivalent to the model derived using singular perturbation approach ?[31]. The techniques to break this algebraic loop will be described in subsection ?4.4.2. 90  4.3     Algebraic Equivalence of Implicit VBR Formulations  An interesting (but not obvious) fact that has contributed to the idea of section ?4.1 of this chapter is that the constant-parameter singular perturbation VBR model ?[31] and the implicit VBR model from subsection ?3.3.2 (herein referred to as IVBR1) are algebraically equivalent to the implicit constant-parameter model derived in section ?4.2 (herein referred to as IVBR2). The interfacing circuits and the remaining state equations of IVBR1 and IVBR2 are identical; however, their sub-transient voltages are formulated differently. Therefore, it remains to show that the q- and d-axis sub-transient voltages of IVBR1 and IVBR2 are in fact equivalent.  For consistency with Chapter 3, the variables altered/approximated by the fictitious )1( ?M th q-axis damper winding are denoted by the triple prime symbol (''' ). The q-axis sub-transient flux of IVBR1 is given by (?3?69). The q-axis magnetizing flux is given by (?3?72), and the sub-transient back EMF voltages are given by (?3?67) and (?3?68). Since )1()1( / ??? MkqMkq rL?,  0)1( ???? ?thenrif Mkq (?4?12) Therefore, the voltage induced in the additional winding given by (?3?66) can be expressed as 91  ????????????????????????????????? ???? ?????MjmqkqjlkqjkqjmdqsmdMlkqmdMkq LrLpiLLLv121)1()1( )(1 ??. (?4?13) Based on (?3?59), the following equalities are derived mqmdMlkqmdLLLL????????????? ???? )1(1 (?4?14) )1( ???????????Mlkqmdmqmdmq LLLLL. (?4?15) By substituting (?3?44) and (?3?74) into (?3?67), after some algebraic manipulations, qe ??? is rewritten as ???? ?????????????????????MjkqjmqlkqjkqjmdMlkqmdMkqdrq LrLLLve12)1()1( )( ????. (?4?16) Substituting (?4?14) into (?4?13), and considering (?3?73), expression for )1( ?Mkqv simplifies to ??????????????????????????? ??? qsmdMjmqkqjlkqjkqjmdmdmqMkq piLLrLLLv12)1( )( ??. (?4?17) Substituting (4?17) into (4?16) and rearranging the resulting terms yields  ???????????????????????? ???MjkqjmqlkqjkqjmqqsMlkqmdmqdrq LrLpiLLLe12)1()( ???? 92  ?? ????????????????????????????MjmdmqMlkqmdmqkqjmqlkqjkqj LLLLLLr1 )1(2 )()( ??. (?4?18) Substituting (?4?15) into (4?18) further yields ??????????????????????????Mjkqjmqlkqjkqjmqqsmdmqdrq LrLpiLLe12 )()( ????. (?4?19) Based on (?3?32) and (?3?33), we have dqmdmq LLLL ??????????? (?4?20) This proves that (?4?19) is identical to (?4?8), i.e., the q-axis sub-transient voltages of IVBR1 and IVBR2 are equal. Next, it is required to show that the d-axis sub-transient voltages of IVBR1 and IVBR2 are equal. Therefore, (?3?74) is substituted into (?3?68) to get ????????????????????????????? ??fdlfdmdNjkdjdlkdjkdjmdqsmdmqrqrd vLLLrLiLLe12 )()( ????? dsmdNj lkdjkdjlfdfdfddlfdfdmd iLLrLrLrL 21222 )( ????????????????? ????. (?4?21) After substituting (?3?45) into (?4?22) and some algebraic manipulations, de ???  simplifies to 93  ?????????????????????????? ??Njkdjmdlkdjkdjmdqsmdmqrqrd LrLiLLe12 )()( ????? )(2 fdmdlfdfdmdfdlfdmdLrLvLL ?? ??????. (?4?22) Considering (?4?20), this result proves that (?4?22) is equal to (?4?9). Comparing (?4?19) and (?4?22) for IVBR1 with (?4?8) and (?4?9) for IVBR2, respectively, it is concluded that these two models are algebraically equivalent. 4.4     Numerically Efficient Explicit Implementation 4.4.1    Approximation of current derivative Exact numerical calculation of current derivative qspi in (?4?8) requires infinitely small time-steps. However, qspi may be approximated with reasonable accuracy using filters or other numerical techniques. Specifically, in continuous-time domain, the derivative of a signal may be approximated using a first-order high-pass filter ?[16] such as 00)( psspsH i ??. (?4?23) 94  which has a pole at 0p? . This filter adds a pole to the overall system. If 0p  is sufficiently large, the filter approximation will have high accuracy. However, the added pole, with magnitude around 0p , could make the overall system numerically stiff.  Alternatively, considering a fixed time-step t? , the derivative may be estimated by using the previous and current values of the stator current (i.e. backward difference), yielding the following discrete-time filter tzzH i ????11)(. (?4?24) where 1?z  denotes a unit delay. For variable time-steps, qspi may be approximated using the first-order backward differentiation formula 11)()()(?????nnnqsnqsnqs tttititpi. (?4?25) Here, the subscripts n  and 1?n  denote the current and previous time-steps, respectively.  The approximations (?4?24) and (?4?25) may add very little computational cost to the VBR model. Moreover, unlike (?4?23), they do not add an extra continuous state variable to the system. The implementation of the proposed constant-parameter VBR model interfaced with an external circuit-system in a typical state-variable-based transient simulation program (e.g. Simulink) is depicted in Figure ?4-1. Here, the interfacing circuit of Figure ?2-2 may be implemented in a circuit program or toolbox (PLECS ?[17], ASMG ?[18],95   Figure ?4-1  Implementation of explicit constant-parameter VBR model using filter Hi in continuous-time (high-pass) or discrete-time (backward difference) to approximate the current derivative.   SimPowerSystems ?[19], etc.) by using basic circuit elements (i.e. resistors, inductors, and controlled voltage-sources). The machine?s rotor state model, the reference frame transformations, and the high-pass filter iH  are implemented using conventional library components (gains, summations, integrators, functions, etc.). 4.4.2    Relaxation of algebraic loop An alternative implementation to current derivative estimation was proposed in section ?4.2, wherein the sub-transient voltage abce ??  evaluated in the rotor subsystem is function of the machine?s terminal voltage qsv that forms an algebraic loop. To break this algebraic96   Figure ?4-2  Implementation of explicit constant-parameter VBR model using filter Hv in continuous time (low-pass) or discrete time (zero- or first-order-hold) to break the algebraic loop.  loop, the voltage qsv can be approximated using a filter vH . A block diagram depicting this implementation is shown in Figure ?4-2. In continuous-time domain, vH  can be implemented using a first-order low-pass filter such as 00)( pspsHv ??, (?4?26) where 0p?  is the filter pole. Similarly to the method of subsection ?4.4.1, a large value of 0p  will make the approximation more accurate at the expense of making the overall system numerically stiffer.  97  A zero-order-hold can also be used to break this algebraic loop. In discrete-time domain, this results in the following filter 1)( ?? zzHv . (?4?27) It is also possible to approximate qsv with higher order filters, for example using a first-order-hold filter 212)( ?? ?? zzzHv . (?4?28) where a constant time-step is assumed. The zero-order-hold simply introduces a one-time-step delay, whereas the first-order-hold implies linear prediction based on the two previous values. An example of using the zero- and first-order-hold filters to approximate a sample rising signal (qsv) is depicted in Figure ?4-3.  Figure ?4-3  Approximation of voltage vqs using zero-order-hold (delay) or first-order-hold (linear prediction). 98  If the time-step is not constant, the discrete-time approximations of qsv by the zero- and first-order-hold filters respectively become )()( 1?? nqsnqs tvtv. (?4?29) )()(1)( 22111211????????????????????? nqsnnnnnqsnnnnnqs tvtttttvtttttv. (?4?30) Using (?4?27) to (?4?30), as opposed to (?4?26), the system numerical stiffness is not changed and no extra continuous state variables are added, which minimizes the additional computational cost. 4.5    Summary of Approximation Techniques Based on the material presented in the last two subsections (?4.4.1 and ?4.4.2), the proposed approximation techniques can be classified into two groups i. Continuous-time approximations ii. Discrete-time approximations A summary of this classification is presented in Table ?4?1. Type (i) approximations include the use of continuous-time filters to obtain either the current derivative, (?4?23), or the voltage algebraic loop relaxation, (?4?26). It can also include the addition of the fictitious damper winding ?[30] (see Chapter 3, subsection ?3.3.1). All these techniques increase the 99  number of continuous state variables by one, and could also make the overall system numerically stiff if the filter pole (or the added winding resistance) is made very large in an attempt to achieve high accuracy. Type (ii) approximation techniques, i.e. (?4?24), (?4?25), and (?4?27) ? (?4?30), are based on discrete-time filters. Such methods therefore do not add continuous state variables and barely increase the computational cost. No discrete-time form is available for the additional damper winding method. Table ?4?1  Summary of Approximation Techniques to Achieve the Interfacing Circuit of Figure ?2-2 Approximation Technique Current Derivative Approximation Voltage Algebraic Loop Relaxation Circuit-Level Approximation Continuous-Time High-Pass Filter Low-Pass Filter Fictitious Damper Winding Discrete-Time Backward Difference Zero-/First- Order-Hold N/A  4.6     Computer Studies In order to validate and assess the accuracy of the proposed explicit formulations (Table ?4?1), the system shown in Figure ?4-4 is considered in this section. This network contains a 555-MVA steam turbine synchronous generator whose neutral point is grounded through a reactor gX for protection purposes ?[53]. The generator is connected to a large network 100  (represented by an infinite bus) through a Wye/Delta unit transformer. The machine and network parameters are summarized in Appendix E. The generator circuit model ?[4] includes two damper windings in the q-axis, and one damper winding and a field winding in the d-axis. The magnetizing reactances mqX and mdX  are very close (1.61 pu and 1.66 pu). However, due to the dominant effect of the damper and field windings in the sub-transient period, the sub-transient reactances mqX ?? and mdX ??  are not equal (0.10 pu vs. 0.08 pu). Approximation techniques are therefore required to obtain the constant-parameter interfacing circuit of Figure ?2-2. In the considered transient study, the generator is assumed to operate initially in steady-state condition with 9.0?eT pu and 35.2?xfde pu, while the voltage of the infinite bus is set to its nominal value of 1 pu. At 0175.0?t s, a single phase-to-ground fault is applied at the machine?s terminals. This study was chosen as to emulate a highly unbalanced condition when qsi and qsv vary greatly to test the proposed approximation techniques in an unfavourable condition.  All considered models and approaches have been implemented using MATLAB/Simulink ?[15], ?[16] and the PLECS toolbox ?[17]. The general-purpose MATLAB solver ode45 (see Appendix B) is used to run all simulation studies with relative and absolute error tolerances both set to 10?3, and the maximum and minimum time-steps set to 10?3 and 10?7 seconds, respectively. Since all models are consistent, the results of all discrete-time approximation methods converge if the time-step is small enough. For the purpose of 101  comparison, the reference solution is obtained by simulating the study with a very small step size of 710???t seconds. All simulation studies are conducted on a personal computer (PC) with a 2.53GHz Intel CPU. To evaluate the accuracy of machine models, the predicted stator currents are compared to a reference. Here, the 2-norm (cumulative) relative error [defined by (?2?64) in section ?2.4] of the stator currents asi , bsi , and csi  for each approximation technique is computed.  As the error may be different in each phase, an average error is considered in this chapter as a more consistent measure of accuracy.   Figure ?4-4  Test system consisting of a grounded steam turbine generator connected to the network via a unit transformer. A single phase-to-ground fault is applied at the machine terminals. 102  4.6.1    Continuous-time approximation techniques The selection of the pole (or cut-off frequency) of the continuous-time filters (?4?23) and (?4?26) (and similarly, the resistance of the additional winding approach) involves a trade-off between accuracy and numerical stiffness (see subsections ?4.4.1 and ?4.4.2). Therefore, to objectively compare the considered approaches, the filter pole and the additional winding resistance values are chosen as to yield the same largest eigenvalue magnitude ( 1000~ ) during the fault condition. The corresponding pole and winding resistance values are summarized in the first column of Table ?4?2. The system eigenvalues will be shown in subsection ?4.6.3. The single phase-to-ground fault study is first executed using the current filter (?4?23), the voltage filter (?4?26), and the additional winding approach. The predicted trajectories of abcsv , asi , ngi , fdi , and eT  are shown in Figure ?4-5. Therein, it is observed that all three approximation techniques give visibly accurate solutions. To see the difference among the three approaches, portions of Figure ?4-5 are enlarged and reproduced in Figure ?4-6 to Figure ?4-8. The magnified views in Figure ?4-6 to Figure ?4-8 demonstrate that in this case all techniques yield acceptable solutions. However, the voltage filter and the additional winding approaches are noticeably more accurate than the current filter technique. The corresponding cumulative errors )( abcsavg i? are summarized in Table ?4?2, where similar observations are made. In particular, the error of the voltage filter (0.68%) is103  Table ?4?2  Comparison of Continuous-Time Approximation Methods (with Added Eigenvalue of 1,000)  Number of Time-Steps Error  ?avg(iabcs)  Simulation Time Additional Winding (rkq3 = 1.247 pu) 126 0.79 % 175 ms Current Filter (p0 = 947 ) 126 2.96 % 186 ms Voltage Filter (p0 = 1031) 126 0.68 % 179 ms   Figure ?4-5  Transient response to a single phase-to-ground fault predicted by continuous-time approximated models. From top to bottom: terminal voltages vabcs, phase a stator current ias, neutral grounding current ing, field current ifd, and electromagnetic torque Te. 104   Figure ?4-6  Detailed view of current ias shown in Figure ?4-5 for continuous-time approximated models.   Figure ?4-7  Detailed view of current ifd shown in Figure ?4-5 for continuous-time approximated models. slightly smaller than that of the additional winding approach (0.79%). Table ?4?2 also contains the number of time-steps and the overall CPU processing time taken by each model which are very similar for all considered models. To further investigate the effect of the added pole on the model accuracy, the cumulative error  )( abcsavg i? has been calculated for different largest eigenvalue magnitudes (from 1,000 to 10,000). The corresponding results are shown in Figure ?4-9 which demonstrate105   Figure ?4-8  Detailed view of electromagnetic torque Te shown in Figure ?4-5 for continuous-time approximated models.   Figure ?4-9  Cumulative error in predicting currents iabcs versus the magnitude of the added eigenvalue for continuous-time approximated models. that all continuous-time approximation techniques are consistent, as their cumulative error decreases when increasing their largest eigenvalue magnitude. Overall, the voltage filter has the best accuracy, closely followed by the additional winding method. 106  4.6.2    Discrete-time approximation techniques Here, the same single phase-to-ground fault study as subsection ?4.6.1 is repeated using the discrete-time approximations: first-order numerical derivative (?4?24); the zero-order-hold (?4?27); and the first-order-hold (?4?28), which are implemented using memory blocks. All models produced results that are visibly very similar to Figure ?4-6. The magnified views of the predicted trajectories of asi , fdi , and eT  are shown in Figure ?4-10 to Figure ?4-12, respectively. Similarly to the continuous-time approaches, all three discrete-time approximations yield acceptable and visibly close results, with some preference given to the first-order-hold of qsv method.  A quantitative assessment of this study is summarized in Table ?4?3 that shows the smallest error, 78.0)( ?abcsavg i?%, is in the same range as those obtained using the most accurate continuous-time techniques (Table ?4?2). The numerical efficiency and the overall computational times of all approaches are almost identical.  Figure ?4-10  Detailed view of current ias shown in Figure ?4-5 for discrete-time approximated models. 107   Figure ?4-11  Detailed view of current ifd shown in Figure ?4-5 for discrete-time approximated models.   Figure ?4-12  Detailed view of electromagnetic torque Te shown in Figure ?4-5 for discrete-time approximated models.  Table ?4?3  Comparison of Discrete-Time Approximation Methods  No. of Time-Steps Error ?avg(iabcs) Simulation Time Backward Difference of  iqs 126 2.64 % 177 ms Zero-Order-Hold of  vqs 126 1.71 % 170 ms First-Order-Hold of  vqs 126 0.78 % 172 ms 108   Figure ?4-13  Cumulative error in predicting currents iabcs versus maximum time-step size for discrete-time approximated models. Next, the study is repeated using different values of maximum time-step ranging from 10 ?s to 2 ms. The resulting cumulative error )( abcsavg i? as a function of time-step is plotted in Figure ?4-13. This figure shows that all three discrete-time approaches are consistent, and that the first-order-hold approximation demonstrates the best accuracy. 4.6.3    Comparison of models and approximation techniques All proposed approximation techniques work with the same VBR model, i.e. they share the same interfacing circuit of Figure ?2-2 and the same rotor state equations. This makes it simple to switch between the approaches. To get more insight into the achieved results, the eigenvalues of the power system of Figure ?4-4 have been calculated for different VBR formulations by linearizing the system around the operating point; the results are summarized in Table ?4?4. As it can be seen in Table ?4?4, unlike discrete-time 109  approximations, the continuous-time approximations add a large eigenvalue (~1,000) to the system, making it stiffer.  Despite the difference in numerical stiffness of models with discrete- and continuous-time approximation techniques, Table ?4?2 and Table ?4?3 show that all models used approximately the same number of time-steps. This is because in this specific case, the step size is predominantly determined by the maximum allowed time-step rather than by numerical stability constraints ?[55].  The filters are placed on the q-axis voltage (or current), which in the rotor reference frame becomes constant in steady-state. Therefore a reasonably good accuracy in general and zero steady-state error in particular is expected since the filters add no error in steady-state. When the original system is not stiff, it appears beneficial to use discrete-time approximations, since such methods have negligible additional computational cost and do not change the number of continuous state variables nor the numerical stiffness. The discrete-time approximation techniques also have the advantage of not requiring any parameter settings. The discrete-time filters might not produce good results in some situations, for example: when using stiff solvers ?[16] and/or if the time-step is too large for a given transient study when the approximated variables have high-frequency components; and if the dynamic saliency is considerable (e.g. 5.1/ ????? mdmq XX). 110  Table ?4?4  Eigenvalues of the System when using Different Approximation Techniques  (a) Continuous-Time Methods (b) Discrete-Time Methods Additional Winding ?[30] (Subsection ?3.3.1) +23.578 ? j39.426  (6.27 Hz) ?0.536 ?3.391 ?6.850 ? j15.790  (2.51 Hz) ?11.401 ?47.124 ?68.871 ? j36.410  (5.79 Hz) ?999.684 N/A Current Derivative Approximation (Subsection ?4.4.1) +23.450 ? j39.762 (6.33 Hz) ?0.544 ?3.480 ?7.096 ? j15.881  (2.53 Hz) ?11.401 ?47.124 ?68.713 ? j36.251 (5.77 Hz) ?999.782 +23.636 ? j40.363 (6.42 Hz) ?0.544 ?3.481 ?7.071 ? j15.863 (2.52 Hz) ?11.401 ?47.124 ?69.822 ? j36.658 (5.83 Hz) Voltage Algebraic Loop Relaxation (Subsection ?4.4.2) +26.231 ? j30.574  (4.87 Hz) ?0.563 ?3.737 ?8.843 ? j15.777 (2.51 Hz) ?11.321 ?47.124 ?69.196 ? j42.156 (6.71Hz) ?999.849 +26.215 ? j30.614 (4.87 Hz) ?0.560 ?3.735 ?8.821 ? j15.757 (2.51 Hz) ?11.321 ?47.124 ?69.214 ? j42.114 (5.70 Hz)  111  The continuous-time approximation approaches are suggested for cases when the original system is stiff from the beginning due to the natural components and time constants present in the system. In those cases, the added eigenvalue is controlled by selecting the filter pole, which then could be increased without making the overall system any stiffer while achieving the needed approximation accuracy or more.           112  CHAPTER 5:  DIRECT INTERFACING OF SYNCHRONOUS MACHINE MODELS FROM STATOR AND ROTOR TERMINALS The VBR models were originally introduced to facilitate the stator interface of synchronous machines to arbitrary networks ?[28]. The rotor interface is neglected in most VBR publications. In this chapter, two new models for synchronous machine with direct interface of rotor and stator terminals are introduced. In the first model, all of the rotor and stator windings are represented by an equivalent circuit. The whole model is represented only with constant-parameter RL impedance and voltage-source branches (VBR form). This model offers direct interfacing of any rotor and/or stator windings. In the second model, the damper windings are represented using standard state-space equations, while the field and stator windings remain in circuit form allowing a direct interface with external ac power network and dc exciter network, respectively. The approximation techniques used in this chapter are similar to those presented in Chapter 4. Similar to the previous chapter, the electromagnetic torque equations will not change and the same generic mechanical system will be considered (see subsection ?3.1.2). With respect to the state-of-the-art models that will be summarized in the following section, this chapter proposes new VBR models that combine the following characteristics: 113  1) Direct stator interfacing; 2) Direct rotor interfacing; and  3) Constant-parameter equivalent interfacing circuit. 5.1     Rotor and Stator Interfacing of Synchronous Machine Models The interfacing properties of various models for synchronous machines transient studies are summarized in Table ?5?1. As depicted in this table, the straightforward implementation of the CC-PD model (CC-PD) ?[1], ?[24], ?[26], ?[27] achieves a direct interface to arbitrary networks from the stator and rotor sides, but due to its rotor-position-dependent inductances, this model has variable parameters and is generally more computationally expensive.  The classical qd0 model possesses constant-parameter elements. Its stator circuit requires an indirect interface with external inductive networks (i.e. it requires snubbers), but its rotor field winding, if represented as a circuit, can be directly interfaced with the external dc exciter circuit. This model is therefore shown in Table ?5?1 as having a direct interface on the rotor side.   The VBR models ?[28] ? ?[31], ?[33], ?[34] have been derived to achieve a direct interface of the stator terminals with an arbitrary network on the ac side. The original/classical VBR formulation ?[28] possesses rotor-position-dependent stator interfacing resistance and114  Table ?5?1  Classification of Synchronous Machine Formulations based on Interfacing with External Inductive Networks Formulation Stator Interface Rotor Interface Circuit Parameters Classical CC-PD ?[1] Direct Direct Variable Classical qd0 ?[1], ?[7] Indirect Direct Constant Classical VBR ?[28], ?[29] Direct Indirect Variable Stator VBR ?[30], ?[31], Chapter 3, & 4 Direct Indirect Constant Stator VBR ?[33] Direct Indirect Variable Stator/Rotor VBR ?[34] Direct Direct Variable Chapter 5?s VBR Direct Direct Constant  inductance matrices. It was later shown in ?[29] how to obtain a constant resistance matrix. Overall, the VBR models presented in ?[28] and ?[29], and also the more recent models proposed in ?[33] and ?[34], have stator interfacing branches with variable parameters (similar to the CC-PD model). Simulation of variable elements is possible in some programs (such as PLECS ?[17] and ASMG ?[18]). However, it is very computationally costly and should be avoided if possible. The variable inductances of the stator interfacing circuit were first made constant by equating the sub-transient inductances (i.e. neglecting dynamic saliency) using: a) an 115  additional fictitious damper winding ?[30]; and b) singular perturbation ?[31] as discussed in subsections ?3.3.1 and ?3.3.2, respectively. A generalized VBR model and a unified interfacing circuit for both induction and synchronous machines was proposed in Chapter 2 and 3. Mathematical methods and several numerically advantageous continuous- and discrete-time approximation techniques are also proposed in Chapter 4 to obtain the same unified stator interfacing equivalent circuit of Chapter 2 and 3. All these models achieve constant-parameter interfacing of the machine's stator terminals as shown in Table ?5?1. The synchronous machine VBR models ?[28] ? ?[31], ?[33] and the ones presented in Chapter 3 and 4 were initially introduced to allow direct network-machine interfacing at the stator terminals, while assuming that a voltage-source representing the exciter is connected to the rotor field winding. In these models, to interface the rotor to an arbitrary network, e.g. an ac excitation systems ?[4], a controlled current-source is used. This also creates an incompatible interface if the field winding terminals are connected to inductive elements or a switching converter, therefore requiring a snubber. To obtain direct rotor interfacing, the model presented in ?[33] was recently extended to also represent the field winding using a VBR formulation ?[34]. Therein, the resulting equivalent interfacing circuit is defined by four coupled and variable RL branches, which achieves a direct interface of both stator and rotor terminals with external circuits as summarized in Table ?5?1.  116  5.2     All-Circuit Formulation The aim of this section is to derive a synchronous machine model consisting solely of basic circuit elements with constant parameters. Such model is easy to implement in traditional simulation packages while offering a direct interface to external circuits. The model can also be used to study special cases, e.g. broken damper windings and machines with excitable auxiliary windings ?[56]. 5.2.1    Stator voltage equations Inserting (?3?44) and (?3?45) into (?3?16) and (?3?17), respectively, and taking the derivative of the result gives qqsqqs piLpp ?? ?????? (?5?1) ddsdds piLpp ?? ??????  (?5?2) where the q- and d-axis sub-transient inductances are defined by (?3?34) and (?3?35). Then, using (?3?13) to (?3?15), the time derivatives of the rotor sub-transient flux linkages are written as ?[28] ????????Mj lkqjkqjkqjkqjmqq LirvLp1? (?5?3) ???????? ???????? ??Nj lkdjkdjkdjkdjlfdfdfdfdmdd LirvLirvLp1? (?5?4) 117  Inserting (?3?16), (?3?17), (?3?22), and (?3?23) into (?3?10) and (?3?11) and substituting (?5?1), (?5?2), (?5?3), and (?5?4) into the resulting equations, after algebraic manipulation, the q- and d-axis stator voltages are rewritten as qsqsqdsdrqssqs eiLpiLirv ???????? ? (?5?5) dsdsdqsqrdssds eiLpiLirv ???????? ? (?5?6) where the sub-transient voltages qse ?? and dse ??  are defined as ????????????????????Mj lkqjkqjkqjkqjmqNjkdjfdmdrqs LirvLiiLe11? (?5?7) ???????? ?????????? ????Nj lkdjkdjkdjkdjlfdfdfdfdmdMjkqjmqrds LirvLirvLiLe11? (?5?8) The q- and d-axis steady-state inductances, similar to (?3?34) and (?3?35), are defined as mqlsq LLL ?? (?5?9)  mdlsd LLL ?? . (?5?10)  As shown in Chapter 4, to obtain constant-parameter stator interfacing branches in stationary abc-phase coordinates, it is convenient to rearrange (?5?5) and (?5?6) as follows qsqsddsdrqssqs eiLpiLirv ??????????? ? (?5?11) dsdsdqsdrdssds eiLpiLirv ??????????? ? (?5?12) 118  where the new sub-transient voltages are  qsdqMj lkqjkqjkqjkqjmqNjkdjfdmddsddrqs iLLpLirvLiiLiLLe )()(11??????????????????????????????? ????? (?5?13) ?????? ????????????????????? ????Nj lkdjkdjkdjkdjlfdfdfdfdmdMjkqjmqqsdqrds LirvLirvLiLiLLe11)(? (?5?14) Using (?5?5), after algebraic manipulation, the derivative term qsdq iLLp )( ????? in (?5?13) can be expressed as )()( qsdsdrqssqsqdqqsdq eiLirvLLLiLLp ?????????????????? ? (?5?15) Equation (?5?13)  can then be rewritten by substituting qse ?? in (5-15) by (?5?7) and inserting the resulting equation into (?5?13), which simplifies to ???????????????????????????????????? ????Mj lkqjkqjkqjkqjqdmqNjkdjfdqdmddsqddrqs LirvLLiiLLiLLe11)1(? )( qssqsqdq irvLLL ????????  (?5?16) Finally, transforming (?5?11), (?5?12), and (?3?12) to abc-phase coordinates yields the constant-parameter stator interfacing circuit presented in Chapter 3 and 4 and shown in Figure ?2-2  and the upper part of Figure ?5-1. The parameters of the stator circuit are given by (?3?81), (?3?82), and (?3?84). 119  5.2.2    Rotor voltage equations Solving for stator currents qsi and dsi in (?3?16) and (?3?17), substituting the result in (?3?22) and (?3?23), after some algebraic manipulations, the magnetizing fluxes md?  are rearranged as ? ? qsmqlsmqMjkqjmqlsmq LLLiLL ?? ??? ??1|| (?5?17)  ? ? dsmdlsmdNjkdjfdmdlsmd LLLiiLL ?????????????? ??1|| (?5?18)  where ? ? 111|| ??? ?? baba LLLL . Inserting (5-17) and (5-18) into (?3?19) to (?3?21) and taking the derivative of the resulting equations gives ? ? MjpLLLiLLiLpp qsmqlsmqMakqamqlskqjlkqjkqj ,,2,1,||1?????????? ?? ???? (?5?19)  ? ? dsmdlsmdNakdafdmdlsfdlfdfd pLLLiiLLiLpp ????????????????????? ??1|| (?5?20)  ? ? NjpLLLiiLLiLpp dsmdlsmdNakdafdmdlskdjlkdjkdj ,,2,1,||1??????????????????? ???? (?5?21)  Incorporating (?3?10), (?3?11), and (5-19) ? (5-21) into (?3?13) ? (?3?15), after algebraic manipulation, yields 120  MjeiLLpipLirv qrMakqamqlskqjlkqjkqjkqjkqj ,,2,1,)||(1???????? ?? (?5?22)  drNakdafdmdlsfdlfdfdfdfd eiiLLpipLirv ????????? ???? ??1)||( (?5?23)  NjeiiLLpipLirv drNakdafdmdlskdjlkdjkdjkdjkdj ,,2,1,)||(1??????????? ???? ?? (?5?24)  where )( dsrqssqsqmqqr irvLLe ??????? (?5?25)  )( qsrdssdsdmddr irvLLe ???????. (?5?26)  In summary, the all-circuit constant-parameter voltage-behind-reactance machine model, herein referred to as AC-CP-VBR, is shown in Figure ?5-1. The constant-parameter circuit representation based on (5-22) to (5-26) is shown in Figure ?5-1 (lower part). The sub-transient algebraic voltage equations are calculated by (?5?16) [or (?5?13)], (?5?14), (?5?25), and (?5?26) in addition to qs? and ds?  given by (?3?16), (?3?17), (?3?22), and (?3?23). As shown in Figure ?5-1, the voltage terminals of the damper windings are available for connection with arbitrary external circuits. The calculation of qse ??? using (?5?13) will require the approximation of a time derivative (see Chapter 4, section ?4.4.1). Alternatively, if qse ??? is calculated using (?5?16), the resulting model 121   Figure ?5-1  All-circuit constant-parameter VBR synchronous machine model (AC-CP-VBR) with direct interfacing to arbitrary external ac and dc networks. 122  is algebraically identical to the classical qd0 model as no approximations are made. However, this may introduce an algebraic constraint (loop) in the q-axis, which needs to be relaxed to achieve a numerically efficient explicit model (see Chapter 4, section ?4.4.2). Several possible discrete- and continuous-time methods to approximate the time derivative and break such algebraic loops have been presented in Chapter 4, section ?4.4. 5.3     Stator-and-Field-Circuit Formulation In most studies, the damper windings of synchronous machines are all short-circuited. As a result, they are not interfaced with external networks, and as such there is no numerical advantage gained by representing them using circuit elements. In fact, representing the dynamics of all damper windings in a standard state-space form with flux linkages as state variables (as opposed to the AC-CP-VBR model presented earlier in section ?5.2) further increases the model efficiency, as will be shown in section ?5.5. 5.3.1    Stator voltage equations After manipulating (?3?19) and (?3?21), the damper winding currents can be written as MjLi lkqjmqkqjkqj ,,2,1, ???? ?? (?5?27)  NjLi lkdjmdkdjkdj ,,2,1, ???? ?? (?5?28)  123  Substituting (?5?27) and (?5?28) into (?3?13) and (?3?15), respectively, the damper windings state equations become  MjvLrp kqjkqjmqlkqjkqjkqj ,,2,1,)( ????? ??? (?5?29)  NjvLrp kdjkdjmdlkdjkdjkdj ,,2,1,)( ????? ??? (?5?30)  Substituting (?5?28) into (?3?23), and after algebraic simplification, the d-axis magnetizing flux is written as ?????????????? ??Nj lkdjkdjfddsmdmd LiiL1?? (?5?31)  where 1111?? ????????????? ?Nj lkdjmdmd LLL. (?5?32)  Inserting (?3?39) and (?5?29) into (?5?1) also gives ?????????Mj lkqjkqjmqkqjmqqsqqs LrLiLpp12)( ???. (?5?33)  Substituting (?3?17) into (?3?10), replacing md?  with (?5?31), and inserting (5?33) into the resulting equation give the q-axis stator voltage equation 124  qsqsqdsdrqssqs eiLpiLirv ??????????? ? (?5?34)  where  mdlsd LLL ????????  (?5?35)  and ???????????????????????Mj lkqjkqjmqkqjmqNj lkdjkdjfdmdrqs LrLLiLe121)( ????. (?5?36)  To find the d-axis stator voltage equation, (?3?40) is first substituted into (?5?2). The flux linkage time derivatives are then replaced by using (?3?14) and (?5?30), which yields ???????? ????????? ??Nj lkdjkdjmdkdjlfdfdfdfdmddsdds LrLirvLiLpp12)( ???. (?5?37)  Inserting (?3?16) and (?3?44) into (?3?11), and substituting the resulting equation into (?5?37), after some algebraic manipulation, the d-axis stator voltage equation becomes dsdsdqsqrdssds eiLpiLirv ?????????? ? (?5?38)  where ????????????????Nj lkdjkdjmdkdjmdlfdfdfdfdmdqrds LrLLirvLe12)()( ????. (?5?39)  125  Following the approach set forth in Chapter 4, section ?4.2 and in this chapter, subsection ?5.2.1, the stator voltage equations (?5?34) and (?5?38) are rearranged as qsqsddsdrqssqs eiLpiLirv ??????????? ? (?5?40)  dsdsdqsdrdssds eiLpiLirv ??????????? ?, (?5?41)  where the new sub-transient voltages are qsdqMj lkqjkqjmqkqjmqNj lkdjkdjfdmddsddrqs iLLpLrLLiLiLLe )()()(121????????????????????????????????????? ???????? (?5?42)  ? ?lfdfdfdfdmdNj lkdjkdjmdkdjmdqqsdqrds LirvLLrLiLLe????????????????????? ??12)()( ???? (?5?43)  Similar to AC-CP-VBR, in order to replace the term qsdq iLLp )( ????? in (?5?42) with state variables and voltages, qsqiLp ?? is first expressed using (?5?34) and (?5?36) and then inserted into (?5?42). After some algebraic manipulations, (?5?42) becomes ???????????????????????????????????????????? ????Mj lkqjkqjmqkqjqdmqNj lkdjkdjfdqdmddsqqddrqs LrLLLiLLiLLLLe121)()(???? ))(( qssqsqdq irvLLL ???????? . (?5?44)  126  When converted to abc-coordinates, (?5?40) and (?5?41) along with (?3?12) form the same constant-parameter stator interfacing circuit as the one in the upper part of Figure ?5-1. The parameters of this interfacing circuit are also defined by (?3?81), (?3?82), and (?3?84). 5.3.2    Rotor voltage equations The state equations for the q- and d-axis rotor damper windings are given by (?5?29) and (?5?30), respectively. Therein, mq? and md?  are computed using (?3?44) and (?3?39), and (?5?31), respectively. The implementation of the damper windings state model is shown in Figure ?5-2 (left side). To allow for a direct interface with external networks, i.e. exciter circuit, the field winding has to be represented as a circuit with input terminals. To obtain the field circuit, (?3?17) is solved for ids, and then the result is substituted into (?5?31) which following algebraic manipulation gives ??????????????? ??Nj lkdjkdjlsdsfdmdmd LLiL1??? (?5?45)  where 11111?? ??????????????? ?Nj lkdjlsmdmd LLLL. (?5?46)  Substituting (?5?45) into (?3?20), taking the derivative of the results, and inserting (?5?30) give 127  ????????????????Nj lkdjkdjmdkdjmddslsmdfdlfdfd LrLLLpiLpp12)( ???? (?5?47)  where mdlfdlfd LLL ????????. (?5?48)  In the next step, (?3?11) is solved for dsp?  and substituted in (?5?47), and the resulting equation is then inserted into (?3?14). This procedure, after some algebraic simplifications, gives  fdfdlfdfdfdfd eiLpirv ??????? (?5?49)  where ???????? ??????????? ??Nj lkdjkdjmdkdjlsqsrdssdsmdfd LrLirvLe12)( ????. (?5?50)  The resulting constant-parameter field winding interfacing circuit is depicted in Figure ?5-2 (right side). The complete model is obtained by replacing the rotor circuit in Figure ?5-1 (bottom) with Figure ?5-2 and updating the sub-transient voltage equations. For this model, the sub-transient voltages are given by (?5?44) [or (?5?42)], (?5?43), and (?5?50). The algebraic equations (?3?44), (?3?39), and (?5?31) are used to compute the magnetizing fluxes. This 128   stator-and-field-circuit constant-parameter VBR model is herein referred to as SFC-CP-VBR. The SFC-CP-VBR model is algebraically equivalent to the classical qd0 model when qse ??? is calculated using (?5?44). However, an algebraic loop may be introduced when the model is interfaced to an external network. This algebraic loop may be relaxed to make an efficient explicit model. Similarly to AC-CP-VBR model, if qse ??? is calculated using (?5?42), a current derivative must be approximated, as will further be explained in section ?5.4.   Figure ?5-2  Rotor subsystem for the stator-and-field-circuit constant-parameter voltage-behind-reactance (SFC-CP-VBR) model wherein the damper windings are represented as a state model and the field winding is made available as an interfacing circuit (the stator interfacing circuit is the same as in Figure ?5-1). 129  5.4     Numerically Efficient Explicit Implementation The AC-CP-VBR and SFC-CP-VBR models have constant parameters. However, when connected to inductive external networks, these models will contain algebraic loops in both the d- and q-axis. The algebraic loops are illustrated in Figure ?5-3 for the AC-CP-VBR formulation. Similar algebraic loops exist for the SFC-CP-VBR model. An algebraic loop of similar nature was also encountered in Chapter 4, section ?4.2 in the q-axis.  Models containing algebraic loops typically require additional iterations within each time-step that considerably add to the computational cost ?[16]. To achieve a numerically efficient solution (with acceptable accuracy but with less computational cost), it is possible to break the algebraic loops using discrete- or continuous-time low-pass filters (see Chapter 4, subsection ?4.4.2). For example, a first-order low-pass filter (LPF) may be considered as given by (?4?26). By increasing the pole?s magnitude, the filter becomes faster, improving the accuracy of the approximation. However, a very fast filter can make the system numerically stiff. Inserting a low-pass filter anywhere in the loop will relax the algebraic constraint. However, a more accurate solution can be obtained when the approximated variable is varying slowly. Herein, qsv and fdv are assumed to be the slowest variables in the q- and d-axis loops, respectively (see Figure ?5-3). Therefore, in the proposed models, the low-pass filters are inserted as to relax these variables accordingly. 130   Figure ?5-3  Algebraic loops in AC-CP-VBR resulting in an implicit formulation. The H blocks indicate where the low-pass filters may be inserted to break the algebraic loops. The algebraic loop in the q-axis of AC-CP-VBR (see Figure ?5-3, left) does not exist if qse ??? is evaluated using (?5?13) [or (?5?42) for the SFC-CP-VBR model], which requires the current derivative qspi. In this case, the derivative can be approximated using a high-pass filter (see Chapter 4, subsection ?4.4.1). 5.5     Computer Studies To verify the proposed models, a test system consisting of a large synchronous generator with a static exciter, as shown in Figure ?5-4, is considered. Without loss of generality, in this case-study, a simplified model of a potential-source controlled-rectifier excitation system [5, see sec. 8.3] is considered, wherein the ac side of the exciter is represented by an inductive Th?venin equivalent. Further details of the excitation system are omitted, as the 131  focus of the studies is on the accuracy of various machine models in predicting the electrical transients on the stator and rotor terminals of the machine. The system parameters are summarized in Appendix F. A single per-unit system based on the machine?s nominal power and stator voltage is used for the generator and all external circuits. The steam turbine generator is assumed to be driven by a constant 0.8 pu mechanical torque in steady-state. The Th?venin voltage-source on the stator side has a frequency of 60 Hz with nominal voltage. On the exciter side, a 12-pulse rectifier fed from a three-phase 60 Hz voltage-source generates a field current of about 1.325 pu. The generator is assumed to be grounded through a grounding resistor (rg = 0.2 pu) as depicted in Figure ?5-4. The test system of Fig. 4 has been implemented using MATLAB/Simulink ?[15], ?[16] and the PLECS toolbox ?[17]. For a numerically efficient and algebraic-loop-free implementation of the AC-CP-VBR and SFC-CP-VBR models, the first-order low-pass filters with poles of ?5,000 and ?1,000, respectively, have been used. For conciseness, only these voltage filters are herein considered. However, similar to Chapter 4, discrete-time filters are applicable here as well. The existing constant-parameter models, namely qd0 ?[1] and the stator constant-parameter VBR (CP-VBR) with voltage loop relaxation  (Chapter 4), were also implemented for comparison. Since PLECS? built-in qd0 model does not include zero-sequence and for consistency with the other models, a custom qd0 model has been created using 132  conventional circuit elements. A circuit-based implementation of the qd0 model (as opposed to the more traditional state-space model) is required for direct connection of its field winding with the exciter. Such implementation is shown in [1, see p. 202]. The stator is also implemented as a circuit which is interfaced to the external power system using controlled current-sources ?[8]. Since the external ac network (its Th?venin equivalent) is inductive, the qd0 model requires a three-phase interfacing snubber, as shown in Figure ?5-4. To have reasonable accuracy, a 50 pu stator snubber is chosen for the qd0 model. This value of snubber results in a system eigenvalue with a magnitude of 2.61?105, as shown in Table ?5?2.   Figure ?5-4  A wye-grounded steam turbine generator with a static 12-pulse rectifier-based exciter system. The stator snubbers are required only for the classical qd0 model and the field snubber is required only for the CP-VBR model from Chapter 4. 133  The CP-VBR model (see Chapter 4, section ?4.2) has a direct interface with the network through its stator terminals, but its field winding is represented as a standard state-space model with a voltage-input current-output formulation. It is therefore interfaced to the external rotor circuit using controlled current-sources. In this case, a resistive (or capacitive) snubber must be connected to the field winding, as shown in Figure ?4-4. To have reasonable accuracy, a 0.030 pu snubber (50 times larger than the field winding resistance) is chosen for the CP-VBR model of Chapter 4. This introduces an eigenvalue with a magnitude of 1.35?106. Without the snubber, the largest eigenvalue magnitude of the system would be smaller than 1,000. Finally, the CC-PD model ?[1] is used to generate a reference solution. This model requires no snubber at either sets of terminals, but it has variable coupled inductances resulting in significant computational cost. To compare the accuracy of all the subject models, similar to the other chapters, the 2-norm cumulative relative error is considered which is given by (?2?64).  Since the CC-PD, AC-CP-VBR, and SFC-CP-VBR models are not very stiff, MATLAB?s default solver ode45 has been used for simulating these models. The qd0 and CP-VBR models require snubbers and therefore are stiff. For these two models, the stiffly-stable solver ode15s is used instead. A list of MATLAB ordinary differential equation (ODE) solvers is given in Appendix B. For consistency among all simulation studies, the relative and absolute error tolerances are set to 10?4, and the maximum and minimum time-steps are set to 10?3 and 10?7 s, respectively. To produce an accurate reference solution, the CC-PD 134  model was used with a very small time-step of 10?6 seconds and both error tolerances were set to 10?6. For consistency, all simulations were run on a personal computer with a 2.53GHz Intel CPU and Windows XP operating system. In all studies, the system is assumed to operate in steady-state prior to a transient event. 5.5.1    Single-phase-to-ground fault in the network To emulate a severe unbalanced condition at the machine stator terminals, a single-phase-to-ground fault is applied at the source. The fault is implemented by setting the phase a voltage to zero at the end of the first cycle ( 60/1?t s). The corresponding source voltages abcv  and currents abci , the machine grounding current ngi , the three-phase exciter current abcexi , the rectified exciter current dci , and the electromagnetic torque eT  predicted by the various models are presented in Figure ?5-5. As it can be seen in Figure ?5-5, all the models give visibly close results. To provide a better comparison, the magnified views of abci , abcexi , and dci  are plotted in Figure ?5-6 to Figure ?5-9.  A summary of the models numerical performance is also given in Table ?5?2. Analyzing Figure ?5-6 to Figure ?5-9 and the second column of Table ?5?2 shows that the stiff models (i.e. qd0 and CP-VBR) need considerably more time-steps (2,384 and 2,739) than the non-stiff models (i.e. CC-PD, AC-CP-VBR, and SFC-CP-VBR) (less than 1,000). The CC-PD model also uses a small number of steps (922) while producing the most accurate results. However, due to its variable parameters, this model is by far the slowest (6.772 s). It was135   Figure ?5-5  Single-phase-to-ground fault case-study transient responses as predicted by the considered models. From top to bottom: bus voltages vabc and currents iabc, machine neutral current ing, three-phase ac exciter current iabcex, rectifier current idc, and electromagnetic torque Te. also determined that using the explicit solver ode45 for the qd0 and CP-VBR models (which are stiff) leads to significantly more time-steps and longer simulation times.  136   Figure ?5-6  Magnified fragment from Figure ?5-5: phase current ic (in steady-state).   Figure ?5-7  Magnified fragment from Figure ?5-5: phase current ic (during transient). Figure ?5-6 shows that as opposed to the other models, qd0 predicts steady-state stator currents with a noticeable error. Similarly, as it can be seen in Figure ?5-8 and Figure ?5-9, the stator CP-VBR model (Chapter 4, section ?4.2) also has a noticeable error in field winding current. Based on these figures and Table ?5?2 (third and fourth columns from left), 137   Figure ?5-8  Magnified fragment from Figure ?5-5: ac exciter phase current icex.   Figure ?5-9  Magnified fragment from Figure ?5-5: rectified current idc. it can be concluded that the presence of snubbers in the stator or rotor side increases the errors in iabc or iabcex, respectively. The two fastest models are AC-CP-VBR and SFC-CP-VBR, both of which do not require snubbers and have direct stator and field windings interfacing. These two models use the138  Table ?5?2  Numerical Performance of the Models for the Single-Phase-to-Ground Fault Study Model Number of Time-Steps Error ?avg Largest Eigenvalue Magnitude max{|?|} Simulation Time Stator Source Current iabc Exciter Source Current iabcex Coupled-Circuit Phase-Domain  (CC-PD) 922 0.03 % 0.01 % 831 6.772 s Classical QD0 With Stator Snubber (qd0) 2,384* 2.63 % 0.08 % 2.61?105 0.642 s Stator CP-VBR (Chapter 4) With Rotor Snubber 2,739* 1.97 % 1.52 % 1.35?106 0.954 s All-Circuit (AC-CP-VBR) 969 2.51 % 0.39 % 5,000 0.619 s Stator-And-Field-Circuit (SFC-CP-VBR) 875 1.95 % 0.23 % 1,000 0.547 s * For the numerically stiff qd0 and CP-VBR models, MATLAB?s ode15s solver is used. The other models are solved with ode45.  139  low-pass filter (?4?26) to relax algebraic loops. As Table ?5?2 shows, the SFC-CP-VBR model offers better accuracy since its damper windings are implemented as a state-space model and do not contribute directly to the sub-transient voltage dre ??  algebraic loop. Therefore, low-pass filter poles with small magnitudes can be chosen (~1,000) which merely affect the system?s stiffness. At the same time, the AC-CP-VBR model requires a filter pole with a much larger magnitude (~5,000) to achieve a comparable accuracy. Having the smallest eigenvalues, in the case studies with lower frequencies, the SFC-CP-VBR model would be able to choose significantly larger time-steps than the other models. 5.5.2    Diode failure in the exciter system When system disturbances (such as faults) in the stator circuit of a synchronous generator happen, the field current fluctuates violently and can become negative and threaten the transient blocking voltage of rectifiers of a static exciter. A computer model can be used for prediction of over-voltages and over-currents and to design the required remedies ?[49]. In case of brushless exciters, testing the diode rectifiers may not be easy. However, detecting diode failure from analyzing exciter waveforms is possible ?[34]. Designing such detection system also depends on computer simulations. To validate the proposed models in these situations, where the field winding and exciter variables are of particular interest, a diode failure scenario is simulated.  A three-leg diode bridge is shown in Figure ?5-10 (left). In this study, the diode D3 in the wye transformer rectifier (see Figure ?5-4) fails and becomes short-circuited when the 140  reverse voltage is at its peak (at 0153.0?t s). This exposes D3 and D1 to excessive currents (more than 100 pu), which burn both of them or their protective fuses ?[50]. It is assumed that the current continues to flow until it goes to zero, after which the diodes become open. Figure ?5-10 (right) shows the simulated diode currents for this scenario. The times of the equivalent short-circuit and open-circuit events are also shown in Figure ?5-10  (right) by the circled numbers 1 and 2, respectively.  Simulation results of the exciter and field variables predicted by various models are presented in Figure ?5-11. This figure shows the delta transformer (secondary) voltages vT? and currents iT?, the wye transformer (tertiary) voltages vTY  and currents iTY, and the rectifier output voltage vfd and current idc. The stator voltages and currents for this study are not affected much and therefore are not shown. A summary of the numerical performance of all considered models is given in Table ?5?3. As seen in the second and fifth columns of Table ?5?3 (from left), the models proposed in this chapter need fewer time-steps and are faster than the stiff models with snubbers. The CC-PD model requires the fewest time-steps (668), but is the slowest due to its variable inductances. Table ?5?3 (third and fourth columns from left) and Figure ?5-11 indicate that all considered models have excellent accuracy except for the CP-VBR model. The resistive snubber of CP-VBR model (see Figure ?5-4) adds a noticeable current to the rectifier current (dci ), as seen in Figure ?5-11 (last subplot). 141   Figure ?5-10  Left: the diode bridge connected to the delta transformer in Figure ?5-4. Right: simulation results of the diode currents in pu: 1) after the initial failure (short-circuit) of D3; and 2) after D1 and D3 (or the corresponding protective diodes) become open due to the ensuing excessive currents. At the same time, the qd0 model with the stator interfacing snubber accurately predicts the field and exciter quantities since its field winding is directly interfaced and the three-phase snubber depicted in Figure ?5-4 mostly affects its stator variables. This study demonstrates that the models with direct interface of the rotor field winding predict the rotor-exciter variables more accurately in steady-state and during transients.  142   Figure ?5-11  Transient responses as predicted by the considered models for the static exciter diode failure case-study. From top to bottom: exciter delta transformer line voltages vT?, delta transformer currents iT?, wye transformer line voltages vTY, wye transformer currents iTY, rectifier output voltage vfd, and rectifier current idc.  143  Table ?5?3  Numerical Performance of the Models for the Exciter Diode Failure Study Model Number of  Time-Steps Error ?avg Simulation Time Stator Source Current iabc Exciter Source Current iabcex Coupled-Circuit Phase-Domain  (CC-PD) 668 0.02 % 0.01 % 5.198 s Classical QD0 With Stator Snubber (qd0) 1,966* 0.02 % 0.09 % 1.069 s Stator CP-VBR (Chapter 4) With Rotor Snubber 2,243* 1.59 % 0.63 % 1.115 s All-Circuit (AC-CP-VBR)  724 0.04 % 0.01 % 0.802 s Stator-And-Field-Circuit (SFC-CP-VBR) 695 0.02 % 0.01 % 0.757 s * For the numerically stiff qd0 and CP-VBR models, MATLAB?s ode15s solver is used. The other models are solved with ode45.  144  CHAPTER 6:  SUMMARY OF CONTRIBUTIONS AND FUTURE WORK 6.1     Conclusions and Contributions In this thesis, initially a general constant-parameter interfacing circuit for synchronous and induction machines models with voltage-behind-reactance (VBR) form was developed. Based on the interfacing circuit, several formulations for ac machines were presented to obtain numerically efficient VBR models. The proposed formulations have direct interface and are demonstrated to be more accurate and numerically efficient than the existing classical qd0 models, and are easier to implement into different simulation programs as well. A new approach for obtaining constant-parameter models of synchronous machines was proposed. This methodology is used to achieve constant-parameter VBR interface for both the rotor and stator terminals of the machine. With respect to the initial objectives of this research, the contributions of this thesis can be summarized as follows: 6.1.1     Objective 1 In Chapter 2, the model presented in ?[32] was extended by including the zero-sequence in the interfacing circuit. The new model has constant-parameters, is explicit, and is easy to implement in commonly used simulation programs. The computer study showed that the 145  proposed model requires the least number of calculations (764) compared to the state-of-the-art VBR model (3,865) and the classical qd0 model (7,070) for the same transient study and similar level of accuracy, which represents a substantial improvement in simulation efficiency.  This contribution addresses Objective 1 of this thesis. 6.1.2     Objective 2 In Chapter 3, the state-of-the-art methods ?[30], ?[31] for obtaining the constant-parameter synchronous machine VBR model were used to formulate a new synchronous machine  model with the same unified interfacing circuit as was developed for the induction machines in Chapter 2. It was demonstrated that the new explicit model has very good accuracy (1.2%, for the case study) and is faster than the classical qd0 model (number of time steps 557 versus 1,871, respectively).  This result addresses the Objective 2 of this thesis. It was also shown that the singular perturbation method ?[31] is actually algebraically equivalent to the classical qd0 formulation. However, the model obtained using the singular perturbation is in fact implicit and therefore may become computationally very expensive to solve and difficult to implement in many commonly used simulation programs. 6.1.3     Objective 3 Chapter 4 extended the work presented in Chapter 3 and presented new approximation approaches to obtain the same constant-parameter VBR interfacing circuit for synchronous 146  machine models. The proposed approach artificially makes the machine's sub-transient reactance equal while maintaining all properties of the original machine model up to the desired level of accuracy specified by the user. This new methodology is based on either continuous or discrete filters that can be used to approximate either current derivatives or the voltage in order to break algebraic loops. The proposed continuous- and discrete-time approximations result in computationally efficient implementation of the proposed models (see section ?4.4). The new approximation methods (while achieving the same goal as the added fictitious damper winding) are very intuitive and easy to implement, do not require complex analysis or fitting procedures, and are shown to give a high degree of accuracy. As can be seen in section ?4.6, the cumulative error remains below 1% at a fairly large time-step of mst 1?? . This accuracy is acceptable for most power systems transient studies. This chapter satisfies the Objective 3 of this thesis. 6.1.4     Objective 4 Chapter 5 presents two new explicit synchronous machine models (AC-CP-VBR, section ?5.2) and (SFC-CP-VBR, section ?5.3) with constant-parameters and direct interfacing of the machine's stator and rotor terminals with external circuits. The proposed models are considerably less stiff than the conventional indirectly-interfaced qd0 and constant-parameter stator VBR models (Chapter 4, section ?4.2). The AC-CP-VBR model goes beyond just the state-variable-based simulation programs since it is represented solely by basic/conventional and constant circuit elements. This model is easy to implement in different simulation programs and is practical for special cases where the equivalent 147  damper windings may not be short-circuited or additional field windings are required. In the second model, SFC-CP-VBR, to increase numerical accuracy and efficiency, the short-circuited damper windings are represented in state-space form (with flux linkages as state variables). Single-phase-to-ground fault and the exciter rectifier diode failure case studies have been conducted to validate and compare the proposed models. It is shown that both models have an excellent combination of accuracy and numerical efficiency. The new models are demonstrated to be about 10 times faster than the CC-PD model and achieve similar accuracy with qd0 models, while having a few orders of magnitude lower stiffness ratio. These results completely address the last objective of this thesis (Objective 4). 6.2     Potential Impacts of Contributions It is further envisioned that the proposed constant-parameter VBR interfacing circuits and approximation techniques will find wide application in many state-variable-based transient simulation programs. Having the same general structure and interfacing circuit for both synchronous and induction machines will make it easier to develop customized electrical machines components, their user interfaces and parameter entry. Software developers and users should be able to easily implement such models in their programs. Therein, depending on the solver selection and the machine/system parameters, the user can select the type of approximation that gives the best numerical results. The proposed fast and accurate machine models will enables the users to model larger systems with more details 148  that is not possible or practical with the present tools and will save countless number of hours spend on running simulations by many power engineers, researchers, and students.  6.3     Future Work The research directions listed below are already being considered by the members of Electric Power and Energy System research group at the University of British Columbia.  6.3.1    Doubly-fed induction machine model with direct interface A stator-fed induction machine (with squirrel cage rotor or when the rotor circuit is connected to voltage-sources) has a simple VBR representation shown in Chapter 2. However, having direct interface for both rotor and stator terminals is required for doubly-fed induction machines. The methods proposed in Chapter 4 and 5 can be investigated for development of a new model for induction machine with constant-parameter decoupled interfacing circuit for both rotor and stator terminals. 6.3.2    Inclusion of magnetic saturation For more accurate prediction of power system transients, the effect of magnetic saturation should be considered. Magnetic saturation results in change of inductances (especially magnetizing inductances of the q- and d-axis) as the magnetizing flux changes. The non-linear (variable) part of the magnetizing inductances can be also incorporated in the sub-149  transient voltage-source using similar approximation approaches as was proposed in Chapter 4. 6.3.3    Application to dynamic phasor solution and shifted frequency analysis  Simulation of power system using dynamic phasors approach and shifted frequency solution has been considered in ?[57] ? ?[60]. More specifically, using the VBR formulation with variable parameters has been investigated recently in ?[61], ?[62] with promising results. However, it should be relatively straightforward to extend the proposed constant-parameter interfacing equivalent circuits to the dynamic phasor and shifted frequency solution approaches. Therein, additional computational saving can be achieved due to constant and decoupled interfacing circuit parameters of the models presented in this thesis.     150  References [1] P. C. Krause, O. Wasynczuk, and S. D. Sudhoff, Analysis of Electric Machinery and Drive Systems, 2nd Edition, Wiley-IEEE Press, Piscataway, NJ, 2002. [2] B. Wu, Y. Lang, N. Zargari, S. Kouro, Power Conversion and Control of Wind Energy Systems, Wiley-IEEE Press, Piscataway, NJ, 2011. [3] L. L. Lai, T. F. Chan, Distributed Generation: Induction and Permanent Magnet Generators, Wiley-IEEE Press, Piscataway, NJ, 2007. [4] P. Kundur, Power System Stability and Control, McGraw-Hill, New York, 1994. [5] E. W. Kimbark, Power System Stability, Vol. III, IEEE Press, 1995. [6] P. M. Anderson, A. A. Fouad, Power System Control and Stability, IEEE Press, Power Systems Engineering Series, 1994. [7] R. H. Park, "Two-reaction theory of synchronous machines generalized method of analysis-part I," Transactions of the American Institute of Electrical Engineers, vol. 48, no. 3, pp. 716-727, Jul. 1929. [8] L. Wang, J. Jatskevich, V. Dinavahi, H. W. Dommel, J. A. Martinez, K.  Strunz, M. Rioual, G. W. Chang, R. Iravani, ?Methods of Interfacing Rotating Machine Models in Transient Simulation Programs,? IEEE Transactions on Power Delivery, vol. 25, no. 2, pp. 891-903, Apr. 2010. [9] H. W. Dommel, EMTP Theory Book, Vancouver, BC, Canada: Microtran Power System Analysis Cooperation, May 1992. [10] ?Alternative Transients Programs, Rule Book,? ATP?EMTP, ATP User Group, OR, 1995. 151  [11]  ?Transients Analysis Program for Power and Power Electronic Circuits, Reference Manual,? MicroTran Power System Analysis Corp., BC, Canada, Jan. 2002. Available: http://www.microtran.com [12] ?Power Systems Computer Aided Design, PSCAD, User?s Guide,? Manitoba HVDC Research Center, MB, Canada, Feb. 2010. Available: http://www.hvdc.ca [13] ?Electromagnetic Transient Program Restructured Version, EMTP-RV,? CEATI Int., Inc., QC, Canada, 2005. Available: http://www.emtp.com [14] ?PSIM User?s Guide,? Version 7.0, Powersim Inc., USA, Jul. 2006. Available: http://powersimtech.com [15] ?MATLAB 7 (R2009a)?, MathWorks, Inc., Natick, MA, 2009. [16] ?Simulink Dynamic System Simulation Software?Users Manual?, MathWorks, Inc., Natick, Massachusetts, 2009. Available: http://www.mathworks.com [17] ?Piecewise Linear Electrical Circuit Simulation (PLECS)?User Manual?, Version. 1.5, Plexim GmbH, 2006. Available: http://www.plexim.com [18] ?Automated State Model Generator (ASMG)?Reference Manual?, Version 2, PC Krause and Associates Inc., IN, 2002. Available: http://www.pcka.com [19] ?Simulink Dynamic System Simulation Software?SimPowerSystems Manual?, The MathWorks, Inc., Natick, MA, The MathWorks, Inc., 2009. [20] ?Modelica - A Unified Object-Oriented Language for Systems Modeling?Language Specification?, Version 3.3, Modelica Association, May 9, 2012. Available: http://www.modelica.org [21] ?acslX?User?s Guide?, Version 2.5, The AEgis Technologies Group, Inc., March 2009. Available: http://www.acslx.com 152  [22] ?EASY5 Engineering Software for the Design, Analysis and Simulation?, MSC SimEnterprise, Inc. Available: http://www.mscsoftware.com [23] ?EUROSTAG: Software for simulation of large electric power systems?, Tractebel Energy Engineering. Available: http://www.eurostag.be and www.tractebel-engineering.com [24] P. Subramaniam and O. P. Malik, ?Digital simulation of a synchronous generator in the direct-phase quantities,? Proceedings of the Institute of Electrical Engineers, vol. 118, no. 1, pp. 153-160, Jan. 1971. [25] G. R. Slemon, ?An equivalent circuit approach to analysis of synchronous machines with saliency and saturation,? IEEE Transactions on Energy Conversion, vol. 5, no. 3, pp. 538?545, Sept. 1990. [26] O. Wasynczuk, S. D. Sudhoff, "Automated state model generation algorithm for power circuits and systems," IEEE Transactions on Power Systems, vol. 11, no. 4, pp. 1951-1956, Nov. 1996. [27] X. Cao, A. Kurita, H. Mitsuma, Y. Tada and H. Okamoto, ?Improvements of numerical stability of electromagnetic transient simulation by use of phase-domain synchronous machine models,? Electrical Engineering in Japan, vol. 128, no. 3, Apr. 1999. [28] S. D. Pekarek, O. Wasynczuk, H. J. Hegner, "An efficient and accurate model for the simulation and analysis of synchronous machine/converter systems," IEEE Transactions on Energy Conversion, vol. 13, no. 1, pp. 42-48, Mar. 1998. [29] W. Zhu, S. D. Pekarek, J. Jatskevich, O. Wasynczuk, D. Delisle, "A model-in-the-loop interface to emulate source dynamics in a zonal DC distribution system," IEEE Transactions on on Power Electronics, vol. 20, no. 2, pp. 438-445, Mar. 2005. [30] S. D. Pekarek, E. A. Walters, "An accurate method of neglecting dynamic saliency of synchronous machines in power electronic based systems," IEEE Transactions on Energy Conversion, vol. 14, no. 4, pp. 1177-1183, Dec. 1999. 153  [31] S. D. Pekarek, M. T. Lemanski, E. A. Walters, "On the use of singular perturbations to neglect the dynamic saliency of synchronous machines," IEEE Transactions on Energy Conversion, vol. 17, no. 3, pp. 385- 391, Sep. 2002. [32] L. Wang, J. Jatskevich, S. D. Pekarek, "Modeling of Induction Machines Using a Voltage-Behind-Reactance Formulation," IEEE Transactions on Energy Conversion, vol. 23, no. 2, pp. 382-392, Jun. 2008. [33] D. C. Aliprantis, O. Wasynczuk, C. D. Rodriguez Valdez, "A Voltage-Behind-Reactance Synchronous Machine Model With Saturation and Arbitrary Rotor Network Representation," IEEE Transactions on  Energy Conversion, vol. 23, no. 2, pp. 499-508, Jun. 2008. [34] A. M. Cramer, B. P. Loop, D. C. Aliprantis, "Synchronous Machine Model With Voltage-Behind-Reactance Formulation of Stator and Field Windings," IEEE Transactions on Energy Conversion, vol. 27, no. 2, pp. 391-402, Jun. 2012. [35] (2012, Sep. 28). Handbook of General Requirements for Electrical Service to Dispersed Generation Customers, Consolidated Edison Company of New York, Inc., New York, NY, Mar. 2006. [Online]. [36] L. Wang and J. Jatskevich, ?A voltage-behind-reactance synchronous machine model for the EMTP-type solution,? IEEE Transactions on Power Systems, vol. 21, no. 4, pp. 1539?1549, Nov. 2006. [37] L. Wang, J. Jatskevich, C. Wang and P. Li, ?A voltage-behind-reactance induction machine model for the EMTP-type solution,? IEEE Transactions on Power Systems vol. 23, no. 3, pp. 1226?1238, Aug. 2008. [38] U. Karaagac, J. Mahseredjian, O. Saad, S. Dennetiere, "Synchronous Machine Modeling Precision and Efficiency in Electromagnetic Transients," IEEE Transactions on Power Delivery, vol. 26, no. 2, pp. 1072-1082, Apr. 2011. 154  [39] U. Karaagac, J. Mahseredjian, O. Saad, "An Efficient Synchronous Machine Model for Electromagnetic Transients," IEEE Transactions on Power Delivery, vol. 26, no. 4, pp. 2456-2465, Oct. 2011. [40] U. Karaagac, J. Mahseredjian, I. Kocar, O. Saad, "An Efficient Voltage-Behind-Reactance Formulation-Based Synchronous Machine Model for Electromagnetic Transients," IEEE Transactions on Power Delivery, vol. 28, no. 3, pp. 1788-1795, Jul. 2013. [41] D. S. Vilchis-Rodriguez, E. Acha, "A Synchronous Generator Internal Fault Model Based on the Voltage-Behind-Reactance Representation," IEEE Transactions on Energy Conversion, vol. 24, no. 1, pp. 184-194, Mar. 2009. [42] D. S. Vilchis-Rodriguez, E. Acha, "Nodal Reduced Induction Machine Modeling for EMTP-Type Simulations," IEEE Transactions on Power Systems, vol. 27, no. 3, pp. 1158-1169, Aug. 2012. [43] L. Wang, J. Jatskevich, "Approximate Voltage-Behind-Reactance Induction Machine Model for Efficient Interface with EMTP Network Solution," IEEE Transactions on Power Systems, vol. 25, no. 2, pp. 1016-1031, May 2010. [44] L. Wang, J. Jatskevich, "Including Magnetic Saturation in Voltage-Behind-Reactance Induction Machine Model for EMTP-Type Solution," IEEE Transactions on Power Systems , vol. 25, no. 2, pp. 975-987, May 2010. [45] L. Wang, J. Jatskevich, "Magnetically-Saturable Voltage-Behind-Reactance Synchronous Machine Model for EMTP-Type Solution," IEEE Transactions on Power Systems, vol. 26, no. 4, pp. 2355-2363, Nov. 2011. [46] F. Therrien, L. Wang, J. Jatskevich, O. Wasynczuk, "Efficient Explicit Representation of AC Machines Main Flux Saturation in State-Variable-Based Transient Simulation Packages," IEEE Transactions on Energy Conversion, vol. 28, no. 2, pp. 380-393, Jun. 2013. 155  [47] L. Wang, J. Jatskevich, "A Phase-Domain Synchronous Machine Model With Constant Equivalent Conductance Matrix for EMTP-Type Solution," IEEE Transactions on Energy Conversion, vol. 28, no. 1, pp. 191-202, Mar. 2013. [48] A. S. Abdel-Khalik, S. Ahmed, A. A. Elserougi, A. M. Massoud, "A Voltage-Behind-Reactance Model of Five-Phase Induction Machines Considering the Effect of Magnetic Saturation," IEEE Transactions on Energy Conversion, vol. 28, no. 3, pp. 576-592, Sep. 2013. [49] M. Goto, A. Isono, K. Okuda, "Transient Behavior of Synchronous Machine with Shunt-Connected Thyristor Exciter Under System Faults," IEEE Transactions on Power Apparatus and Systems, vol. PAS-90, no. 5, pp. 2218-2227, Sep. 1971. [50] T. Zouaghi and M. Poloujadoff, "Modeling of polyphase brushless exciter behavior for failing diode operation," IEEE Transactions on Energy Conversion, vol. 13, no. 3, pp. 214-220, Sep. 1998. [51] J. J. Cathey, R. K. Cavin, A. K. Ayoub, "Transient Load Model of an Induction Motor," IEEE Transactions on Power Apparatus and Systems, vol. PAS-92, no. 4, pp. 1399-1406, Jul. 1973. [52] W. Gautchi, Numerical Analysis: An Introduction, Birkhauser, Boston, Massachusetts, 1997. [53] IEEE Guide for the Application of Neutral Grounding in Electric Utility Systems, IEEE C62.92-2, Part II?Grounding of Synchronous Generator Systems, IEEE Std. CG2.92.2 ?1989 (R 2005). [54] L. F. Shampine, M. W. Reichelt, and J.A. Kierzenka, "Solving Index-1 DAEs in MATLAB and Simulink," SIAM Review, vol. 41, pp. 538-552, 1999. [55] U. M. Ascher and L. R. Petzold, Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, Philadelphia, PA: SIAM, 1998. 156  [56] A. M. El-Serafi and N. C. Kar, ?Methods for determining the intermediate-axis saturation characteristics of salient-pole synchronous machines from the measured d-axis characteristics,? IEEE Transactions on Energy Conversion, vol. 20, no. 1, pp. 88-97, Mar. 2005. [57] S. Henschel, Analysis of electromagnetic and electromechanical power system transients with dynamic phasors, Doctoral dissertation, the University of British Columbia, Vancouver, BC, 1999. [58] V. Venkatasubramanian, ?Tools for dynamic analysis of the general large power system using time varying phasors,? International Journal of Electrical Power & Energy Systems, vol. 16, no. 6, pp. 365?376, Dec. 1994. [59] P. Zhang, J. R. Marti, H. W. Dommel, "Shifted-Frequency Analysis for EMTP Simulation of Power-System Dynamics," IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 57, no. 9, pp. 2564-2574, Sep. 2010. [60] K. Strunz, R. Shintaku, F. Gao, "Frequency-Adaptive Network Modeling for Integrative Simulation of Natural and Envelope Waveforms in Power Systems and Circuits," IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 53, no. 12, pp. 2788-2803, Dec. 2006. [61] P. Zhang, J. R. Marti, H. W. Dommel, "Synchronous Machine Modeling Based on Shifted Frequency Analysis," IEEE Transactions on  Power Systems, vol. 22, no. 3, pp. 1139-1147, Aug. 2007. [62] P. Zhang, J. R. Marti, H. W. Dommel, "Induction Machine Modeling Based on Shifted Frequency Analysis," IEEE Transactions on  Power Systems, vol. 24, no. 1, pp. 157-164, Feb. 2009. [63] L. F. Shampine, M. W. Reichelt, ?The MATLAB ODE Suite,? SIAM Journal on Scientific Computing, vol. 18, no. 1, pp. 1-22, 1997.  157  Appendices:      Appendix A:    Parameters of the Induction Machine Model in Section ?2.4 50-hp, 60 Hz induction machine ?[51] 4-pole, 3-phase, 460 V, 1705 rpm, rs = 0.087 ?, Xls = 0.302 ?, Xm = 13.08 ?, rr = 0.228 ?, Xlr = 0.302 ?. Appendix B:    MATLAB Ordinary Differential Equation (ODE) Solvers ?[15], ?[63] Solver Type of Problem Accuracy Method MATLAB?s Suggestion  ode45 Non-stiff Medium One-step explicit Runge-Kutta (4,5) (the Dormand-Prince pair) The first solver to try ode23 Non-stiff  Low  One-step explicit Runge-Kutta (2,3) (the Bogacki-Shampine pair) May be more efficient than ode45 at crude tolerances and for moderately stiff systems ode113 Non-stiff Low to high Multi-step variable order Adams-Bashforth-Moulton PECE May be more efficient than ode45 at stringent tolerances and when the evaluation of the ODE file function is expensive 158  Solver Type of Problem Accuracy Method MATLAB?s Suggestion  ode15s Stiff  Low to medium  Multi-step variable order based on the numerical differentiation formulas NDFs; optionally uses the backward differentiation formulas (BDFs) Try whenever ode45 fails or is very inefficient, or whenever the system is stiff; can solve DAEs ode23s Stiff Low One-step modified Rosenbrock formula of order 2 May be more efficient than ode15s at crude tolerances; can be more effective than ode15s for some stiff problems ode23t Moderately stiff  Low  One-step Trapezoidal Use for moderately stiff problem where numerical damping is not needed; can solve DAEs ode23tb Stiff  Low Implicit, TR-BDF2 May be more efficient than ode15s at crude tolerances ode15i Fully implicit   BDFs To solve fully implicit differential equations  159  Appendix C:    Parameters of the Synchronous Machine in Subsections ?3.7.1 to ?3.7.3 125-kW, 60 Hz synchronous machine (in pu) ?[30] rs = 0.00515, Xls = 0.0800, Xmq = 1.00, Xmd = 1.77, rkq1 = 0.0610, rfd = 0.00111, rkd1 = 0.024, Xlkq1 = 0.330, Xlfd = 0.137, Xlkd1 = 0.334, H = 2.5 s, 248.0???qX  and 0921.0???dX . Appendix D:    Parameters of the Power System in Subsection ?3.7.4      555-MVA steam turbine generator (in pu) [4, sec 8.3]: 24 kV, 2-pole, 3600 r/min, 0.9 p.f., 60 Hz, rs = 0.003, Xls = 0.15, Xmq = 1.61, Xmd = 1.66, rkq1 = 0.00619, rkq2 = 0.02368, rfd = 0.0006, rkd1 = 0.0284, Xlkq1 = 0.7252, Xlkq2 = 0.125, Xlfd = 0.165, Xlkd1 = 0.1713, H = 5.6 s, 25.0???qX , and 23.0???dX . Appendix E:    Parameters of the Power System in Section ?4.6 555-MVA steam turbine generator (in pu)  [4, sec. 8.3]: The machine is the same as subsection ?3.7.4 which is given in Appendix D. Grounding reactance (in pu, on the machine base): Xg = 0.3. Unit transformer impedance (in pu, on the machine base): Xt = 0.16 and rt = 0.02. 160  Appendix F:    Parameters of the Power System in Section ?5.5 555-MVA steam turbine generator (in pu)  [4, sec. 8.3]: The machine is the same as subsection ?3.7.4 (and section ?4.6) which is given in Appendix D. Grounding resistance (in pu, on the machine base): Rg = 0.2. Th?venin impedance of the network (in pu, on the machine base): Xs = 14% and X/R = 14. Rotor-exciter Ydy transformer impedance (in pu, on the transformer base): 10 MVA, 24:0.144:0.144 kV, XPS = XPT = XST = 8% and XPS /RPS = XPT/RPT = XST /RST =13.33.  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0165681/manifest

Comment

Related Items