"Applied Science, Faculty of"@en .
"Electrical and Computer Engineering, Department of"@en .
"DSpace"@en .
"UBCV"@en .
"Zhang, Peng"@en .
"2009-04-30T15:35:21Z"@en .
"2009"@en .
"Doctor of Philosophy - PhD"@en .
"University of British Columbia"@en .
"Electromagnetic Transients Program (EMTP) simulators are being widely used in power\nsystem dynamics studies. However, their capability in real time simulation of power systems is\ncompromised due to the small time step required resulting in slow simulation speeds.\nThis thesis proposes a Shifted Frequency Analysis (SFA) theory to accelerate EMTP\nsolutions for simulation of power system operational dynamics. A main advantage of the SFA is\nthat it allows the use of large time steps in the EMTP solution environment to accurately\nsimulate dynamic frequencies within a band centered around the fundamental frequency.\nThe thesis presents a new synchronous machine model based on the SFA theory, which\nuses dynamic phasor variables rather than instantaneous time domain variables. Apart from using\ncomplex numbers, discrete-time SFA synchronous machine models have the same form as the\nstandard EMTP models. Dynamic phasors provide envelopes of the time domain waveforms and\ncan be accurately transformed back to instantaneous time values. When the frequency spectra of\nthe signals are close to the fundamental power frequency, the SFA model allows the use of large\ntime steps without sacrificing accuracy. Speedups of more than fifty times over the traditional\nEMTP synchronous machine model were obtained for a case of mechanical torque step changes.\nThis thesis also extends the SFA method to model induction machines in the EMTP. By\nanalyzing the relationship between rotor and stator physical variables, a phase-coordinate model\nwith lower number of equations is first derived. Based on this, a SFA model is proposed as a\ngeneral purpose model capable of simulating both fast transients and slow dynamics in induction\nmachines. Case study results show that the SFA model is in excess of seventy times faster than\nthe phase-coordinate EMTP model when simulating the slow dynamics.\nIn order to realize the advantage of SFA models in the context of the simulation of the\ncomplete electrical network, a dynamic-phasor-based EMTP simulation tool has been developed."@en .
"https://circle.library.ubc.ca/rest/handle/2429/7727?expand=metadata"@en .
"3788957 bytes"@en .
"application/pdf"@en .
"SHIFTED FREQUENCY ANALYSIS FOR EMTP SIMULATION OF POWER SYSTEM DYNAMICS by Peng Zhang B.Sc., Shandong University, China, 1996 M.Sc., Shandong University, China, 1999 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in The Faculty of Graduate Studies (Electrical and Computer Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) March 2009 \u00C2\u00A9 Peng Zhang, 2009 Abstract Electromagnetic Transients Program (EMTP) simulators are being widely used in power system dynamics studies. However, their capability in real time simulation of power systems is compromised due to the small time step required resulting in slow simulation speeds. This thesis proposes a Shifted Frequency Analysis (SFA) theory to accelerate EMTP solutions for simulation of power system operational dynamics. A main advantage of the SFA is that it allows the use of large time steps in the EMTP solution environment to accurately simulate dynamic frequencies within a band centered around the fundamental frequency. The thesis presents a new synchronous machine model based on the SFA theory, which uses dynamic phasor variables rather than instantaneous time domain variables. Apart from using complex numbers, discrete-time SFA synchronous machine models have the same form as the standard EMTP models. Dynamic phasors provide envelopes of the time domain waveforms and can be accurately transformed back to instantaneous time values. When the frequency spectra of the signals are close to the fundamental power frequency, the SFA model allows the use of large time steps without sacrificing accuracy. Speedups of more than fifty times over the traditional EMTP synchronous machine model were obtained for a case of mechanical torque step changes. This thesis also extends the SFA method to model induction machines in the EMTP. By analyzing the relationship between rotor and stator physical variables, a phase-coordinate model with lower number of equations is first derived. Based on this, a SFA model is proposed as a general purpose model capable of simulating both fast transients and slow dynamics in induction machines. Case study results show that the SFA model is in excess of seventy times faster than the phase-coordinate EMTP model when simulating the slow dynamics. In order to realize the advantage of SFA models in the context of the simulation of the complete electrical network, a dynamic-phasor-based EMTP simulation tool has been developed. 11 Table of Contents Abstract ii Table of Contents iii List of Tables vi List of Figures vii Acknowledgements xi Chapter 1 Introduction 1 1.1 Background 1 1.2 Motivation 3 1.3 Contributions 5 Chapter 2 Shifted Frequency Analysis 6 2.1 Analytic Signal and Hubert Transform 6 2.1 .1 Shifted Frequency Analysis 7 2.2 SFA-Based Network Component Models 10 2.2.1 Equivalent Circuits for RLC in the Shifted Frequency Domain 11 2.2.2 Options between Complex Arithmetic and Real Arithmetic 15 2.2.3 Transformer Model in the Shifted Frequency Domain 17 2.2.4 Load Models 20 2.3 Numerical Accuracy Analysis 27 Chapter 3 Synchronous Machine Modelling Based on Shifted Frequency Analysis 33 3.1 Introduction 33 3.2 Voltage-behind-Reactance Synchronous Machine Model 33 111 Table of Contents 3.3 Synchronous Machine Modelling with SFA .38 3.3.1 Synchronous Machine Model Based on SFA 38 3.3.2 Discrete Time Model 39 3.3.3 Note on the Cylindrical-Rotor Machine Model 45 3.4 Simulation Results 47 3.4 Summary 57 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis 58 4.1 Introduction 58 4.2 Equivalent-Reduction Approach to Induction Machine Modelling in EMTP. 59 4.2.1 Equivalent-Reduction (ER) to Stator Quantities 60 4.2.2 Induction Machine Modelling Based on ER Technique 60 4.2.3 Simulation Results 65 4.3 Induction Machine Modelling with SFA 71 4.3.1 Induction Machine Modelling Based on SFA 71 4.3.2 Discrete Time Model 72 4.3.3 Simulation Results 73 4.4 Summary 80 Chapter 5 EMTP Implementation 81 5.1 Introduction 81 5.2 Program Structure 82 5.3 Test Cases 86 Chapter 6 Conclusions 102 6.1 Summary of Contributions 102 6.2 Future Research 103 Bibliography 105 Appendix A Machine Parameters 111 A. 1 Synchronous Machine Parameters 111 iv Table of Contents A.2 Induction Machine Parameters.112 Appendix B Voltage behind Reactance Induction Machine Model 113 B. 1 Voltage Behind Reactance Model of Induction Machine 113 B. 2 Discrete Time VBR Model 116 B. 3 Free Acceleration Simulation of a 3-hp Induction Machine 118 B. 4 Dynamic Performance of Induction Machine during Mechanical Torque Changes 118 Appendix C Test Case Data 123 V List of Tables Table 2.1 Load Increased with the Time 23 Table 3.1 CPU Times for 4s Simulation in Case C 50 Table 3.2 Accuracy Comparison between SFA Model and VBR EMTP Model in Case C 56 Table 4.1 CPU Times for Simulations 71 Table 4.2 CPU Times for Simulations 75 Table 4.3 CPU Times for Simulations 79 vi List of Figures Figure 2.1(a) M Coupled Inductances; (b) EMTP Equivalent Circuit; (c) SFA Equivalent Circuit 9 Figure 2.2 (a) Series Connection of M-phase R and L; (b) EMTP Equivalent Circuit; (c) SFA Equivalent Circuit 13 Figure 2.3 Linear Test Case 13 Figure 2.4 Current Flowing through Branch 1-2 (At = ims) 14 Figure 2.5 Linear Time Varying Test Case 14 Figure 2.6 Voltage across the Time Varying Inductance (At 0.5 ms) 15 Figure 2.7 Single-Phase Two-Winding Transformer 18 Figure 2.8 Secondary Winding Voltage (At = 50 ms) 20 Figure 2.9 Equivalent Circuit for the Exponential Load 21 Figure 2.10 Two Node Test Case 21 Figure 2.11 Simulation Results of the Node Voltage (a) SFA Solution (At = 10 ms) (b) EMTP Solution (At = 0.5 ms) 22 Figure 2.12 Simulation Results of Voltage Collapse (a) SFA Solution (At = 10 ms) (b) EMTP Solution (At = 0.5 ms) 22 Figure 2.13 SFA Simulation Results of the Voltage Collapse Process (At = 10 ms) 23 Figure 2.14 Steady State Equivalent Circuit 24 Figure 2.15 Discrete Time Equivalent Circuit 24 Figure 2.16 Motor Test Case 26 Figure 2.17 Induction Motor Terminal Voltage (At = 5 ms) 27 Figure 2.18 Real Power Absorbed by the Induction Motor (At = 5 ms) 27 Figure 2.19 Slip of Induction Motor (At = 5 ms) 27 Figure 2.20 Equivalent Circuit Le for Backward Euler and Forward Euler (L=1; & 1 cycle, 3 cycles, 5 cycles) 30 Figure 2.21 Accuracy of integration rules (a) Magnitude (b) Phase 31 vii List of Figures Figure 3.1 Salient-Pole Synchronous Machine and Its Windings 34 Figure 3.2 SFA Equivalent Circuit of a Synchronous Generator 43 Figure 3.3 Simulation Results with the SFA Model for Field Voltage and Mechanical Torque Changes(At=7ms) 51 Figure 3.4 Simulation Results with the SFA Model for a Three Phase Fault (At 0.5 ms): A Stable Case (fault removed at 1 .4s) 52 Figure 3.5 Simulation Results with the SFA Model for a Three Phase Fault (At = 0.5 ms): An Unstable Case (fault removed at 1 .5s) 53 Figure 3.6 Simulation Results with the SFA Model for a Three Phase Fault (At = 0.5 ms): Zoomed-in View of a Portion of Results in Figure 3. 5 54 Figure 3.7 Simulation Results with the SFA Model for a Single-Phase-to-Ground Fault (At = 0.5 ms) 54 Figure 3.8 Simulation Results by SFA: Phase A Stator Current 55 Figure 3.9 Simulation Results by SFA: Zoomed-in View of Phase A Stator Current 55 Figure 3.10 Simulation Results by SFA: Field Current 55 Figure 3.11 Simulation Results by SFA: Electromagnetic Torque 56 Figure 3.12 Simulation Results by SFA: Rotor Speed 56 Figure 3.13 Time Domain Results Using the VBR Model (At = 10 ms) 57 Figure 4.1 Torque-Speed Characteristics during Free Acceleration of a 2250-hp Induction Machine (At = 500 us) 66 Figure 4.2 Dynamic Performance of a 3-hp Induction Machine during Free Acceleration (At =500ts) 67 Figure 4.3 Dynamic Performance of a 3-hp Induction Machine during Step Changes in Load Torque (At = 500 jis) 68 Figure 4.4 Simulation Results for a 3-Phase Fault at the Terminals of a 2250-hp Induction Machine (At = 100 us) 69 Figure 4.5 Simulation Results for a 3-Phase Fault with Different Time Steps 70 Figure 4.6 SFA Equivalent Circuit of an Induction Machine 73 Figure 4.7 Stator Current during Free Acceleration of a 2250-hp Induction Machine (At = 500 us) 74 viii List of Figures Figure 4.8 Dynamic Performance of a 3-hp Induction Machine during Free Acceleration (At =500.is) 75 Figure 4.9 Dynamic Performance of a 3-hp Induction Machine during Step Changes in Load Torque (RE At = 0.5 ms, SFA At =5 ms) 76 Figure 4.10 Simulation Results for a 3-Phase Fault at the Terminals of a 2250-hp Induction Machine (At = 500 ps) 77 Figure 4.11 Simulation Results by SFA 78 Figure 4.12 Time Domain Results from an EMTP Algorithm Implemented in MATLAB 79 Figure 5.1 Schematic Structure for the EMTP Dynamic Phasor Simulation Tool 83 Figure 5.2 The H Transmission Line Model 84 Figure 5.3 Contributions of the H-Circuit Transmission Line Model to the Nodal Admittance Matrix 85 Figure 5.4 The MATLAB Code for Inserting H-Circuit Transmission Line Model into G Matrix 86 Figure 5.5 Flow Chart for Dynamic-Phasor-Based EMTP Simulator 87 Figure 5.6 One-line Diagram of a Radial Test System 87 Figure 5.7 Phase A Voltage at the Load Node (At = 5 ms) 88 Figure 5.8 Zoomed-in View of Phase A Voltage at the Load Node (At = 5 ms) 88 Figure 5.9 Phase B Voltage at the Load Node (At = 5 ms) 88 Figure 5.10 Zoomed-in View of Phase B Voltage at the Load Node (At = 5 ms) 89 Figure 5.11 Phase C Voltage at the Load Node (At = 5 ms) 89 Figure 5.12 Zoomed-in View of Phase C Voltage at the Load Node (At = 5 ms) 89 Figure 5.13 One Line Diagram for the Test Feeder 90 Figure 5.14 Phase A Voltage at the Induction Machine Node (At = 1 ms) 91 Figure 5.15 Zoomed-in View of Phase A Voltage at the Induction Machine Node (At = 1 ms) ..91 Figure 5.16 Phase B Voltage at Load Node 6 (At = 1 ms) 91 Figure 5.17 Zoomed-in View of Phase B Voltage at Load Node 6 (At = 1 ms) 92 Figure 5.18 Phase C Voltage at Node 1 (At = 1 ms) 92 Figure 5.19 The First Benchmark Network for Subsynchronous Resonance Studies 93 Figure 5.20 Generator Terminal Voltage: Phase A (At = 1 ms) 94 Figure 5.21 Transformer High Side (Bus A) Voltage: Phase A (At = 1 ms) 95 ix List of Figures Figure 5.22 Voltage across Series Capacitor: Phase A (At = 1 ms) 95 Figure 5.23 Infinite Bus (Bus B) Voltage: Phase A (At = 1 ms) 96 Figure 5.24 The Second Benchmark Network for Subsynchronous Resonance Studies 98 Figure 5.25 Generator Terminal Voltages: Phase A (At = 1 ms) 99 Figure 5.26 Voltage across Series Capacitor: Phase A (At = 1 ms) 99 Figure 5.27 Transformer High Side (Bus 2) Voltage: Phase A (At = 1 ms) 100 Figure 5.28 Infinite Bus (Bus 1) Voltage: Phase A (At = 1 ms) 100 Figure B. 1 Free Acceleration Characteristics in Stationary Reference Frame 119 Figure B.2 Free Acceleration Characteristics in Rotor Reference Frame 120 Figure B.3 Free Acceleration Characteristics in Synchronous Reference Frame 121 Figure B.4 Dynamic Performance of a 3-hps Induction Machine during Step Changes in Load Torque 122 x Acknowledgements I would like to express my deep gratitude to my supervisors Dr. Jos\u00C3\u00A9 R. MartI and Dr. Hermann W. Dommel. My admiration for their great achievements and contributions was the reason why I decided to study with them in Canada. I want to thank Dr. MartI for his consistent guidance and financial support during my study at UBC. The inspiration and warm personality of Dr. MartI and Dr. Dommel have won the author\u00E2\u0080\u0099s highest respect and love. I sincerely thank my friends and colleagues at UBC and British Columbia Transmission Corporation for their encouragement and understanding. My wife Helen, my parents Yihua and Guangfu, are sources of unconditional love and support since the beginning of this project. They are entitled to more thanks than can be expressed here. My parents-in-law Guier and Jingxiang, my brother Kun, my brother-in-law Hailu, my sister-in-law Haiyang are also important members of the team that supports me permanently. To all my family, from the depths of my heart, thank you. xi Chapter 1 Introduction Chapter 1 Introduction 1.1 Background Although a power system blackout is a low-probability event, whenever it happens, its impact on the power system and the society is catastrophic. A typical example is the Aug. 14, 2003 Blackout in the USA and Canada, which caused a total loss of 62,000 MW of loads and tripping of 531 generators [1]. Fifty million people were affected by this system collapse. During the next two years, twelve major blackouts happened in Europe, Asia and North America. A recent blackout happened in Los Angeles on September 12, 2005 affecting millions in California. In today\u00E2\u0080\u0099s deregulated electricity market environment, power system planning and operations are largely driven by economic motivations, without sufficient investments in new generation and transmission equipments. This might have contributed to the higher frequency of blackouts in recent years [2]. Because of aging infrastructures operating under stressed conditions, power system stability, including transient stability, voltage stability, and inter-area oscillations have become a major concern in the North American power grid. Traditionally in power system studies a number of simplifying assumptions are made to analyze different types of stability problems with specific tools [3]-[5]. For instance, the static techniques for long-term voltage stability analysis [6] only solve the steady-state (algebraic) response of the system. Another example is quasi-steady-state (QSS) analysis [7], which neglects the fast network transients and only considers electromechanical dynamics. Because the time varying electrical quantities in QSS are represented with phasors, these tools can only capture snapshots of the system operation ignoring the dynamics between states. Problems such as voltage stability depend strongly not only on the machine-network interaction, but also on the interaction between the electrical network and the loads. Moreover, 1 Chapter 1 Introduction severe disturbances in the system often involve frequency excursions in addition to voltage excursions [8]. Especially when the power system is approaching critical collapse points, the QSS or steady-state assumptions can result in inaccurate predictions. Tools such as the EMTP [9][10] can most closely simulate the real power system dynamics by continuously tracing the evolution of the system state in arbitrary multi-phase networks with lumped or distributed parameters. Therefore EMTP-type simulators [11] -[15] are very appealing to be used for power system stability assessment. These tools, however, require small discretization steps dictated by the need to trace the instantaneous values of all waveforms. This makes the EMTP unnecessarily slow to trace phenomena around the 60 Hz fundamental frequency. As a matter of fact, nowadays one can only simulate small-scale equivalent power systems on the real-time EMTP simulators [10] [16]. Therefore, efforts should be taken to bridge the gap between the EMTP-type simulators and a unified power system analysis tool. Reference [18] aims at combining electromagnetic and electromechanical simulations [19] by considering both the network and machine differential equations without making QSS assumptions. However, no efficient solution algorithms are provided in this reference and the presented machine models are still reduced order models in phasor form. Reference [20] provides the first systematic method for a unified steady-state and transient power system analysis tool, by combining dynamic phasor concepts with EMTP simulation. Dynamic phasor models for transmission lines and linear elements are presented in this reference. However, the machine model is still based on the dq0 model, which can present numerical stability problems for large time steps and may not be the most appropriate for a general EMTP-based dynamic phasor unified theory. Important contributions to dynamic phasor modeling of electric machines have been presented in [21], [22], [23] based on the approach of the generalized averaging method [23] which, however, is also an approximate method with higher frequency components truncated. In general, it is a challenging topic to correctly, efficiently obtain the time-domain simulation results in the neighborhood of the fundamental frequency without making QSS or other simplif4ng assumptions. Hence, the intention of this thesis is to explore an effective way to combine the EMTP solution and the dynamic phasor concept, to build dynamic phasor models for power system components, and to develop a general-purpose simulation tool to obtain dynamic phasor results from EMTP solutions. 2 Chapter 1 Introduction 1.2 Motivation It is important to have some insight into the power system signals before we can further develop an efficient and accurate method for power system dynamic simulation. A conmion observation is that, when a contingency happens, the frequencies in a power system are usually close to fundamental frequency co3 (60 Hz). That is to say, in the dynamic situation, the electrical signals in power systems have a fundamental frequency component modulated by slower events. This fact is analogous to the situation in communication systems where the carrier frequency is very high and carries on its sidebands the much lower frequency (e.g., audio or video) information. In communications theory, the carrier frequency is modulated to send out the information and then demodulated to recover the information at the receiving end [24]. Typically the carrier in a communication system is a sinusoidally varying signal at some frequency which is much higher than the frequencies contained in the information or message. The process of modulation gives the signal a nonzero bandwidth that is usually much smaller than the carrier frequency. Thus the signal can be regarded as a narrowband process, and can be accurately modeled as the product of a bandlimited \u00E2\u0080\u0098message\u00E2\u0080\u0099 waveform and a sinusoidal carrier. By this analogy, one can directly draw the conclusion that usually a power system signal u(t) can be modeled by a narrowband signal which has a small bandwidth around a. In steady-state AC circuit theory, one uses complex exponentials (i.e. phasors) to represent real sinusoidal signals with the understanding that the real part of the complex exponential gives the physical time-domain quantity of the AC signal. Going beyond the steady-state analysis, the analytic signal of a real function plays the similar role for a dynamic waveform. In fact, the analytic signal z(t) of a bandpass signal u(t) can be derived as: (1.1) where UQ) is called the dynamic phasor of u(t), zQ) = uQ) + 1H[u(t)], and H[.] denotes the Hilbert transform [25][26]. The dynamic phasor UQ) is the complex envelope of u(t), and will degenerate to the phasor U if uQ) becomes a sinusoidal signal in steady state. In fact, the analytic signal of a steady-state sinusoidal function Acos(a3t+\u00C3\u00A7o) is Be\u00C2\u00B0\u00E2\u0080\u0099 where B = Here the dynamic phasor 3 Chapter 1 Introduction becomes which is the phasor ofAcos(cot+q). Similar to the phasor, the dynamic phasors can be transformed back to time domain waveforms by taking the real part of the analytic signal in (1.1) (proof is shown in Chapter 2) It can be seen in (1.1) that the spectrum of U(t) is the spectrum of the analytic signal of u(t) shifted by synchronous angular frequency \u00E2\u0080\u0094 o. This can be called the \u00E2\u0080\u009Cfrequency shifting\u00E2\u0080\u009D property. If the power network is formulated by using dynamic phasors instead of the instantaneous time quantities, it is said to be modeled in the shifted frequency domain, as opposed to the model in the time domain. Similar idea is proposed in reference [27]. The way to build power system component models and solve the network equations in the shifted frequency domain is called the Shifted Frequency Analysis (SFA). A major advantage of the Shifted Frequency Analysis is that it allows the use of large time steps in the EMTP solution environment to accurately simulate dynamic frequencies within a band centered around a fundamental frequency. The original system is transformed into a shifted frequency system where the frequencies around the power frequency (e.g. 60Hz) become frequencies around dc (0 Hz) [28][29]. Because the time step in EMTP simulations is limited by the Nyquist frequency, low frequencies in the shifted frequency system imply the feasibility of using large integration time steps. By coupling SFA models into an EMTP-type simulator, the shifted system is then numerically integrated to obtain dynamic phasor solutions, which are more easily understood by power system operators and planners than instantaneous time domain results. At the same time, the dynamic phasor results can then be transformed back to time domain waveforms using the inverse transformation. EMTP uses the trapezoidal integration rule that is A-stable, reasonably accurate and simple. The discrete-time nodal equations formed in EMTP are elegant and can preserve the sparsity of the power system structure. All these features makes EMTP a standard simulation tool in the power industry that can model the power system at the device level. The goal of this thesis is to implement the SFA method into EMTP, which is the first practical step for building a unified power system analysis tool based on the EMTP solution. 4 Chapter 1 Introduction 1.3 Contributions A two-step strategy is employed in the research of Shifted Frequency Analysis. First, the SFA-based models for power network components are developed and tested. Second, a general- purpose SFA-based EMTP tool is built to simulate a power network. The main contributions are presented in Chapter 2 to Chapter 5 and are briefly summarized as follows. I. The definition of dyhamic phasor and the theoretical basis of the SFA method are introduced. The basic procedures for analyzing time functions (signals) in the SFA domain are discussed, and the numerical accuracy of the discrete-time SFA solution is analyzed. Chapter 2 summarizes these theoretical works. II. A few types of power system components are modeled in the SFA domain. Chapter 2 describes the SFA-based models for the linear circuit components, transformer and some load models such as exponential load and steady state induction motor load. III. A new efficient synchronous machine model for the simulation of slow system dynamics using EMTP modelling in the Shifted Frequency Analysis framework is presented in Chapter 3. Discrete-time SFA synchronous machine models are derived which have a similar form as the EMTP component models. IV. In Chapter 4, the SFA method is extended to model the induction machines in EMTP. By analyzing the relationship between rotor and stator physical variables, a phase-coordinate model with lower number of equations is first derived. Based on this, a SFA model is proposed as a general purpose model capable of simulating both fast transients and slow dynamics. Case study results have confirmed that the SFA induction machine model is a valuable component for real time EMTP simulations. It is observed that the SFA model is in excess of 70 times faster than the EMTP phase-coordinate model when simulating dynamics with frequency spectra close to the fundamental power frequency. V. A final contribution in this thesis is a dynamic-phasor-based EMTP simulation tool. The structure of the tool and the test system results are reported in Chapter 5. 5 Chapter 2 Shifted Frequency Analysis Chapter 2 Shifted Frequency Analysis 2.1 Analytic Signal and Hubert Transform The Fourier transform F(co) of a real-valued functionf(t) is a Hermitian function, which means F(\u00E2\u0080\u0094co) = F*(w) in the frequency domain. Thus the F for negative frequencies can be expressed by F* (complex conjugate of F) for positive ones. This means that the positive frequency spectra is adequate to represent a real signal, and the negative frequency components of the Fourier transform can be discarded without loss of information. Now one can construct a function Fa(0) that contains only the non-negative frequency components of F(co) by defining Fa(O)) F((D) + sgn(aF(\u00C3\u00BC) (2.1) where 1 a)>O sgn(o)= 0 o=0 \u00E2\u0080\u00941 o<0 The analytic signalfa(t) is defined to be the inverse Fourier transform Of Fa(CO) L(t)=f(o0)td = -- f F(o)e\u00C2\u00B0td (2.2) It is clearly seen thatLQ) is a complex function. If we define the Hubert transform off(t) as the imaginary part Offa(t), and denote it by H[fQ)] = f(t), thenL(t) can be expressed as fQ) =JQ) +jfQ) (2.3) 6 Chapter 2 Shifted Frequency Analysis From equation (2.1) and (2.3) we can get the following Fourier pair JQ) +jfQ) -> F(o) + sgn(F(w) F Because we haveJ(t) F(w) by definition, we can derive that f(t) -jsgn(co)F(co) iF Using the transform pair \u00E2\u0080\u0094 -+ \u00E2\u0080\u0094jsgn(co) and taking the inverse Fourier transform of \u00E2\u0080\u0094jsgn(o).F(a), we can get H[fQ)]=j(t) = \u00C2\u00B1*fQ)= --PVfdr (2.4) Note that the Hubert transform is an improper integral, thus it is calculated using the Cauchy principal value (see the \u00E2\u0080\u0098PV\u00E2\u0080\u0099 in (2.4)), i.e. f(t) iimi[j dr + A dr].-A t\u00E2\u0080\u0094r Obviously, H[f(t)] is real. The analytic signal has no negative frequency components; moreover, the original real signal can be converted back from it by simply dropping the imaginary part. The analytic signal is thus a generalization of the phasor concept. This important implication leads to the shifted frequency analysis method as is explained in the following section. 2.1.1 Shifted Frequency Analysis If a power system signal u(t) is a bandpass one, it can then be represented [201 as. u(t) = u1 (t) cos oit \u00E2\u0080\u0094 UQ (t) sin ot (2.5) where the lowpass signals u1 (t) and UQ (t) are the in-phase and quadrature components of the bandpass signal, respectively. The dynamic phasor U(t) for the signal u(t) is defined as U(t)=uJ(t)+juQ(t) which is the complex envelope of u(t). Assume that z(t) is the analytic signal ofu(t), then z(t)e\u00C2\u00B0\u00E2\u0080\u0099 = {uQ) + jH[u(t)Je0t 7 Chapter 2 Shifted Frequency Analysis = {u1 (t) cos ot \u00E2\u0080\u0094 UQ (t) sin ot + j[1 (t) sin cost + UQ (t) cos u1(t)+fu0t) (2.6) where H[.] denotes the Hubert transform. Therefore, U(t) z(t)e\u00C2\u00B0\u00E2\u0080\u0099 (2.7) Note that the dynamic phasor U(t) will degenerate to the phasor U if u(t) becomes a sinusoidal signal in steady state. Equation (2.7) shows that the spectrum of U(t) is the spectrum of the analytic signal of u(t) shifted by synchronous angular frequency \u00E2\u0080\u0094 co,. This is called \u00E2\u0080\u009Cfrequency shifting\u00E2\u0080\u009D. Shifted frequency analysis (SFA) [24] allows the exact simulation of frequencies within a band centered around a fundamental frequency using large time steps in a discrete-time EMTP type of solution environment. The original system is transformed into a shifted frequency system where the frequencies around the power frequency (60Hz or 50Hz) become frequencies around dc (0 Hz). The shifted system is then solved using an EMTP solution. The equivalent circuit for network components in the shifted frequency domain can be derived in three steps: (i) Create the phase-coordinate differential equations of the component in the normal unshifted domain; (ii) Transform phase quantities into dynamic phasor variables according to Equation (2.7); (iii) Discretize the dynamic phasor equations using an integration method, e.g. the trapezoidal rule, and build the equivalent circuit suitable for an EMTP solution. For instance, the time domain voltage equations for M coupled inductances, shown in Figure 2.1(a), can be written as v(t)=tLt) (2.8) As illustrated in Figure 2.1(b), the standard EMTP equivalent circuit derived with the trapezoidal rule can be written as 1(t) =R1v(t) \u00E2\u0080\u0094 hL (t) (2.9) where 8 Chapter 2 Shifted Frequency Analysis RL =- (2.10) is an MxM matrix of resistances, and hL(t)=L\u00E2\u0080\u0099v(t \u00E2\u0080\u0094 At)At + hL (t \u00E2\u0080\u0094 At) (2.11) is a vector of past histories. To obtain the SFA equivalent circuit (Figure 2.1(c)), we transform (2.8) using (2.7) to get V(t)= L dIQ) +jo3L1(t) (2.12) where V(t) and 1(t) are the dynamic phasor vectors corresponding to the physical time vectors v(t) and i(t), respectively. We can now transform (2.12) to discrete time using, for example, the trapezoidal rule, 1(t) =R1V(t) \u00E2\u0080\u0094 HL (t) (2.13) where the equivalent resistances RL and history terms HL (t) are expressed as 2L (2 RL =\u00E2\u0080\u0094+KL (2.14) 4 2. \u00E2\u0080\u0094-\u00E2\u0080\u0094+JO) At 2 L1V(t \u00E2\u0080\u0094 At) \u00E2\u0080\u0094 HL (t \u00E2\u0080\u0094 At) (2.15) -+jcO The corresponding SFA equivalent circuit is shown in Figure 2.1(c). L p p (a RL = 2L/At RL = (2/At+ jco)L hL(t) HLQ) (b) (c) Figure 2.1 (a) M Coupled Inductances; (b) EMTP Equivalent Circuit; (c) SFA Equivalent Circuit A major advantage of SFA modelling is that since the equivalent circuits are centered 9 Chapter 2 Shifted Frequency Analysis around 60 Hz (ja)L for the inductances of Figure 2.1), deviations from 60 Hz correspond numerically to very low frequencies (like deviations from 0 Hz in the unshifted domain) and a large integration step can be used in the solution. For example, respecting the accuracy limits imposed by the Nyquist frequency and the distortion introduced by the discretization rule [30], allowing a 3% error with the trapezoidal rule, the integration step can go, for example, from about three 60 Hz cycles, 5Oms, for frequency deviations of\u00C2\u00B1 2Hz to as large as is for frequency deviations of\u00C2\u00B1 0.1Hz. With these time steps, the SFA method provides an effective way to make the EMTP program a general purpose simulator for power system stability problems. By applying SFA in an EMTP-type simulator, time varying phasor solutions are obtained, which can be easily understood by power system operators and planners. At the same time, detailed waveform results can also be traced back from the SFA results by using the inverse transformation u(t) = Re[U(t)e.t 1 (2.16) Unless specifically noted, in this paper uppercase letters are used to represent dynamic phasors in the shifted frequency domain, while lowercase letters denote real variables in the time domain. 2.2 SFA-Based Network Component Models The component models in the SFA domain are building blocks for a SFA-based network simulator. In this section, we implement the SFA modeling technique described in Section 2.1.1 and build component models for the linear RLC components, transformer, exponential load and steady-state induction machine. SFA equivalent circuits for other network components can also be derived following the procedures in Section 2.1.1. Particularly, two important dynamic components in the power system, i.e. synchronous machine and induction machine, are modeled in SFA domain in Chapter 3 and Chapter 4. 10 Chapter 2 Shifted Frequency Analysis 2.2.1 Equivalent Circuits for RLC in the Shifted Frequency Domain The equivalent circuits for network components in the shifted frequency domain can be derived by first writing the component equations in the phase domain and then relating phase and dynamic phasor variables according to the frequency shifting transformation. A. M-Phase Resistances The time domain voltage equations for M-Phase resistances are expressed by vQ)=Ri(t) (2.17) By using (2.7), we can obtain the SFA form of equation (2.17) VQ)=RI(t) (2.18) Series resistance matrix sometimes appears as a part of the it-circuit representation of the transmission line. if an M-phase transmission line is modeled as an M-phase it-circuit, the series resistance matrix will be a full MxM matrix. The off-diagonal elements come about because the earth return is eliminated as the (M+1)th conductor [9]. B. M-Phase Coupled Inductances Equations (2.12)-(2.15) give us the shifted frequency domain equations and the discrete- time equivalent circuit for L. C. M-Phase Coupled Capacitances The time domain voltage equations for the M-Phase capacitances are written as follows (2.19) The SFA form of equation (2.8) is obtained by employing the shifted frequency transformation, I (t) = C dV(t) +K V (t) (2.20) where K jwC. By discretizing (2.20) by the trapezoidal integration rule, the difference equations are obtained I(t)GV(t)h(t) (2.21) where 11 Chapter 2 Shifted Frequency Analysis Rc=(Gc)\u00E2\u0080\u0099 = hc(t)= Equation (2.21) is the EMTP equivalent circuit for C in the shifted frequency domain. D. M-Phase Coupled R-L Branches The time domain voltage equations for a series connection of M-Phase resistances and coupled inductances, shown in Figure 2.2(a), can be written as v(t) R 1(t) + L diQ) (2.22) The standard EMTP equivalent circuit (Figure 2.2(b)) for the series connection of R and L are reproduced as follows i(t) =rv(t) \u00E2\u0080\u0094h(t) where rL =R+--L=g is a matrix of resistances, and _I]v(t_At) is a vector of history terms. We now transform (2.22) using (2.7) to obtain the dynamic phasor equations V(t)= R I (t) + L dI(t) +KL 1(t) (2.23) where KL =jco3L, V(t) and 1(t) are the dynamic phasor vectors corresponding to the time-domain real variables v(t) and 1(t), respectively. Equation (2.23) can be transformed to discrete time by using the trapezoidal rule, 1(t) R;V(t) \u00E2\u0080\u0094 H (t) (2.24) where the equivalent resistances R and past histories H (t) are expressed as R= 12 Chapter 2 Shifted Frequency Analysis _I1V(t_At)_G(R_\u00E2\u0080\u0094-L +KLJ H(t-At) The SFA equivalent circuit described in (2.24) is illustrated in Figure 2.2(c). R L r=R+(2/At)L (a) R=R+ [(2/At) +jcoJL H(t) (c) Figure 2.2 (a) Series Connection of M-phase R and L; (b) EMTP Equivalent Circuit; (c) SFA Equivalent Circuit E. Case Study Two test cases are simulated to explore the feasibility of the SFA modeling method. The first test case is to simulate the switching operations in a typical linear time invariant (LTI) circuit that is shown in Figure 2.3. A voltage source is switched into the circuit at t = 0 s. After 8 ms, the switch Si opens once the current flowing through the switch crosses zero. In this circuit, R1 = 3 2, R2=50 2, L1 = 300 mH, L2= 1000 mH, C1 = 20 jiF, C2 = 6 1iF. The voltage source has a rms value of 230 kV and a frequency of 60 Hz. Figure 2.3 Linear Test Case The current following through the branch connecting node 1 and node 2 is illustrated in Figure 2.4. From the simulation result we can find (b) 1 Ri Li 2 3 L2R2 13 Chapter 2 Shifted Frequency Analysis The second test case is a linear time varying (LTV) circuit. A sinusoidal current source 1(t) = A cos(c#) is applied to a time varying inductance L cos(o0t), as is shown in Figure 2.5. The frequency of the current source is co = 60Hz, and co =10Hz. The current source has a magnitude of 1 A, and we choose an integration step At = 500ps. The unknown variable v(t) is the voltage across the inductance. Figure 2.5 Linear Time Varying Test Case The numerical solution to the dynamic phasor V(t) can be readily obtained by using the trapezoidal rule, which is given by (1) During the transient state, the SFA solutions are the envelop of the time domain solutions; (2) During the steady state, the SFA solutions degenerate to the phasor, which is the concept being widely used in power industry; (3) Time domain simulation results can be accurately traced back from the SFA solutions. We can see from Figure 2.4 that the time domain curve is identical to that obtained from the EMTP time-domain simulation. Shifted Frequency Analysis: Test Case -I- 9 Time (s) Figure 2.4 Current Flowing through Branch 1-2 (At = ims) 1(t) = Acos(o0 v(t) 14 Chapter 2 Shifted Frequency Analysis V(t) = \u00E2\u0080\u0094V(t \u00E2\u0080\u0094 At) + [f0JL cos(co0t)+ o)0Lsin(a0t)+ ..cos(aot)]I(t) + [JoL cos(o.0(t \u00E2\u0080\u0094 At)) + L sin(a0(t \u00E2\u0080\u0094 At)) \u00E2\u0080\u0094 cos(a0(t \u00E2\u0080\u0094 At))]I(t \u00E2\u0080\u0094 At) The time domain solution of voltage across the inductance can be transformed back from V(t) by v(t) = Re[V(t)e3]. Figure 2.6 Voltage across the Time Varying Inductance (At = 0.5 ms) The simulation result is shown in Figure 2.6. It is clearly seen that the dynamic phasor solution is the envelope of the time domain solution. This test case validates that the SFA method is suitable for the LTV circuit analysis. 2.2.2 Options between Complex Arithmetic and Real Arithmetic As can be seen in Section 2.2.1, the discretized power systems equations based on SFA models will be a system of linear complex equations, which can be solved by using either the complex arithmetic or the real arithmetic. Some references [311 [321 concluded that, for the complex matrices, the complex inversion of them might be up to twice as fast as the real inversion, and that the rounding error bound for complex inversion is tighter than that for real inversion in terms of Gauss elimination. Fortunately, this may no longer be a notable issue with the modem computer architectures. Therefore, the real arithmetic method can be used as an alternative to the complex-arithmetic-based system solver. This section gives some examples about building real valued equivalent circuit for the RLC elements. time (s) 15 Chapter 2 Shifted Frequency Analysis A. M-Phase Resistances The discretized equations (2.18) can be decomposed into two equation sets corresponding to the real and imaginary parts of the dynamic phasors for voltages and currents. This leads to the real valued equivalent circuit of R. fIre(t)1 = rR\u00E2\u0080\u0099 lEVre(t)1 (2.25) [\u00E2\u0080\u0098im (t)J L R\u00E2\u0080\u0099 ]LYm (t)] B. M-Phase Coupled Inductances The decomposed difference equations for the M-Phase coupled inductances can be derived as Ire (t)1 = raL1 \u00E2\u0080\u0094 bL\u00E2\u0080\u0099 lrvre (t)1 \u00E2\u0080\u0094 rhLre (t)1 (2.26)[\u00E2\u0080\u0098im (t)] LbL1 aL1 ]L\u00E2\u0080\u009Dim(t)J L1\u00E2\u0080\u0099L,Im (t)j where P\u00E2\u0080\u0099 (t)1 \u00E2\u0080\u0094 rcL-\u00E2\u0080\u0099 \u00E2\u0080\u0094 \u00E2\u0080\u0098-\u00E2\u0080\u0098 (t \u00E2\u0080\u0094 At)1 + r e \u00E2\u0080\u0094 fThL,re (t \u00E2\u0080\u0094 At) [hLIrn (t)] \u00E2\u0080\u0094 LdL1 cL1 ][V. (t \u00E2\u0080\u0094 At)] Lf e ][hLjm (t \u00E2\u0080\u0094 At) 2 a=y bCzD2 4E(22 21 (42 -\u00E2\u0080\u0094-il\u00E2\u0080\u0094-i-a I = Atk At) d = Lzxt 2 r(2 21 E(2 2II\u00E2\u0080\u0094I+(L) I Il\u00E2\u0080\u0094I +0)LAt) SJ LAt) S (22 2 4 -I\u00E2\u0080\u0094-I+a) 0)\u00E2\u0080\u0094 \u00E2\u0080\u0094 VAt) S e\u00E2\u0080\u0094 (22 2 (2 2I\u00E2\u0080\u0094I+w \u00E2\u0080\u00941+0) VAt) VAt) C. M-Phase Coupled Capacitances We can get the real valued equivalent circuit for C as the following 16 Chapter 2 Shifted Frequency Analysis rlre(t)1 = C O)C rvre(t)1_rhc,re(t)1 (2.27) LIirn(t)J c --c LVirn(t)J Lhcim(t)J At where 4 r hCre(t)1 \u00E2\u0080\u0094 Ac 0 rvre(t \u00E2\u0080\u0094 At)1 \u00E2\u0080\u0094 Ehc,re(t \u00E2\u0080\u0094 At) Lh17 0\u00E2\u0080\u0099)] \u00E2\u0080\u0094 c L\u00E2\u0080\u0099m(t \u00E2\u0080\u0094 At)] [hcim (t \u00E2\u0080\u0094 At) At D. M-Phase Coupled R-L Branch By decomposing the equivalent circuit in (2.24), we can obtain the real valued equivalent for a series RL branch rlre(t)1 \u00E2\u0080\u0094 rG,re _Gpj,jmjVre(t)][hp,re(t) 228LI1,(t)] \u00E2\u0080\u0094 [Gjm Gre iLVim 0\u00E2\u0080\u0099)] [h1,0\u00E2\u0080\u0099) where Eh,re (t)1 \u00E2\u0080\u0094 rc \u00E2\u0080\u0094 Dlrvre 0\u00E2\u0080\u0099 \u00E2\u0080\u0094 At) \u00E2\u0080\u0094 \u00E2\u0080\u0094 FlEhc reQ \u00E2\u0080\u0094 At) 0\u00E2\u0080\u0099)] \u00E2\u0080\u0094 LD c ]LV,rn 0\u00E2\u0080\u0099 \u00E2\u0080\u0094 At)] [F E ][hcim (i\u00E2\u0080\u0099 \u00E2\u0080\u0094 At) A=R---L At B=o3L C = \u00E2\u0080\u0094 G RL,re + G RL re AG jLre \u00E2\u0080\u0094 GpjBG \u00E2\u0080\u0094 G jL,reBGi?J,,im \u00E2\u0080\u0094 G jm D = \u00E2\u0080\u0094G jL,im + G BG + GpjAG + G RL,re AG jm \u00E2\u0080\u0094 G jm BGpjm E G RL,re A \u00E2\u0080\u0094 GpjjB F 2.2.3 Transformer Model in the Shifted Frequency Domain Since the SFA modeling focuses on the low frequency dynamics as opposed tO the standard EMTP modeling for simulating the fast transients, here only the low frequency transformer model is built in the SFA domain. For transients with highest frequency less than 2-3 17 Chapter 2 Shifted Frequency Analysis kllz, the transformer can be modeled as a series connection of multi-phase coupled R and L branches, in which each RL branch can represent one transformer winding [9]. Therefore, the SFA domain model of a transformer follows directly the model described in equation (2.24). Here a major issue is how to obtain the G1 (thus the history terms) in (2.24). A traditional way for single-phase transformers is the [Z] matrix method, which first builds the coupled impedances from the transformer parameters. For example, the [Z] matrix of a two winding transformer illustrated in Figure 2.7 can be obtained as below [z]= Then, [Zj is inverted to get the admittance matrix [Y] that will be used to calculate G1. This method, however, has limitations such as: (i) Magnetizing impedance Zm cannot be neglected or set to be co in this model. Otherwise the model cannot work. (ii) Zm usually dominates [Z] matrix because Zm >> leakage impedances Zi and Z2. However, it is the leakage impedances that largely determine the simulation results. This means all data inputs are required to be very precise. Moreover, [Zf1 is ill-conditioned because of the dominating Zm, which may negatively affect the accuracy in the simulation results. [Zj Z1 N1:N2 Z2 I\u00E2\u0080\u0099 VI Zm j I Figure 2.7 Single-Phase Two-Winding Transformer In this section the [Lf1 model [9] adopted in MicroTran, rather than the [Zj matrix model, is used. This model can directly formulate the [Lf\u00E2\u0080\u0099 matrix without doing inversion on the [Z] matrix. It works for any number of windings, and for single-phase as well as three-phase transformers. For instance, the [Lf1 matrix of a two winding transformer in Figure 2.7 can be built by 18 Chapter 2 Shifted Frequency Analysis 1 N1 [L]1= Ni (N2i (2.29) N2L N2) L where L=Li+JL2 An alternative to the [Lf\u00E2\u0080\u0099 model is to use an inductance matrix in series with an ideal transformer [14]. The theory behind the [Lf\u00E2\u0080\u0099 model and the [Lf\u00E2\u0080\u0099 matrices for other types of transformers can be found in [91. The equivalent conductance and history terms matrices can be calculated by using (2.24). Note that some modifications in (2.24) are needed for the usage of [Lf\u00E2\u0080\u0099, as expressed below G =R+--L + KU =[LR+(_+J)I]1L-\u00E2\u0080\u0099 HRL(t)={[L\u00E2\u0080\u0099R + + frD JI] [LR + + JDsJI1G \u00E2\u0080\u0094 G}V(t \u00E2\u0080\u0094 At)\u00E2\u0080\u0094 [L-\u00E2\u0080\u0099R + +. [L-\u00E2\u0080\u0099R + + joi] H(t\u00E2\u0080\u0094At) The magnetizing branch is not required in this model. It can be added at the terminals when it is needed. This will allow one to model the saturation effect of the core by adding the nonlinear Lm at the terminals. As a test case, a two winding transformer with the on-load tap changer (OLTC) is simulated using the SFA model. The rated line-to-line voltage of the primary winding is 110kV. The winding ratio is N1: N2 = 110 : 28.4. R = 10 , X 33.2174 2. The transformer is first connected to the rated voltages and operating until the primary winding voltage suddenly drops to 0.95 p.u. at t = lOs. Five seconds later, the secondary side OLTC changes the tap ratios to 1+ 0.0125, attempting to bring the voltage back. The total simulation time is 20s. Figure 2.8 illustrates the dynamic phasor result together with the time domain result transformed back from SFA solution. In this simulation, a very large time step At = 50 ms is used. On the other hand, if 19 Chapter 2 Shifted Frequency Analysis we simulate this case using the EMTP algorithm, the time step for an accurate solution should at least satisfr At < 11(5 x 60) =3.3 ms to respect the Nyquist frequency limit. Figure 2.8 Secondary Winding Voltage (At = 50 ms) 2.2.4 Load Models A. Exponential Load Model An exponential load means the power consumed by the load depends exponentially on the load voltage, which can be expressed as ( \flp PQ) = PoLJ (2.30) 0 ( \,nq QQ) = (2.31) 0 where Po, Qo, and Vo are the pre-disturbance conditions of the load. For special cases such that load parameters np or nq becomes 0, 1, or 2, the load model will represent the constant power, constant current, or constant impedance load, respectively. The following equations express the discrete-time equivalent circuit (see Figure 2.9) for the exponential load. They can be directly derived from (2.30) and (2.31) by using the backward Euler rule. The backward Euler rule avoids predictions in the simulation with load models, though it may introduce some damping into the equivalent circuit. The trapezoidal rule can also Transformer Test Case 02 \u00E2\u0080\u0094 SFA Dynamic Phasor Solution zsrV_ - V 0 2 4 6 8 10 12 14. 16 18 20 Time (s) 20 Chapter 2 Shifted Frequency Analysis IV(t - At)12 RlOOd (t) = P(t) V(t - At)12 kOd(t) = oQ(t) RL1d(t) = + JO)0 L10 (t) \u00E2\u0080\u00944 2 hLld(t) = At L1 (t)V(t \u00E2\u0080\u0094 At) + At \u00E2\u0080\u0094 2 hL,d (t \u00E2\u0080\u0094 At)(2 At load I + Joi At P(t) Iv(t\u00E2\u0080\u0094At)= (VQ\u00E2\u0080\u0094At) jiq vc P+jQ \u00E2\u0080\u0094b Roa\u00E2\u0080\u0099 Lload \u00E2\u0080\u0094, RlaRLl0 hLload Figure 2.9 Equivalent Circuit for the Exponential Load In this section, a two node test case shown in Figure 2.10 is simulated to demonstrate the load modeling by the SFA method. This test case for voltage stability studies is taken from [33]. Figure 2.10 Two Node Test Case be used but iterations have to be applied at each time step due to the nonlinearity of the exponential load model. (2.33) (2.34) (2.32) (2.35) (2.36) (2.37) where I V(t)I is the magnitude of dynamic phasor V(t). A B ZjO.25 P+jQ 21 Chapter 2 Shifted Frequency Analysis First, the system is operating with the load P= O.8p.u., and Q=O.4p.u. The load voltage at node B is illustrated in Figure 2.11. It can be seen that the SFA curve is the complex envelope of the time domain curve, to which the SFA result can be accurately transformed. Two Node Test Csse Two Node Test Case 92 90 88 > 86 84 -\u00E2\u0080\u00940 1 2 3 4 t (a) (a) (b) Figure 2.11 Simulation Results of the Node Voltage (a) SFA Solution (At = 10 ms) (b) EMTP Solution (At = 0.5 ms) Second, the system is operating with a heavier load of P= 1 .85p.u., Q=O.1 5p.u. From the nodal voltage shown in Figure 2.12, it can be seen that the system collapses in this case. These results in Figure 2.12 are almost the same as the results obtained in [33], where the system is found to collapse when P=1.843p.u and Q=O.15p.u. Two Node Test Case Figure 2.12 Simulation Results of Voltage Collapse (a) SFA Solution (At 10 ms) (b) EMTP Solution (At = 0.5 Now the system performance is simulated when the system load keeps increasing with time. The system load is assumed to change in the way shown in Table 2.1. The SFA results of 22 - 0 0.2 0.4 0.6 t (s) 0.8 Two Node Test Case t (s) ms) Chapter 2 Shifted Frequency Analysis the load node voltage are illustrated in Figure 2.13. It can be seen that the SFA model can accurately trace the voltage stability behavior of the system when the exponential dynamics of the load are taken into consideration. Table 2.1 Load Increased with the Time Time interval (s) [0, 2] [3, 4) [4, 5] [5, 6] [6, 7] [7, 8] ActiveLoad(p.u.) 0.7 0.8 0.8. 0.9 1.0 1.85 Reactive Load (p.u.) 0.3 0.4 0.5 0.5 0.6 0.15 Two Node Test Case 12C 110 100 ______ 90 180 70 , 60 50 40 t (s) Figure 2.13 SFA Simulation Results of the Voltage Collapse Process (At = 10 ms) B. Steady State Induction Machine Load The steady state induction motor load is modeled in this section. The machine may be first initialized at the beginning of the simulation. Usually, the induction machine parameters such as R, X, XM, XR, RR are known quantities as input data. In addition, the real power P and the terminal voltage are also known, where V can be determined by the steady state power flow calculations. Hence the unknowns for initialization are the slip s, and reactive power Q. With the core loss neglected, the equivalent circuit of the induction motor is shown in Figure 2.14. 23 Chapter 2 Shifted Frequency Analysis P + jQ, X -- XM :RIs Figure 2.14 Steady State Equivalent Circuit The formula to calculate the initial slip s and the reactive power Q can be derived [34] based on the equivalent circuit in Figure 2.14, as expressed below. RR (2.38) _g2 -4AC Q=Im(V2Y) (2.39) where A=A\u00E2\u0080\u0099---R V2 B=B\u00E2\u0080\u0099---X C=C\u00E2\u0080\u0099---RS(XS+XM)V2 A\u00E2\u0080\u0099 = + (x + XM )2 B\u00E2\u0080\u0099 = 2RX C\u00E2\u0080\u0099 = (XSXM + XSXR + XRXM )2 + R:(xR + XM)2 Then the initial slip can be validated by checking whether the power factor falls into a reasonable range, for instance, 0.7 ___________ O 1855 1805 1755 1702 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3 Ume (s) Figure 4.3 Dynamic Performance of a 3-hp Induction Machine during Step Changes in Load Torque (At = 500 is) C. Three Phase Fault Test In this simulation of fast transients, the 2250-hp induction machine is first operating under rated conditions with a load torque equal to TB. Here TB is the base torque defined by TB = where PB is the rated power output of the machine and 0b is the rated speed of the (2/p)cob machine. A 3-phase fault is then applied at the terminals at t = 6.1 s. After 6 cycles, the fault is cleared. The simulation time is 10 s with a time step of 100 ,us. Figure 4.4 illustrates that the proposed ER model can simulate fast dynamics of the machine as accurately as the VBR model. Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis \u00E2\u0080\u009400 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 300 2000 1005 -1005 -2005 11* ER Result -- VBRResuH I I 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3 68 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis E Figure 4.4 Simulation Results for a 3-Phase Fault at the Terminals of a 2250-hp Induction Machine (At = 100 its) Cases A to C show that the ER model is as accurate as the VBR model for the same At step. It is also interesting to explore the numerical accuracy and robustness of the ER model by using different At steps in the fast transients simulation of Case C. Results obtained from the VBR simulation [53] with a very small At = 10 ps are assumed to be the reference accurate time domain results. A zoomed-in view of phase A stator current is shown in Figure 4.5. Figure 4.5 indicates that, if we run an EMTP simulation with the ER model using a time step of 1 ms, relatively accurate results are still obtained. Even with a time step of 1.5 ms, the numerical solution remains stable. 5000 -5000 ER Reso1 - - VBR Res,1t 6.1 6.2 6.3 6.4 6.5 6.6 6.7 6.8 6.9 69 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis Table 4.1 summarizes the CPU times needed for the VBR and ER simulations in cases A to C. The timing in Table 4.1 was obtained on a 1.83G Hz PC. It can be calculated that the ER is usually 33%-36% faster than the VBR model for the same accuracy. The efficiency of the ER model is due to its concise form and the low computer burden in updating electrical variables and flux linkages. The ER model has a further advantage in addition to the faster speed. Instead of being specifically designed for the machines to achieve efficient phase domain models, the ER technique used in the induction machine modeling is a general-purpose approach. ER can also be used for building more efficient models for other types of machines, transformers and power system elements, as long as they can be mathematically represented by differential equations with mutually coupled (or decoupled) state variables. Since the ER-based model proposed in this chapter is actually a full-order model, it is different from the Model-Order-Reduction-based models. Model Order Reduction (MOR) approaches have been proposed, either for constructing a reduced-order model of the machine itself by neglecting stator flux derivatives, or for constructing a reduced-order model of the external system by using techniques such as Krylov-subspace or Fourier methods [54]-[57]. MOR methods improve the simulation speed; however, undesirable inaccurate step-by-step simulation results will be produced, especially from the latter one. Figure 4.5 Simulation Results for a 3-Phase Fault with Different Time Steps 70 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis Table 4.1 CPU Times for Simulations Test Case VBR Model ER Model Case A (3-hp induction machine, * 0.875 s 0.556 sT= is, IXt=500,us) Case A (2250-hp induction machine, 2 601 . 1 698T=3s,At=5004u ) Case B (3-hp induction machine, 2 606 . 696T=3s,At=SOOus) Case C (2250-hp induction machine, 8 391 . 5 562T=iOs,At=500ps) Case C (2250-hp induction machine, 82.796 . 55.031 .T lOs,At504us) * T stands for the total simulation time. 4.3 Induction Machine Modelling with SFA 4.3.1 Induction Machine Modelling Based on SFA To build the dynamic-phasor-based stator equations, the stator variables are first expressed as \u00E2\u0080\u0098abcs (t) = abcsl (t) C05 \u00E2\u0080\u0094 \u00E2\u0080\u0098abcsQ (t) sin w,t (4.12) V abcs (t) = v abcsl (t) cos o t \u00E2\u0080\u0094 v abcso (t) sin OJJ t (4.13) Substituting (4.12) and (4.13) into (4.6), we get 1abcsI (t) cos (J)5t \u00E2\u0080\u0094 VbQ (t) sin co5t = r5 [Iabcd (t) cos cot \u00E2\u0080\u0094 I abcsQ (t) sin th5]+ \u00E2\u0080\u0098\u00E2\u0080\u0098abcs [pi (t) cos cost \u00E2\u0080\u0094 abcs! (t)o5 sin oi5t \u00E2\u0080\u0094 P\u00E2\u0080\u0099abcso (t) sin ot \u00E2\u0080\u0094 abcsQ (t)co5 cos oJ5t] (4.14) Then we apply the Hilbert transform to both sides of (4.14) and construct the analytic signal, which leads to Vabcs (t)e\u00C2\u00B0\u00E2\u0080\u0099 = r5 \u00E2\u0080\u0098abcs (t)e\u00C2\u00B0\u00E2\u0080\u0099 + Labcs [Plabcs (t) + fcOSI abcs (t)] e\u00C2\u00B0\u00E2\u0080\u0099 + V (t)e\u00C2\u00B0 (4.15) where VabL,s (t) and \u00E2\u0080\u0098abcs (t) are dynamic phasors of stator voltages Vabcs (t) and currents \u00E2\u0080\u0098abcs (t), respectively, and Vre, (t) is the dynamic phasor of v (t). 71 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis Applying frequency shifting as described in (2.7), we can obtain the dynamic-phasor equations for the stator part as Vb(t) r5 Iabcs(t)+Labcs[PIab( )+ J\u00C3\u00BCSIabcS(t)] + Vres(t) (4.16) 4.3.2 Discrete Time Model Discretizing (4.16) with the trapezoidal rule of integration gives us Vb(t) = [r +--+foJLb ] Is(t) +V(t)+Ehl(t) (4.17) where Ehl (t) [i. + + js Jtabcs] \u00E2\u0080\u0098abcs (t \u00E2\u0080\u0094 At) \u00E2\u0080\u0094 Vb (t \u00E2\u0080\u0094 At) + \u00E2\u0080\u0098res (t \u00E2\u0080\u0094 At) From the discrete-time equation of Vres (t) (see Section 4.2), the corresponding discrete- time equation of r8S(r) can be derived as Vres = k(t)Iabcs (t) + E$h (t) e1(00t) e1(0__2) 1 \u00E2\u0080\u0094 k(t)Iabcs(t) + eo23) ej(ot232) 1 . [k4(t)] (4 18) e j(9\u00E2\u0080\u0094eot+22r/3) e j(O\u00E2\u0080\u0094aI+2,r/3\u00E2\u0080\u0094nr/2) 1 The factors k(t) and k4(t) are explained in Section 4.2. By substituting (4.18) into (4.17), we get the equivalent circuit of the induction machine (4.11) in the SFA domain as follows Vabcs(t) _ Req Iabcs(t)+EhQ) (4.19) where Req = r5 + + + k(t) (Equivalent resistances) Eh (t) = Ehl (t) + Eh 0\u00E2\u0080\u0099) (History term) The SFA equivalent circuit for the induction machine is illustrated in Figure 4.6. 72 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis Figure 4.6 SFA Equivalent Circuit of an Induction Machine Discretizing the rotor mechanical part equations [9], we can get the discrete time equations of the mechanical part O)r(t) = cr(t\u00E2\u0080\u0094At) (4.20) I (t) = f[. (t)lq (t) + Aq3 (t)idS (t)] (4.21) Or (t) = Or (t \u00E2\u0080\u0094 At) + [C0r (t) + cUr (t \u00E2\u0080\u0094 At)] (4.22) where p is the number ofpoles, \u00E2\u0080\u0098qdOs K(Or(t))Re[Iabcs(t)e\u00E2\u0080\u0099] Lilq \u00C2\u00B1 mq = Lislqs + LM (1qs + qr) Aas = LIS1d + kd = LiSidS + LM (Ids + dr) and Iqar(t) can be updated by using (4.9). The more complex mechanical part equations based on the multi-mass model can be found in [9]. 4.3.3 Simulation Results The same three test cases described in Section 4.2.3 are simulated on a 1.83 GHz PC to illustrate the efficiency and accuracy of the SFA model. The first case illustrates free acceleration characteristics of a 3-hp induction machine and a 2250-hp induction machine. The second case is a load torque change test which shows the slow dynamics of the 3-hp machine. The third case simulates the fast transients caused by a three phase fault. The following subsections A, B, and C explain that the SFA-based model is an accurate and general-purpose one. Subsection D explores the very high efficiency of the SFA model. The author has developed 73 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis a MATLAB program based on the EMTP algorithm and the ER machine model. All SFA results in the following subsections have been compared with the results from this EMTP program in the same MATLAB environment. The induction machine parameters for the three cases can be found in Appendix A. Rotor parameters have been converted to the stator side. A. Simulation of Free Acceleration The start-up transients of the 3-hp and 2250-hp induction machines during free acceleration from stall are simulated. The total simulation times are 1 s for the 3-hp induction machine and 3 s for the 2250-hp induction machine, respectively. In this simulation, a time step of & = 500 ps is used. Figure 4.7 illustrates the stator current of the 2250-hp induction machine during free acceleration. Note that the SFA model can produce both dynamic phasor and time domain results. The simulation results for the 3-hp induction machine are shown in Figure 4.8. The variables plotted in this figure are \u00E2\u0080\u0098as (phase A stator current), Vm (magnetizing flux linkage), Te (electromagnetic torque) and mechanical rotor speed. These machine variables are observed in the rotor reference frame; results viewed in other reference frames are omitted here. As is shown in Figure 4.8, the time domain values reconstructed from the dynamic phasor results are identical to those obtained from the EMTP algorithm using a ER model. The CPU times used for simulations are listed in Table 4.2. Figure 4.7 Stator Current during Free Acceleration of a 2250-hp Induction Machine (At 500 jts) 74 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis Table 4.2 CPU Times for Simulations Free Acceleration Case SFA Model ER Model 3-hp induction machine, T = 1 s, tt = 500 4us 0.728 s 0.556 s 2250-hp induction machine, T 3 s, & = 500 ,us 2.43 9 s 1.698 s * T stands for the total simulation time. B. Load Torque Change Test Slow dynamics during step changes in load torque are simulated in this case. The 3-hp machine originally operates in steady state with no load. A mechanical torque Tm of 12 Nm is applied at t = 2.05 s. Then Tm is reversed to \u00E2\u0080\u009412Nm at t = 2.5 s. The simulation time is 3 s with a Time (s) Figure 4.8 Dynamic Performance of a 3-hp Induction Machine during Free Acceleration (At = 500 jts) 75 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis time step of 0.5 ms for the ER model, and a time step of 5 ms for the SFA model. Figure 4.9 shows the phase A stator current and the slip during load torque changes. Figure 4.9 shows that simulation results from the SFA model with At = 5 ms are as accurate as those from an EMTP solution with At = 0.5 ms. The CPU time for a 3 s simulation with the SFA model is 0.266 s, while the simulation with the ER model uses 1.696 s. Therefore, to achieve the same numerical accuracy, the SFA induction machine model can be more than 5 times faster than its corresponding ER model in this case. Figure 4.9 Dynamic Performance of a 3-hp Induction Machine during Step Changes in Load Torque (RE At = 0.5 ms, SFA At =5 ms) C. Three Phase Fault Test In this simulation of fast transients, the 2250-hp induction machine is first operating 76 2.5 Time (s) Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis under rated conditions with a load torque equal to TB. A 3-phase fault is then applied at the terminals at t = 6.1 s. After 6 cycles, the fault is cleared. The simulation time is 10 s with a time step of 500 ps. Figure 4.9 illustrates that the proposed SFA model can simulate fast dynamics of the machine as accurately as the ER model. The CPU times are 7.258s and 5.562 s for SFA simulation and ER simulation, respectively. D. Comparison between SFA and an EMTP Phase-Domain Model Cases A to C show that the SFA model is as accurate as EMTP phase-domain models for the same time step i.t. Now we will carefully explore the numerical accuracy and robustness of Time (s) Figure 4.10 Simulation Results for a 3-Phase Fault at the Terminals of a 2250-hp Induction Machine (At 500 J4s) 77 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis the SFA model, by using different time steps At in the slow dynamics simulation. The 2250-hp machine originally operates in steady state with rated torque of 8900 Nm. The machine dynamics are recorded when a sudden reverse of mechanical torque to -8900 Nm is applied at t = 6 s. The total simulation time is 10 s. The SFA model was run with different time steps At = 0.lms, At = lms, At = 5ms, and At = lOms. The results obtained from an original phase-domain model (refer to p. 142-147 in [43]) with a small At = 0.1 ms are assumed to be the reference accurate time domain results. A zoomed-in view of phase A stator current is shown in Figure 4.11, which shows that the SFA model can accurately simulate the system dynamics with a large time step of 10 ms. In fact, even with a time step of 15 ms, the numerical solution with SFA remains stable. On the contrary, an EMTP simulation in the time domain with a time step of more than 4 ms became numerically unstable, as is shown in Figure 4.12. 5.8 6 6.2 6.4 6.6 Time (s) Figure 4.11 Simulation Results by SFA 6.8 / 7.2 78 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis $ Figure 4.12 Time Domain Results from an EMTP Algorithm Implemented in MATLAB Table 4.3 summarizes the CPU times needed for the SFA and ER simulations in this case. It can be seen that the SFA model was about 10 times faster (0.370 s vs 2.80 1 s) than the EMTP phase-domain model (ER model) to achieve a comparable accuracy. In certain situations, for instance, when simulating a system with power electronics components, an EMTP-type simulator may have to use very small time steps, for example, less than 0.1 ms. As can be seen from Table 4.3 that a SFA model with zt = lOms can be about 75 times faster than an EMTP model with & = 0.lms, for reasonable accuracy. Therefore, it is very promising to apply the SFA model in hybrid simulators, where the power electronics components as well as the network can be simulated with a small time step, while the SFA can be used for machine models with a larger time step. Further, the new SFA model has the capability to be extended to model the saturation and deep bar effects in induction machines. The SFA modeling framework is shown to be successful for machine modeling as a new method, which opens the door for more future research such as applying the SFA to wind turbine generator modeling, and developing SFA based EMTP simulators or hybrid simulators, and so on. Table 4.3 CPU Times for Simulations 2250-hp induction machine SFA Model ER ModelTorque Change Case T=lOs,&0.lms 36.108 s 27.328s T=lOs,&=lms 3.725s 2.801s T=lOs,&5ms 0.738s -- T= 10 s, & = 10 ms 0.370 s * Numerically unstable. 79 Chapter 4 Induction Machine Modelling Based on Shifted Frequency Analysis 4.4 Summary (1) An Equivalent-Reduction-based induction machine model is proposed. With a lower number of model equations, it is actually a full-order phase-domain model but is simpler and faster. Park\u00E2\u0080\u0099s transformation is implicitly used in the ER model in order to maintain the elegance and simplicity in the discrete-time equivalent circuit. As a non-ad hoc approach, the ER can also be used for building more efficient models for other power system elements. (2) Based on the ER model, a SFA model is proposed as a general-purpose model capable of simulating both fast transients and slow dynamics. Case study results have confirmed the SFA induction machine model is a valuable component for real-time EMTP simulations. It is observed that the SFA model is in excess of 70 times faster than the EMTP solution with the RE model when simulating dynamics with frequency spectra close to the fundamental power frequency. 80 Chapter 5 EMTP Implementation Chapter 5 EMTP Implementation 5.1 Introduction Since a number of power system component models in the shifted frequency domain have been built in the previous chapters, it is now time to ask how they can be used in a general purpose simulation tool based on the Shifted Frequency Analysis. EMTP, which was originally developed for calculating the transient overvoltages in transmission systems, has been significantly expanded to tracing the evolution of the system states in arbitrary multi-phase power networks consisting of all types of components. With improved functionality, accuracy and numerical stability, EMTP has become a standard tool being widely used in the power industry for system planning and designing purposes. Now EMTP is seeing broader applications in power system steady state [58] and dynamics studies and will remain one of the mainstreams in power system research. The intention of this thesis is to expand EMTP to efficiently simulate the slow dynamics in power systems, and to bridge the gap between the EMTP and a unified power system analysis tool. This means that the EMTP algorithm will be adopted in the SFA-based simulation tool. What makes a difference here is that the electrical variables are described by dynamic phasors instead of instantaneous time values. 81 Chapter 5 EMTP Implementation 5.2 Program Structure The implicit trapezoidal rule of integration, which has attractive characteristics in terms of accuracy and numerical stability, is used in EMTP and is followed in the SFA-based modelling of system components. By using the trapezoidal rule, the differential equations representing all network components are converted into algebraic relationships which can be interpreted as equivalent resistances or admittances, voltages, currents, and known history terms. Then the nodal equations of the system can be written directly as follows. [G][V(t)] = [I(t)]+ [H(t)] (5.1) where [G] is the matrix of the nodal equivalent admittances [VQ)] is the vector of the nodal voltages, which are dynamic phasors [1(t)] is the vector of the nodal currents, which are dynamic phasors [11(t)] is the vector of the nodal history terms The voltages of the nodes connected to voltage sources are known quantities, therefore the corresponding equations can be eliminated. Suppose A is the index set denoting the nodes with unknown voltages, and B is the index set for nodes with known voltages, then the nodal equations can be written in a block matrix form rGM G1EvA(t)1 = EIA(tHAt)1 (5.2) [GBA GBB][VB(t)J [IB(t)+HB(t)] Thus the unknown voltages can be obtained at time t by solving the following equations GVA(t) [IA(t) + HA(t)]\u00E2\u0080\u0094GVB(t) (5.3) In this thesis, a toolbox for simulating power system dynamics in the shifted frequency domain has been developed with MATLAB. This dynamic phasor tool consists of several files, which are depicted in the schematic structure in Figure 5.1. A. Input Data File Input data file provides the data needed for the simulation of a multi-phase power network. These include transmission line data, load data, machine data, switching operations 82 Chapter 5 EMTP Implementation data, voltage sources data, etc. All these data are documented in a MATLAB script file and most of the data are provided as the MATLAB arrays. Refer to Appendix C for the input data format. System Solver Main.m PiLineParameters.m MacParameters.m Data Processing LoadParameters.m RLCParameters.m EMTP Dynamic Phasor Simulation Tool (EMTDP) Equivalent Admittance Matrix Gsys.m Forming Updating UpdateG.m Admittance Matrix and History Terms UpdateHist.m Data Input File Case_i 3Bus.m Figure 5.1 Schematic Structure for the EMTP Dynamic Phasor Simulation Tool B. Data Processing Files When the input data are read by the main program, the data processing functions are called to establish the equivalent circuit for different power system components, and initialize their history terms for the simulation at the first time step. The equivalent circuit for different components can be found in previous chapters. The initialization in the SFA simulations is based on the snapshot method. First, let the dynamic phasor program run with zero initial conditions and reaches the steady state. Then, a snapshot of.the system is taken by saving system variables and history terms at a particular time step. These system variables and history terms from the snapshot file are fed to the SFA program as the initial conditions. With system variables and history terms initialized, the SFA simulation 83 Chapter 5 EMTP Implementation will run with a \u00E2\u0080\u0098flat\u00E2\u0080\u0099 start, and the advantage of using large time step in simulating 60 HZ dynamics is achieved. C. File for Building G Matrix The file \u00E2\u0080\u0098Gsys.m\u00E2\u0080\u0099 builds the system [G] matrix in equation (5.1) and returns the [Gj and [Gp.i] to the system solver for calculating the unknown voltages. The following explains how this simulation tool inserts power system components into equivalent admittance matrix [G] by taking the it-circuit model of a transmission line as an example. Suppose there is a three-phase it-circuit for a transmission line connecting node set J (=[ji, j2, j3}) and node set K (=[k1,k2,k3]) as depicted in Figure 5.2. As can be seen in Figure 5.2, the it-circuit model consists of a coupled RE branch and two coupled capacitances. When discretized by the trapezoidal rule, the coupled RE branch will contribute a 3 x3 equivalent admittance matrix to [G], and contributes a 3x1 history term hpj to [H]. Similarly, the coupled capacitances will also contribute a 3 x 3 equivalent admittance matrix Gc and a 3 xl history term hc. From equations (2.21) and (2.24), it can be found that GRL will be entered into [G]j,j and [G]j\u00C3\u00A7K in the system admittance matrix [G], and \u00E2\u0080\u0094G1j. into [G]J,K and [G]j\u00C3\u00A7j. For instance, G(ll) will be added to the element Ggi(1,2) will be added to [GJJ1i2, GII(l,2)to [G]13, GRI(2,1) to [G]j21, and so on. This procedure is illustrated in Figure 5.3. The MATLAB implementation of this process is shown in Figure 5.4. This routine also moves the nodes connecting to the voltage sources to the bottom of the system equation set, and then extract the sub-blocks [Gj and [Gpj] from [G]. C13 Figure 5.2 The fl Transmission Line Model 84 Chapter 5 EMTP Implementation D. File for Updating [G] and [H] The MATLAB files UpdateG.m will update the G matrix when changes happen in the network configuration, e.g. line tripping or faults. The file UpdateHist.m will calculate the history terms at each time step and return them to the solver for calculating unknown variables. E. System Solver The system solver Main.m is the core function in the simulation tool. In each time step, this routine assembles the right hand side of (5.3), performs downward operations on it, and does the back-substitution to obtain the unknown nodal voltages. After the unknowns are found, the routines for updating history terms and/or admittance matrix are called, in preparation for the calculations at the next step. The simulation run will continue until the total simulation time is reached. A schematic flow chart of the system solver is presented in Figure 5.5. 111213 k123 \ +Gc +G -G \___ --- / / -G +G + G1 [Gj [hj Figure 5.3 Contributions of the fl-Circuit Transmission Line Model to the Nodal Admittance Matrix 32 33 k1 Ic2 k3 jl 32 33 k1 Ic2 Ic3 + hcj /hRL + hc K 85 Chapter 5 EMTP Implementation % Pi Lines Contribution to G Matrix for kk = 1: nLine ii = (LineFromBus(kk)- 1 )*3 + 1; jj = (LineToBus(kk)1)*3 + 1; Gsystem( ii: ii + 2,11: ii + 2) = Gsystem( ii: ii + 2,ii: ii + 2) + GRL(:,:, kk) + GC(:,:, kk); Gsystem( ii: ii + 2,jj: jj + 2) = Gsystem( ii: ii + 2,jj: jj + 2) - GRL(:,:, kk); Gsystem( jj: jj + 2,11: ii + 2) = Gsystem( jj: jj + 2,ii: ii + 2) - GRL(:,:, kk); Gsystem(jj: jj + 2,jj:jj + 2) = Gsystem(jj:jj + 2,jj:jj + 2) + GRL(:,:, kic) + GC(:,:, kk); end Figure 5.4 The MATLAB Code for Inserting H-Circuit Transmission Line Model into G Matrix 5.3 Test Cases Four test cases are tested in the SFA-based EMTP simulator developed in this chapter. A. Radial Transmission Line Case The first case is a radial network consisting of a voltage source, a three phase load and a double-circuit overhead line between the source and the load. The one-line diagram is shown in Figure 5.6. The transmission line is energized at t = 0. At t = 3s, one of the parallel circuits is tripped. The total simulation time is 5 seconds. A time step of 5 ms is used in this simulation. Detailed data for this case are included in Appendix C. Figure 5.12-Figure 5.12 illustrates the three phase voltages at the load node, including the zoomed-in view for these voltages in the chosen time interval from 2.9s to 3.ls. A MATLAB program based on the EMTP pi-circuit model is used to generate the time domain results for comparison. It can be seen that the EMTP solutions in the shifted frequency domain and those in the time domain are almost identical, and the dynamic phasor result is the envelop of the time domain curve. 86 Chapter 5 EMTP Implementation ( Input Data Find Steady State Solution to Initialize History Terms and Variables at I = Jr Build GAA, GAB for Transient Solution 1. Factor G into LU Form Jr Start Time Step Simulation = Evaluate Current and Voltage Sources, __________ Update the Right Hand Side of(5.3) J,N I=t+Ai Perform Downward Operations on RHS of(5.3), Do Back-substitution to Solve Nodal Voltages \u00E2\u0080\u0098I, Update History Terms R,X Load Any Switching \u00E2\u0080\u0098\u00E2\u0080\u0098 Modify and Output Results Figure 5.5 Flow Chart for Dynamic-Phasor-Based EMTP Simulator Source Busi Bus2 Figure 5.6 One-line Diagram of a Radial Test System 87 Chapter 5 EMTP Implementation \u00E2\u0080\u0098C I 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 Time (s) Figure 5.9 Phase B Voltage at the Load Node (At = 5 ms) 88 SFA Solution \u00E2\u0080\u0098nfl \u00E2\u0080\u0094l 0.5 1 1.5 2 2.5 3 3.5 4 45 5 Time (s) Figure 5.7 Phase A Voltage at the Load Node (At = 5 ms) Figure 5.8 Zoomed-in View of Phase A Voltage at the Load Node (At = S ms) \u00E2\u0080\u0098 r 1 Chapter 5 EMTP Implementation 0.5 > > -0.5 Li \u00E2\u0080\u00941 ,i i 2.9 2.92 2.94 2.96 2.98 Time(s) (b) 3.02 3.04 3.06 3.08 3 1 0 0.5 I 1.5 2 2.5 3 3.5 4 4.5 5 Time (s) Figure 5.11 Phase C Voltage at the Load Node (At =5 ms) x I I I SFA Solution EMTP Solution 0.5 0 -0.5 \u00E2\u0080\u0094l 2.9 - 2.92 2.94 2.96 2.98 3 302 3.04 3.06 108 3.1 Time (s) Figure 5.12 Zoomed-in View of Phase C Voltage at the Load Node (At = 5 ms) 89 I I I I I I .1 SFA Solution Figure 5.10 Zoomed-in View of Phase B Voltage at the Load Node (At = 5 ms) x ___________ 1 IVVVVV\IVVVVV ChapterS EMTP Implementation B. Distribution Network Case The second case is a 13 bus network with the configuration adapted from the IEEE test feeder [59]. The network configuration is shown in Figure 5.13. The detailed network data are documented in Appendix C. The network is energized at t = 0. The induction machine load connected to node 11 is initially operating at no-load. At t = 3s, a rated mechanical torque is applied on the induction machine. Then one circuit of the double-circuit lines connecting node 0 and node 1 is tripped at t = 6 s. The total simulation time is 8 s. A time step of 1 ms is used in this simulation. Figure 5.14-Figure 5.18 illustrate part of the three phase voltages at three different nodes. A MATLAB program using the EMTP algorithm is used to produce the time domain results. The time domain simulation results and the dynamic phasor results are both shown in each figure to illustrate the correctness of the dynamic phasor results. From the figures, it can be seen that the increase of the induction motor load causes drops in the system voltage. Loss of one of the double circuit lines further weakens the systems and causes voltage dip problems. 5 4 0 1 2 3 11\u00E2\u0080\u00A2 Figure 5.13 One Line Diagram for the Test Feeder 90 ChapterS EMTP Implementation > Figure 5.15 Zoomed-in View of Phase A Voltage at the Induction Machine Node (At = 1 ms) 3 3.5 4 45 5 5.5 6 6.5 7 7.5 8 Time (s) Figure 5.16 Phase B Voltage at Load Node 6 (At = 1 ms) 91 \u00E2\u0080\u0094 SFA Solution EMTP Solution 3 3.5 4 4.5 5 5.5 6 6.5 7 7.5 Time (s) Figure 5.14 Phase A Voltage at the Induction Machine Node (At = 1 ms) I SFA Solution EMTP Solution I I I I I I I I Chapter 5 EMTP Implementation > I I I I I 3.5 4 4.5 5 5.5 6 Time (s) Figure 5.18 Phase C Voltage at Node I (At = 1 ms) C. First Benchmark System for the Subsynchronous Resonance Studies This test system was prepared by an IEEE Subsynebronous Resonance Task Force [601 as a standard test case for computer programs to simulate subsynchronous resonance phenomena. The test system consists of an 892.4 MVA turbine-generator connected through a step-up transformer to a 500 kV transmission line with series capacitor compensation. The power system at the receiving end is represented by a Thevenin equivalent circuit (infinite bus behind reactance) [61]. Figure 5.19 shows the one-line diagram for the first benchmark system. In Figure 5.17 Zoomed-in View of Phase B Voltage at Load Node 6 (At = 1 ms) SFA Solution EMTP Solution .4OQO 6.5 92 Chapter 5 EMTP Implementation Figure 5.19, all data are represented in p.u. based on 892.4 MVA and 500 kV. The generator parameters of the electrical part are listed as below Xci = 1.79 p.U. Xq 1.71 p.U. Ra= 0 X\u00E2\u0080\u0099d = 0.169 p.U. X\u00E2\u0080\u0099q = 0.228 p.U. X1 0.13 p.U. X\u00E2\u0080\u009Dd = 0.135 p.U. X\u00E2\u0080\u009Dq = 0.2 p.u. f 60 Hz T\u00E2\u0080\u0099cio = 4.3 S T\u00E2\u0080\u0099qo 0.85 S it(O) = 1 f.U. T\u00E2\u0080\u009Ddo 0.032 S T\u00E2\u0080\u009Dqo = 0.05 5 The original purpose of this case was to simulate the interaction between the mechanical torque placed on the generator turbines and the electrical torque related to the power network, and the resulting shaft torsional oscillations. Accordingly, a detailed multi-mass modcl of the mechanical shaft is adopted in [601, [611, where the generator shaft system is modeled by 6 masses including 4 turbine sections HP, IP, LPA and LPB, 1 generator and 1 exciter. This thesis focuses on the feasibility of applying the SFA method, and thus uses a simpler single-mass representation of the mechanical part. The total inertia constant of the turbine and generator is 2.894 seconds. The SFA synchronous machine model in Chapter 3 is used in this test case. It should be noted that the SFA synchronous machine model of Chapter 3 can be extended to model the multi-mass shaft system by adding differential equations for individual spring-masses to the differential equations of the shaft system. The self and mutual damping effect can also be easily included in the shaft system equations. The detailed network parameters can be found in [601 and [611. 26kV/500kV A R1 0.02 X1 = 0.50 B X1 = 0 060=0.50 0=l.56 \u00E2\u0080\u0094 I Generator = 0.14 \u00E2\u0080\u0094 Infinite Bus X0=0.14 Xc-O.371 X=0.04 1 0 = 0.04 Figure 5.19 The First Benchmark Network for Subsynchronous Resonance Studies 93 Chapter 5 EMTP Implementation In this test, a three-phase fault occurs at bus B at t = 3 s. After 4 cycles, the fault is cleared. The total simulation time is 3.5 s. A time step of 1 ms is used in this simulation. The whole system is represented using the SFA models proposed in the previous chapters. The generator terminal voltages, the voltages at bus A and bus B, and the voltages across the series capacitors are monitored during the SFA simulation. Figure 5.20 to Figure 5.23 illustrate the dynamic phasor results for the system voltages at different locations as well as the time domain results transformed back from the corresponding dynamic phasors. The dynamic phasor results show that, after the fault is applied and cleared, there are low frequency oscillations happening at the generator terminals, on the transformer sides, and across the series capacitor banks. \u00E2\u0080\u0098U \u00E2\u0080\u0098U C C \u00E2\u0080\u0098U F C C \u00E2\u0080\u0098U C \u00E2\u0080\u0098U time (s) Figure 5.20 Generator Terminal Voltage: Phase A (At = 1 ms) 94 Chapter 5 EMTP Implementation 0 0 .0 I Figure 5.21 Transformer High Side (Bus A) Voltage: Phase A (At = 1 ms) Figure 5.22 Voltage across Series Capacitor: Phase A (At = 1 ms) These oscillations are an electrical phenomenon because the shaft system is modeled as a single mass. The slow oscillations are caused by resonances with the series capacitor, which are excited by the fault and the switching operations. The resonance mode(s) is determined by the inherent characteristics of the power network. In fact, the modes or natural frequencies can be quantitatively found by either performing a frequency scan (steady-state solutions over a frequency range) or by calculating the eigenvalues of the admittance matrix of the electrical 95 3 time (s) Chapter 5 EMTP Implementation network. Similar to the network resonance analysis, if the shaft system is represented by the detailed multi-mass model, we can also determine the torsional natural frequencies (eigenvalues) and mode shapes (eigenvectors) by applying the modal analysis to the shaft differential equations. If the complement of the natural frequency of the network is close to one of the torsional frequencies of the shaft system, torsional oscillations will be excited [331. This is the mechanism of the subsynchronous resonance. The dynamic phasor results obtained from the SFA simulations show that the oscillations with frequencies around 30 Hz have been excited in this benchmark system. On the other hand, the natural frequencies of the 6-mass shaft system were found to be 15.71 Hz, 20.21 Hz, 25.55 Hz, and 32.28 Hz [601. That means the complement of the natural frequency of the series-compensated network is close to the torsional natural frequencies, which would result in a subsynchronous resonance in this system and would build up torsional oscillations on the shaft. This has been verified by EMTP simulations [611. In summary, the SFA simulations are able to capture the slow oscillations in the benchmark system and can be used for subsynchronous resonance studies once the detailed shaft model is incorporated into the SFA synchronous machine models. U, 0 > U, 0 0 to D. Second Benchmark System for the Subsynchronous Resonance Studies Eight years after the first subsynchronous resonance benchmark was published, the IEEE Subsynchronous Resonance Working Group proposed a second benchmark [621. The second 96 time (s) Figure 5.23 Infinite Bus (Bus B) Voltage: Phase A (At = 1 ms) Chapter 5 EMTP Implementation benchmark is a case consisting of a 600 MVA generator connected through a step-up transformer to two parallel lines, one of which is series compensated. The compensation rate of the line is 55%. The power system at the receiving end is represented by a Thevenin equivalent circuit (infinite bus behind impedance). Figure 5.24 shows the one-line diagram for the second benchmark system. All data in Figure 5.24 are represented in p.u. on a 100 MVA, 500kV base. The mechanical part of the round rotor turbine generator is again represented as a single-mass model with a total inertia of 2.683 seconds in the SFA simulation. In [621, a more detailed multi- mass representation is used with 4 masses. The parameters of the electrical part are based on the machine ratings, which are given below. Xd = 1.65 p.U. Xq = 1.59 p.u. Ra = 0.0045 3.U. X\u00E2\u0080\u0099d = 0.25 p.u. X\u00E2\u0080\u0099q = 0.46 p.U. Xi = 0.14 p.U X\u00E2\u0080\u009Dd = 0.20 p.U. X\u00E2\u0080\u009Dq = 0.20 p.u. f 60 Hz T\u00E2\u0080\u0099do = 4.5 S T\u00E2\u0080\u0099qo = 0.55 S if(O) = 1 3.U. Tdo = 0.04 S T\u00E2\u0080\u009Dqo = 0.09 s number of poles = 2; The detailed benchmark network parameters can be found in [62] and [61]. In this test, a three-phase fault is applied at the high voltage side of the transformer connecting to bus 2 at t = 3 s. The fault clearing time is 1 cycle, and the total simulation time is 3.5 s. The time step used in this simulation is 1 ms. Figure 5.25 to Figure 5.28 show the dynamic phasor results for the system voltages at different locations together with the time domain results transformed back from the corresponding dynamic phasors. Obviously, oscillations slower than the fundamental frequency are excited in the system after the fault occurs, which can be seen from the dynamic phasors of the generator terminal voltage, the transformer high side voltage and the voltage across the series capacitor. Note that bus 1 does not see large deviations in the voltage because it is electrically close to the infinite external system, and also because the series capacitor is acting as a highpass filter that blocks the slow oscillations in the system. The test results on the second benchmark system further verify that the SFA method can capture the so called \u00E2\u0080\u0098parallel resonance\u00E2\u0080\u0099 phenomenon [62] in the meshed power system. Note that relatively small time-step is used for the simulation of the subsynchronous resonance because large frequency deviation occurs after the fault is applied. The reason why a small time step has to be used is that large frequency deviation defines a wide bandwidth in the shifted frequency domain, which in turn requires smaller time steps to respect the Nyquist 97 Chapter 5 EMTP Implementation frequency limit in the shifted frequency domain. This indicates that a variable time step scheme is a future research direction to realize full potential of the SFA method. The time step can be reduced when system states are changing rapidly in order to achieve better accuracy in the SFA simulation. On the other hand, when system dynamics slowdown, large time steps can be used to avoid unnecessarily long computational time while still achieving reasonable accuracy. Generator R=0.0002 X=0.02 Bus 2 R1 = 0.0067 0=0.0186 R1 = 0.0074 0.022 Xi = 0.0739 0=0.21 Bus 1 Xi = 0.08 0=0.24 Bus C Xiijne Ri =R00.0014 X1 = X0= 0.03 Infinite Bus Figure 5.24 The Second Benchmark Network for Subsynchronous Resonance Studies 98 Chapter 5 EMTP Implementation V C: -C V C: C > I V C: C: -C V C C U Ct U C: V V C\u00E2\u0080\u0099) Figure 5.25 Generator Termiual Voltages: Phase A (At = 1 ms) Figure 5.26 Voltage across Series Capacitor: Phase A (At = 1 ms) time (s) 3 time (s) 99 Chapter 5 EMTP Implementation 0 \u00E2\u0080\u00A20 CD I I > 0 ci, 0 Figure 5.27 Transformer High Side (Bus 2) Voltage: Phase A (At = 1 ms) Figure 5.28 Infinite Bus (Bus 1) Voltage: Phase A (At = 1 ms) The test results from all the above cases indicate that the dynamic phasor is a generalization of the phasor concept, which can represent the dynamic waveform in power systems, without loss of important information. The SFA method with its implementation in the EMTP environment can integrate the differential equations of the power system in the SFA domain and can produce dynamic phasors for electrical variables, which are visually clear and 100 2.5 time (s) Chapter 5 EMTP Implementation easy to follow for power engineers. With the SFA method, power engineers may gain better insight into the EMTP simulation results, and we would expect broader applications of the EMTP in power system steady state and dynamics studies, beyond the fast transient simulations. 101 Chapter 6 Conclusions Chapter 6 Conclusions 6.1 Summary of Contributions The goal of this thesis is to extend EMTP functionality for power system dynamic simulation, especially for simulating dynamics with frequency spectra close to the fundamental power frequency. This has been accomplished by developing the Shifted Frequency Analysis (SFA) method, modeling system components with SFA, and accelerating the EMTP simulations for dynamics around 60Hz. A series of contributions made in this thesis are the following. I. The theory of Shifted Frequency Analysis is proposed with the help of Hubert transform and analytic signal concept. Numerical accuracy analysis is performed for the discrete- time SFA simulation. II. Linear circuit components, transformer, exponential load and steady-state induction motor are modeled in the shifted frequency domain. III. An efficient SFA-based synchronous machine model for the simulation of slow system dynamics is developed. This model is a general-purpose model that can be used for evaluating the dynamic performance of both the salient-pole and the cylindrical-rotor machines. IV. The SFA method is extended to model the induction machines. This thesis proposes a new phase-domain induction machine model based on the equivalent-reduction (ER) approach. The ER model has a concise discrete-time equivalent circuit that can be directly incorporated into EMTP-type simulators. Based on the ER model, a SFA induction machine model is proposed as a general purpose model capable of simulating both fast transients and slow dynamics. V. An EMTP simulation tool based on the SFA is developed. Simulation results 102 Chapter 6 Conclusions validate that the SFA method is capable to efficiently simulate power system fundamental frequency dynamics. This is the first practical accomplishment to build a unified power system analysis tool based on the EMTP solution. 6.2 Future Research The accomplishment achieved in this thesis will inspire researchers to apply SFA and associated techniques to power system analysis and other potential areas. The success of the applications of SFA may be achieved in different aspects such as: Apply the SFA to the modeling and analysis of renewable energy resources, dispersed generation plants, and industrial power systems. In the foreseeable future, more and more independent-power-producer-owned (IPP-owned) generation units will be connected into the power transmission and distribution system due to a deregulated electricity market. Many of these plants will be cleaner, and many of them will be renewable energy sources such as biomass, solar, wind, geothermal, small hydro, ocean energy, and so on. The IPP interconnection impact studies, which identify the system constraints brought about by the integration of IPPs, and determine network upgrades and remedial action schemes, are therefore critical for the reliable and safe operation of the TPPs and associated power systems. By applying the SFA technique in the IPP modeling, it may lead to accurate and efficient solutions in the IPP impact study, which is hardly achievable with the current phasor tools. II. Develop new component models and controller models, and improve the SFA-based EMTP solver. In addition, the future SFA simulator will adopt more object-oriented design such that each module for new components can simply \u00E2\u0080\u0098plug into\u00E2\u0080\u0099 the simulation engine without modifying the core codes. The robustness and efflciency of the core codes will also be improved by using up-to-date sparsity techniques. III. The equivalent-reduction idea can be used for building more efficient models for various types of machines, transformers and power system equipments, as long as they can be mathematically represented by differential equations with mutually coupled state variables. These new models can then be further explored in the shifted frequency domain. 103 Chapter 6 Conclusions IV. Transform the SFA-based simulation tool into parallel programs. The UBC Object Virtual Network Integrator (OVNI) [631 is a real time parallel simulator, which uses PC clusters [64] as hardware and Muti-Area Thevenin Equivalent (MATE) algorithm as solution engine [65]-[69]. MATE partitions the power system into subsystems and solves them in parallel. By changing the real-valued component models to the dynamic phasor models, a new OVNI simulator based on SFA can be implemented. In the future, one may expect this new simulator to serve as a distributed simulation tool for the supervision and control of self-healing power infrastructures. V. Further investigate the basic theory of the SFA method. This will still be an interesting field in future research. There are two limitations in the current SFA method. First, SFA simulation is computationally more expensive than EMTP simulation for the same integration step At. This is due to the higher computational cost of operating with complex numbers in SFA as compared to real numbers in EMTP and also to the computational cost involved in transferring between shifted frequency domain and time domain. Second, the aliasing effect may occur when simulating very fast transients in a system [24]. This may introduce error or distortion in SFA simulation result. Some new theories such as discrete-time analytic signal, Hilbert-Huang transform, and new antialiasing techniques are likely to bring theoretical breakthroughs and may lead to the next generation SFA method, which may be more flexible and accurate than the current SFA. VI. It may be worthwhile to look into transformations to d,q,0-quantities as an alternative to dynamic phasors [70], particularly for cases where the three-phase impedances are balanced, and where the faults are symmetrical three-phase faults. VII. Last but not least, the applications of the SFA method is not limited to the EMTP solution. The Shifted Frequency Analysis can also be introduced into any other circuit simulators, for instances, some state-space-based circuit and/or electric machinery solvers, SPICE, and so on. It is equally suitable to be implemented in the variable time step simulators such as Simulink\u00C2\u00AE. It can also be used in hybrid simulations, and hardware-in-loop simulators. All in all, the SFA is believed to have opened new doors of opportunity for research relating to dynamic simulations, and further theoretical work and applications are justified. 104 Bibliography Bibliography [1] U.S.-Canada Power System Outage Task Force, Final Report on the August 14, 2003 Blackout in the United States and Canada: Causes and Recommendations, Apr. 2004. [2] G. T. Heydt, C. C. Liu, A. G. Phadke, V. Vittal, \u00E2\u0080\u009CSolution for the Crisis in Electric Power Supply,\u00E2\u0080\u009D IEEE Computer Applications in Power, vol. 14, no. 3, pp. 22-30, Jul. 2001. [3] IEEE/CIGRE Joint Task Force on Stability Terms and Definitions, \u00E2\u0080\u009CDefinition and Classification of Power System Stability,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 19, no. 3,pp. 1387- 1401,Aug.2004. [4] W. Li, Risk Assessment of Power Systems: Models, Methods, and Applications, New York: IEEE Press and Wiley & Sons, 2005. [5] S. R. Aumuri, L. R. Malone, V. Burtnyk, \u00E2\u0080\u009CRepresentation of Single-Pole Open Conditions in Stability Studies\u00E2\u0080\u009D, IEEE Transactions on Power Systems, vol. 6, no. 1, pp 9-15, Feb. 1991. [6] V. Ajjarapu, C. Christy, \u00E2\u0080\u009CThe Continuation Power Flow: a Tool for Steady State Voltage Stability Analysis\u00E2\u0080\u009D, IEEE Transactions on Power Systems, vol. 7, no. 1, pp 416-423, February 1992. [7] T. Van Cutsem, C.Voumas, Voltage Stability of Electric Power Systems, Kluwer Academic Publishers, 1998. [8] IEEE/PES Power System Stability Subcommittee, Voltage Stability Assessment: Concepts, Practices and Tools, Final Document, Aug. 2002. [9] H. W. Dommel, EMTP Theory Book, 2nd edition, Vancouver, BC: MicroTran Power System Analysis Corporation, 1996. [10] H. W. Dommel, \u00E2\u0080\u009CDigital Computer Solution of Electromagnetic Transients in Single- and Muti-Phase Networks,\u00E2\u0080\u009D IEEE Transactions on Power Apparatus and Systems, vol. 88, no. 2, pp. 734-771, Apr. 1969. 105 Bibliography [11] MicroTran Power System Analysis Corporation, MicroTran Reference Manual: Transients Analysis Program for Power and Power Electronic Circuits, Vancouver, BC, Canada, 2002. [121 Manitoba HVDC Research Center, User\u00E2\u0080\u0099s Guide of PSCAD, Winnipeg, Manitoba, Canada, 2003. [13] Canada/America EMTP Users Group, ATP Rule Book, 1998. [14] J. Mahseredjian, S. Dennetiere, L. Dube, B. Khodabakhchian, and L. Gerin-Lajoie, \u00E2\u0080\u009COn a New Approach for the Simulation of Transients in Power Systems,\u00E2\u0080\u009D Electric Power Systems Research, vol. 77, no. 11, pp. 15 14-1520, 2007. [15] DIgSILENT GmbH, DlgSlLENTPowerFactory V13 User Manual, 2005. [16] RTDS Technologies, Real Time Digital Simulator User Manual, 2003 [17] C. Henville, R. Folkers, A. Hiebert and R. Wierckx, \u00E2\u0080\u009CDynamic Simulations Challenge Protection Perfonnance,\u00E2\u0080\u009D Proceedings of the 30th Western Protective Relay Conference, Spokane, WA, Oct. 2003. [18] V. Venkatasubramanian, \u00E2\u0080\u009CTools for Dynamic Analysis of the General Large Power System Using Time-Varying Phasors,\u00E2\u0080\u009D International Journal ofElectrical Power and Energy Systems, vol. 16, no. 6, pp. 365-376, Dec. 1994. [19] M. Stubbe, A. Bihain, J. Deuse, J. C. Baader, \u00E2\u0080\u009CSTAG-A New Unified Software Program for the Study of the Dynamic Behaviour of Electrical Power Systems\u00E2\u0080\u009D, IEEE Transactions on Power Systems, vol. 4, no. 1, pp 129-138, Feb. 1989. [20] S. Henschel, \u00E2\u0080\u009CAnalysis of Electromagnetic and Electromechanical Power System Transients with Dynamic Phasors,\u00E2\u0080\u009D Ph. D. dissertation, The University of British Columbia, Vancouver, Canada, 1999. [21] A. M. Stankovi\u00C3\u00A1, T. Aydin, \u00E2\u0080\u009CAnalysis of Asymmetrical Faults in Power Systems Using Dynamic Phasors,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 15, no. 3, pp. 1062-1068, Aug. 2000. [22] A. M. Stankovi\u00C3\u00A1, S. R. Sanders, T. Aydin, \u00E2\u0080\u009CDynamic Phasors in Modeling and Analysis of Unbalanced Polyphase AC Machines,\u00E2\u0080\u009D IEEE Transactions on Energy Conversion, vol. 17, no. l,pp. 107-113, Mar. 2002. 106 Bibliography [23] A. M. Stankovi\u00C3\u00A9, B. C. Lesieutre, T. Aydin, \u00E2\u0080\u009CModeling and Analysis of Single-Phase Induction Machines with Dynamic Phasors,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 14, no. 1, pp. 9-14, Feb. 1999. [24] J. R. Marti, \u00E2\u0080\u009CShifted Frequency Analysis (SFA) for EMTP Simulation of Fundamental Frequency Power System Dynamics\u00E2\u0080\u009D, Internal Report, Power System Laboratory, The University of British Columbia, Vancouver, Canada, Apr. 2005. [25] 5. L. Hahn, Hubert Transforms in Signal Processing, Norwood, MA: Artech House, 1996. [26] A.V. Oppenheim and A. S. Wilisky, Signal and Systems, Upper Saddle River, NJ: Prentice-Hall, 1996. [27] R. Shintaku and K. Strunz, \u00E2\u0080\u009CBranch Companion Modeling for Diverse Simulation of Electromagnetic and Electromechanical Transients,\u00E2\u0080\u009D IPST\u00E2\u0080\u0099OS, Montreal, Canada, Jun. 2005. [28] P. Zhang, J. R. Marti, H. W. Dommel, \u00E2\u0080\u009CSynchronous Machine Modelling Based on Shifted Frequency Analysis,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 22, no. 3, pp. 1139- 1147, Aug. 2007. [29] P. Zhang, J. R. Marti, H. W. Dommel, \u00E2\u0080\u009CInduction Machine Modelling Based on Shifted Frequency Analysis,\u00E2\u0080\u009D IEEE Transactions on Power Systems, Accepted for Publication. [30] J. R. MartI, J. Lin, \u00E2\u0080\u009CSuppression of Numerical Oscillations in the EMTP,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 4, no. 2, pp. 739-747, May 1989. [31] L. W. Ehrlich, \u00E2\u0080\u009CComplex Matrix Inversion Versus Real,\u00E2\u0080\u009D Communications of the ACM vol. 13, no. 9, pp. 561-562, Sept. 1970. [32] L. Tornheim, \u00E2\u0080\u009CInversion of a Complex Matrix,\u00E2\u0080\u009D Communications of the ACM, vol. 4, issue 9, pp. 398, Sept. 1961. [33] P. Kundur, Power System Stability and Control, New York: McGraw-Hill, 1994. [34] BPA Transient Stability Program User\u00E2\u0080\u0099s Guide, Bonneville Power Administration, 1987. [35] J. R. MartI and K. W. Louie, \u00E2\u0080\u009CPhase-Domain Synchronous Generator Model Including Saturation Effects,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 12, no. 1, pp. 222\u00E2\u0080\u0094229, Feb. 1997. [36] S. D. Pekarek, 0. Wasynczuk, and H. J. Hegner, \u00E2\u0080\u009CAn Efficient and Accurate Model for the Simulation and Analysis of Synchronous Machine/Converter Systems,\u00E2\u0080\u009D IEEE Transactions on Energy Conversion, vol.13, no. 1, pp. 42-48, Mar. 1998. 107 Bibliography [37] L. Wang and J. Jatskevich, \u00E2\u0080\u009CA Voltage-behind-Reactance Synchronous Machine Model for the EMTP-Type Solution\u00E2\u0080\u009D, IEEE Transactions on Power Systems, vol. 21, no. 4, pp 1539\u00E2\u0080\u00941549, Nov. 2006. [38] J. Machowski, J. Bialek, and J. Bumby, Power System Dynamics and Stability, J. Wiley & Sons, 1997. [39] J. R. Marti, L. R. Linares, J. Calviflo, H. W. Dommel, J. Lin, \u00E2\u0080\u009COVNI: An Object Approach to Real-Time Power System Simulators, \u00E2\u0080\u009CProceedings of the 1998 International Conference on Power System Technology, Powercon\u00E2\u0080\u009998, Beijing, China, August 18-21, 1998, pp. 977-981. [40] F. Moreira, J. R. MartI, \u00E2\u0080\u009CLatency Techniques for Time-Domain Power System Transients Simulation,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 20, no. 1, pp. 246-253, Feb. 2005. [41] V. Brandwajn, \u00E2\u0080\u009CSynchronous Generator Models for the Analysis of Electromagnetic Transients\u00E2\u0080\u009D, Ph.D. Thesis, The University of British Columbia, Vancouver, BC, Canada, 1977. [42] 5. D. Pekarek, E. A. Walters, \u00E2\u0080\u009CAn Accurate Method of Neglecting the Dynamic Saliency of Synchronous Machines in Power Electronic Based Systems,\u00E2\u0080\u009D IEEE Transactions on Energy Conversion, vol. 14, no.4, pp. l1771184, Dec. 1999. [43] P. Krause, 0. Wasynczuk, and S. Sudhoff, Analysis of Electric Machinery and Drive Systems, 2nd ed., New York: IEEE Press, 2002. [44] IEEE Recommended Practice for Power Systems Analysis, IEEE Standard 399-1997, Sept. 1997. [45] J. Lin, J. Mahseredjian and S. Lefebvre, \u00E2\u0080\u009CImprovement of Synchronous Machine Saturation Simulation in the EMTP,\u00E2\u0080\u009D Proceedings of IEEE Region 10 Conference, vol 5, pp. 127-132, Oct. 1993. [46] R. Hung, H. W. Dommel, \u00E2\u0080\u009CSynchronous Machine Models for Simulation of Induction Motor Transients,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 11, no. 2, pp 833\u00E2\u0080\u0094838, May 1996. [47] IEEE Tutorial on Digital Simulation of Electrical Transient Phenomena, IEEE Catalog No. 81 EHO173-5-PWR, 1981. 108 Bibliography [48] G. J. Rogers, D. Shirmohammadi, \u00E2\u0080\u009CInduction Machine Modelling for Electromagnetic Transient Program,\u00E2\u0080\u009D IEEE Transactions on Energy Conversion, vol. 2, no. 4, pp 622\u00E2\u0080\u0094628, Dec. 1987. [49] R. H. Park, \u00E2\u0080\u009CTwo-Reaction Theory of Synchronous Machines\u00E2\u0080\u0094Generalized Methods of Analysis\u00E2\u0080\u0094Part I,\u00E2\u0080\u009D AIEE Trans., vol. 48, pp. 716\u00E2\u0080\u0094727, July 1929. [50] H. K. Lauw and W. S. Meyer, \u00E2\u0080\u009CUniversal Machine Modeling for the Presentation of Rotating Electric Machinery in an Electromagnetic Transients Program,\u00E2\u0080\u009D IEEE Transactions on Power Apparatus and Systems, vol. PAS-lOl, no. 6, pp 1342\u00E2\u0080\u00941351, Jun. 1982. [51] J. R. MartI and T. 0. Meyers, \u00E2\u0080\u009CPhase-Domain Induction Motor Model for Power System Simulators,\u00E2\u0080\u009D in Proc. IEEE Western Canada Conf Exhibition, pp. 276-282, 1995. [52] G. Kron, Tensor Analysis ofNetworks, New York: J. Wiley & Sons; London: Chapman &Hall, 1939. [53] P. Zhang, \u00E2\u0080\u009CSymmetrical Induction Machine Simulation Based on the Voltage-behind- Reactance Model\u00E2\u0080\u009D, Project Report, Power System Laboratory, The University of British Columbia, Vancouver, Canada, Dec. 2005. [54] D. Chaniotis and M. A. Pai, \u00E2\u0080\u009CModel Reduction in Power Systems Using Krylov Subspace Methods,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 20, no. 2, pp 888\u00E2\u0080\u0094894, May 2005. [55] J. Roychowdhury, \u00E2\u0080\u009CReduced-Order Modeling of Time-Varying Systems\u00E2\u0080\u009D, IEEE Transactions on Circuits and Systems\u00E2\u0080\u0094Il: Analog and Digital Signal Processing, vol. 46, no. 10, pp. 1273\u00E2\u0080\u00941288, Oct. 1999. [56] R. W. Freund, \u00E2\u0080\u009CKrylov-subspace Methods for Reduced-order Modeling in Circuit Simulation\u00E2\u0080\u009D Journal of Computational and Applied Mathematics, vol. 123, pp. 395\u00E2\u0080\u0094421, 2000. [57] A. C. Antoulas, D. C. Sorensen, and S. Gugercin, \u00E2\u0080\u009CA Survey of Model Reduction Methods for Large-Scale Systems,\u00E2\u0080\u009D Contemporary Mathematics, vol. 280, pp 193\u00E2\u0080\u0094219, 2001. [58] J. C. Fraga, \u00E2\u0080\u009CTunable Difference Equations for the Time Domain Simulation of Power System Operating States\u00E2\u0080\u009D, Ph.D. Thesis, The University of British Columbia, Vancouver, BC, Canada, 2003. [59] W. H. Kersting, \u00E2\u0080\u009CRadial Distribution Test Feeders,\u00E2\u0080\u009D IEEE Transactions On Power Systems, vol. 6, no. 3, pp. 975-985, Aug. 1991. 109 Bibliography [60] IEEE Subsyncbronous Resonance Task Force, \u00E2\u0080\u009CFirst Benchmark Model for Computer Simulation of Subsynchronous Resonance,\u00E2\u0080\u009D IEEE Transactions on Power Apparatus and Systems, vol. PAS-96, no. 5, pp. 1565-1572, September/October 1977. [611 H. W. Dommel, Case Studies for Electromagnetic Transients, 3rd edition, Vancouver, BC, Mar. 2004. [62] IEEE Subsynchronous Resonance Working Group, \u00E2\u0080\u009CSecond Benchmark Model for Computer Simulation of Subsynchronous Resonance,\u00E2\u0080\u009D IEEE Transactions on Power Apparatus and Systems, vol. PAS-104, no. 5, pp. 1057-1066, May 1985. [63] J. R. MartI, L. R. Linares, J. A. Hollman, F. A. Moreira, \u00E2\u0080\u009COVNI: Integrated software/Hardware Solution for Real-time Simulation of Large Power Systems,\u00E2\u0080\u009D in Proceedings of the PSCCO2, Sevilla, Spain, June, 2002 [64] J. A. Hollman, J. R. MartI, \u00E2\u0080\u009CReal Time Network Simulation with PC-Clusters,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 18, no. 2, pp. 563-569, May 2003. [65] M. Armstrong, J. R. MartI, L. R. Linares, and P. Kundur, \u00E2\u0080\u009CMultilevel MATE for Efficient Simultaneous Solution of Control Systems and Nonlinearities in the OVNI Simulator,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 21,110. 3, pp. 1250-1259, Aug. 2006. [66] P. Zhang, J. R. MartI and H. W. Dommel, \u00E2\u0080\u009CNetwork Partitioning for Real-time Power System Simulation,\u00E2\u0080\u009D IPST\u00E2\u0080\u009905, Montreal, Canada, Jun. 2005. [67] J. R. MartI, L. R. Linares, \u00E2\u0080\u009CReal-time EMTP-Based Transients Simulation,\u00E2\u0080\u009D IEEE Transactions on Power Systems, vol. 9, no. 3, pp. 1309-13 17, Aug. 1994. [68] T. De Rybel, J. Hoilman, J. R. MartI, \u00E2\u0080\u009COVNI-NET: a Flexible Cluster Interconnect for the New OVNI Real-Time Simulator,\u00E2\u0080\u009D 15th Power Systems Computation Conference, Liege, Belgium, Aug. 22-26, 2005. [69] F. A. Moreira, J. R. MartI, L. C. Zanetta, Jr., L. R. Linares, \u00E2\u0080\u009CMultirate Simulations With Simultaneous-Solution Using Direct Integration Methods in a Partitioned Network Environment,\u00E2\u0080\u009D IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 53, no. 12, pp. 2765-2778, Dec. 2006. [70] Note from H. W. Dommel, \u00E2\u0080\u0098Dynamic Phasors Versus dqo-Quantities,\u00E2\u0080\u0099 Oct. 25, 2008. 110 Appendix A Machine Parameters Appendix A Machine Parameters A.1 Synchronous Machine Parameters TABLE A-I SmcHRoNous MACHINE PARAM1TERs Parameters Hydro Turbine Generator Steam Turbine Generator r3 0.00234 2 0.00243 2 X, 0.1478(2 0.1538(2 X 0.5911(2 1.457(2 Xd 1.0467(2 1.457(2 rfd 0.0005 (2 0.00075 (2 X 0.2523(2 0.1145(2 r1 - 0.00144(2 Xlkql - 0.6578 (2 r12 0.01675 (2 0.0068 1 (2 X1 0.1267(2 0.07602(2 rij 0.01736(2 0.0108(2 X1 0.1970(2 0.06577(2 J 3.51x10Js2 6.58x104Ps2 Rating 325 MVA 835 MVA Line-to-line 20 KV 26 KV voltage Poles 64 2 111 Appendix A Machine Parameters A.2 Induction Machine Parameters TABLE A-Il INDUCTION MACHINE PARAMi\u00E2\u0080\u0099rERs 3-hp Induction Machine 0.435 2 0.754 2 2250-hp Induction Machine 0.0292 0.2262 Parameters r x\u00E2\u0080\u0099s XM 26.132 13.04 Xir 0.754 2 0.226 rr 0.816Q 0.022Q J 0.089kgm2 63.87 kgm2 TB 11.9Nm 8.9x10Nm Line-to-line 220 V 2.3 kV voltage Poles 4 4 112 Appendix B Voltage behind Reactance Induction Machine Model Appendix B Voltage behind Reactance Induction Machine Model This appendix faithfully reproduces reference [53], a project report for UBC course EECE 549. The only difference here is that all equation numbers and figure numbers are added with a \u00E2\u0080\u0098B.\u00E2\u0080\u0099. This report implemented the Voltage behind Reactance (VBR) model for the symmetrical induction machine simulations. The discrete-time VBR model has the similar form as the EMTP-type equivalent circuits. The simulation results validate the accuracy and efficiency of this new induction machine model. B. 1 Voltage Behind Reactance Model of Induction Machine The equations of the induction machine can be expressed in the arbitrary reference frame as Vqs rslqs + + P)Lq (B.1) VdS 1ds CO)Lqs +pAd$ (B.2) = + pA0, (B.3) o = rriqr + (o \u00E2\u0080\u0094 OJr )adr + Pqr (B.4) o = rridr ( \u00E2\u0080\u0094 0)r )Aqr + P21dr (B.5) o = rrior + pA,0,. (B.6) where each variable and parameter have been converted.to stator side, and co = 0 for stationary reference frame, co = co,. for rotor reference frame and co = coe for 113 Appendix B Voltage behind Reactance Induction Machine Model synchronous reference frame. Flux linkages equations are expressed as = Lislqs + Amq = Lis \u00E2\u0080\u0098qs + LM (\u00E2\u0080\u0098qs + qr) (B.7) = L1SIdS + Ad = L,SidS + LM (\u00E2\u0080\u0098di + ld) (B.8) A05 = L,5i05 (B.9) Aqr LirIqr + Amq (B. 10) Adr = Lirl& + A (B. 11) \u00E2\u0080\u00980r =L1i0 (B.12) Rearrange and manipulate (B.7) - (B.12), we can obtain A,q = L [qs + (B. 13) Amd1M[1ds +.z\u00E2\u0080\u0099,J (B.14) where LL =(+J\u00E2\u0080\u0099 Substitute (B.13), (B.14) into (B.7), (B.8), the q-axis and d-axis flux linkages can be rewritten as Aq5 = L\u00E2\u0080\u009D1q5 + (B. 15) A =L\u00E2\u0080\u009DidS +A\u00E2\u0080\u0099 (B.16) where L\u00E2\u0080\u009D\u00E2\u0080\u0094L -\u00E2\u0080\u0098-L\u00E2\u0080\u009D \u00E2\u0080\u0094 is M (B.17) (B.18) Therefore (B.1) and (B.2) can be written as Vqs = + o(L\u00E2\u0080\u009DidS + A) + P(1\u00E2\u0080\u009Dlqs + A:) (B.19) 114 Appendix B Voltage behind Reactance Induction Machine Model = rSidS \u00E2\u0080\u0094 (D(L\u00E2\u0080\u009Dqs + + p(L\u00E2\u0080\u009Di + (B.20) By manipulating (B.1)-(B.18), we can obtain = LZIqs )\u00E2\u0080\u0094(Co _Cor)2] (B.21) 1 =r[ +LZidS _ar)+(0_C0r)2qr] (B.22) Substitute (B.21), (B.22) into (B.19), (B.20), the voltage behind reactance form of the voltage equations can be expressed Vq = rsiqs + OJL\u00E2\u0080\u009DidS + P(L\u00E2\u0080\u009Dlqs) + (B.23) VdS = rSidS \u00E2\u0080\u0094 (DL\u00E2\u0080\u009Dlqs + P(L\u00E2\u0080\u009Dlqs) + v (B.24) where = + \u00E2\u0080\u0098\u00E2\u0080\u0098 \u00E2\u0080\u0094 lqg \u00E2\u0080\u0094 \u00E2\u0080\u0094(co \u00E2\u0080\u0094 (B.25) v=\u00E2\u0080\u0094coA:+(;\u00E2\u0080\u0094Adr)+i+(co\u00E2\u0080\u0094cor), (B.26) Applying the inverse transformation to (B.23), (B.24), we can get the phase domain equations V (t) = (t) + p[L\u00E2\u0080\u009Di, (t)] + v(t) (B.27) where v(t) = [K(o)]\u00E2\u0080\u0099 v 0 L\u00E2\u0080\u009D L +\u00C2\u00B1L\u00E2\u0080\u009D _LIs3M 3 3 LbCS = \u00E2\u0080\u0094 \u00E2\u0080\u0094- L1+\u00E2\u0080\u0094L L,+LZ 115 Based on the trapezoidal rule, equation (B.27) can be discretized and rearranged as \u00E2\u0080\u00987abcs (t) = + --L5)iabcs(t) + vbC$(t) + ek (t) (B.28) Ch (t) = \u00E2\u0080\u0094 L5} abcs (t \u00E2\u0080\u0094 At) + VbCS (t \u00E2\u0080\u0094 At) \u00E2\u0080\u0094 Vabcs (t \u00E2\u0080\u0094 At) (B.29) Rewrite (B.25) and (B.26), we can get the matrix form ofthesubtransient voltages = fa11 a12 (t)1 + Ea13 0 ji (t)1 (B.30) [v (t)J [a21 a22 ]Ldr (t)] [ 0 a3 ][i (t)] a12 =Wr(t)jj a22 = a11 a23 = a13 manipulating (B.4),( B.5), (B.10) and (B.11), the rotor flux linkage (B.31) Appendix B Voltage behind Reactance Induction Machine Model B. 2 Discrete Time VBR Model where where a =\u00E2\u0080\u0098\u00E2\u0080\u0098I\u00E2\u0080\u0094i L LL1r a21 = a12 a Lrr 13 L By rearranging and equations can be obtained pH (t) E11 b12 (t)1 + rb!3 0 jIqs (t)[A (t)J [b21 b22 ][A.. (t)] [0 b23 ][i (t) where b =--1-\u00E2\u0080\u00941 \u00E2\u0080\u0098 L1, b21 = \u00E2\u0080\u0094b12 b22 = b11 Discretizing (B.3 1) by using the Trapezoidal rule, then E \u00C2\u00B0\u00E2\u0080\u0098)1 r 2 \u00E2\u0080\u0094b11At \u00E2\u0080\u0094 b12 (t)At1 2 +b11At b12 (t \u00E2\u0080\u0094 At)At1E (t \u00E2\u0080\u0094 At) L\u00E2\u0080\u0099dr (t)J [\u00E2\u0080\u0094 b21 (t)At 2 \u00E2\u0080\u0094b22At J [b2, (t \u00E2\u0080\u0094 At)At 2 +b22At J[2ar (t \u00E2\u0080\u0094 At) + 2 \u00E2\u0080\u0094b11At \u00E2\u0080\u0094 b12 (t)At1\u00E2\u0080\u0099 Ebi3At 0 Jj\u00E2\u0080\u00A2\u00E2\u0080\u00A21qS (t)1 + 1q5 (t \u00E2\u0080\u0094 At) 2 \u00E2\u0080\u0094 b22 At] [ 0 b23At] (t)J [i (t \u00E2\u0080\u0094 At) 116 b12 = \u00E2\u0080\u0094[a(t) \u00E2\u0080\u0094 0)r (t)] b13 =--LZ b23 = b13 (B.32) Appendix B Voltage behind Reactance Induction Machine Model Substitute (B.32) into (B.30), we can get the discrete time equations r1 (i\u00E2\u0080\u0099)1[ j =K(t)[ j+K2(t)vd(t) ldS(t) where K (t) = Hii a12 2\u00E2\u0080\u0094b11At \u00E2\u0080\u0094 b12 (t)AtlEbi3A 0 1 + E\u00E2\u0080\u0099u13 [a21 a22 J[\u00E2\u0080\u0094 b2 (t)At 2 \u00E2\u0080\u0094b22At] L 0 b23 At] [ 0 a23 K (t) \u00E2\u0080\u0094 raii a12 ir 2 \u00E2\u0080\u0094b11At \u00E2\u0080\u0094 b12 (t)Atlr 2 +b11At b12 (t \u00E2\u0080\u0094 At)&1r (t \u00E2\u0080\u00942 \u00E2\u0080\u0094 [ a22 ]L b21 (t)At 2 \u00E2\u0080\u0094b22At J[b2,(t \u00E2\u0080\u0094 At)At 2 +b22At JL\u00E2\u0080\u0099dr (t \u00E2\u0080\u0094 At) + Eaii a12 j 2 \u00E2\u0080\u0094 b,1At \u00E2\u0080\u0094 b,2 (t)At1\u00E2\u0080\u0099 P13 At 0 jqs 0\u00E2\u0080\u0099 \u00E2\u0080\u0094 At)[a2, a22 ][\u00E2\u0080\u0094 b21 (t)At 2 \u00E2\u0080\u0094b22At] [ 0 b23At][idS (t \u00E2\u0080\u0094 At) Applying the inverse transformation to (B.33), we can transfer the qd0 variables to the phase domain VbCS (t) = K (t)Iabcs (t) + ChS (t) (B.34) where [K (t) 01 K(t) = [K5 (0(t))]-\u00E2\u0080\u0099 [ j[K. (0(t))] ,EK 0\u00E2\u0080\u0099) ehS(t)=[KS(0(t))] L Finally, by substituting (B.34) into (B.28), the equivalent circuit for the stator voltage equations can be written as \u00E2\u0080\u0098abcs (t) = G eqVabcs (t) \u00E2\u0080\u0094 h(t) (B.35) where Geq [r5 + + K(t)] h(t) =Geq[Ch(t) + ehS(t)] 117 Appendix B Voltage behind Reactance Induction Machine Model B. 3 Free Acceleration Simulation of a 3-hp Induction Machine The first test case is to simulate the start-up transients of a 3-hp induction machine. The parameters of the machine are P=4 V1 220V r 0.4352 X, 0.754 2 XM= 26.13 2 X1r O.754 rr = 0.816c2 J 0.089 kg.m2. The simulation results viewed from different reference frames are shown in Fig. B. 1- B.3. (1) Figure B.1: Stationary Reference Frame Because the initially position of the stationary qd0 reference frame is zero, thef = fqs. Therefore the Iqs and Vqs are the same as i and v, respectively. The rotor variables are varying at 60 HZ after referred to the stationary reference frame. (2) Figure B.2: Rotor Reference Frame At the very beginning, when the speed of the machine is slow, the variables in the rotor reference frame are similar to those in the stationary reference frame. Then the stator variables will vary at slip frequency when the rotor speed is increasing. After the speed reach the synchronous speed, all variables in rotor reference frames becomes constant. (3) Figure B.3: Synchronous Reference Frame At stall the variables in the synchronous reference are varying at 60 HZ. The stator and rotor variables will become constant when the rotor speed reaches synchronous speed. Note that the zero position of synchronous reference frame that we chose is zero, therefore the magnitude of Vqs is identical to this ofv and Vds is zero. B. 4 Dynamic Performance of Induction Machine during Mechanical Torque Changes In this case, the machine originally operates in the steady state. At t = 0.05s, a constant mechanical torque of 12 N.m2 is applied. Then the mechanical torque is reversed to -12N.m at t = 0.5s. The simulation time is ls. Figure B.4 illustrates the simulation results with the VBR model. 118 Appendix B Voltage behind Reactance Induction Machine Model 0\u00E2\u0080\u0099 V 200 .0 E 100 0 0 >0\u00E2\u0080\u00A2 S 0 0S -100 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Figure B.1 Free Acceleration Characteristics in Stationary Reference Frame 119 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Appendix B Voltage behind Reactance Induction Machine Model I i.J.. I .: 0 06 07 08 09 1 OC - I 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 30C I \\45 0.6 0.7 0.8 0.9 ,nn ______________________________________________________ I I I 2O0OjOi2Oi3 04 05 06 07 08 09 Figure B.2 Free Acceleration Characteristics in Rotor Reference Frame 120 Appendix B Voltage behind Reactance Induction Machine Model I I I I I I I 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 2C I C\u00E2\u0080\u0099 -4C I I I I I I I 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 -I (VI J I0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 50 0 -50 100 > -20 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 05 06 07 08 09 100O 1 0I2 0I3 0I4 0i5 08 0I7 08 09 300 I I I I I I I I I Figure B.3 Free Acceleration Characteristics in Synchronous Reference Frame 121 Appendix B Voltage behind Reactance Induction Machine Model 21 22 22.12.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3 j 500C -500C I I I I I I 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3 220C I I iGoc 500C \u00E2\u0080\u0094--- \u00E2\u0080\u0094 I I I I I I I 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 Figure BA Dynamic Performance of a 3-hps Induction Machine during Step Changes in Load Torque 122 Appendix C Test Case Data Appendix C Test Case Data %*********** ************ *** *** ********** ************ ******** *************** % Test Case 1 % AnElectromagnetic Transient Program with Dynamic Phasor Solution % % % % % % nBus=2; Version 1.0 = = Last revised: 24 Feb.2007 = Peng Zhang = = Copyright 2008 \u00C2\u00A9 Peng Zhang = % P1 Line data nLine=1; % Three phase line LineFromBus(1)=1; LineToBus(1)2; I Ri ohm I Xl ohm I Cl muF 1R21 X2lC2l.... LineParameters =[0.0868455 0.0298305 0.0288883 0.2025449 0.0847210 0.0719161 0.00274 -0.0007 -0.00034; 0.0298305 0.0887966 0.0298305 0.0847210 0.1961452 0.0847210 -0.0007 0.00296 -0.00071 0.0288883 0.0298305 0.0868455 0.0719161 0.0847210 0.2025449 -0.00034 -0.00071 0.00274 1*10.56 % Load data nLoadl;% Three phase load LoadBus(1)=1; rLoad(1,1)=(12.47*1000/sqrt(3))P2/(2000*1000); % phase A resistance xLoad(1,1)=(12.47*1000/sqrt(3))\u00E2\u0080\u00992/(0.328684105 17886306346562595373367*2000*1000); % phase A reactance rLoad(1,2)=(12.47*1000/sqrt(3)y2/(2000* 1000); % phase B resistance xLoad(1,2)=(12.47*1000/sqrt(3)y2/(0.32868410517886306346562595373367*2000*1000); % phase B reactance rLoad(1,3)(12.47*1000/sqrt(3))t2/ 20 *1000);% phase C resistance xLoad(1,3)(12.47* 1000/sqrt(3)y2/(0.32868410517886306346562595373367*2000*1000);% phase C reactance % Voltage source data nVSourcel; VSourceBus(1)2; Vabcs(1,1) 10182;%Vamp Vabcs(1,2) = 10182 * exp(-jay * 2 * p / 3); 123 Appendix C Test Case Data Vabcs(1,3) 10182 * exp(jay* 2 * pi/3); T10= 3; LineTripped = 1; return; %************************************************************************** % Test Case 2 % An Electromagnetic Transient Program with Dynamic Phasor Solution % % % % % % nBusl2; Version 1.0 = = Last revised: 24 Feb.2007 = Peng Zhang = = Copyright 2008 \u00C2\u00A9 Peng Zhang = simT = 8; % total simulation time ja\u00E2\u0080\u0099sqrt(-l); ws2*pi*60; % PT Line data nLine=1 1; % Three phase line LineFromBus(l)12; LineToBus(1)=1; LineFromBus(2)1; LineToBus(2)4; LineFromBus(3)4; LineToBus(3)5; LineFromBus(4)=1; LineToBus(4)=2; LineFromBus(5)2; LineToBus(5)=3; LineFromBus(6)=l; LineToBus(6)6; LineFromBus(7)6; LineToBus(7)8; LineFromBus(8)8; LineToBus(8)9; LineFromBus(9)8; LineToBus(9)l0; LineFromBus(10)=6; LineToBus(10)1 1; LineFromBus(1 1)6; 124 Appendix C Test Case Data LineToBus(1 1)7; I Ri I Xl I ClmuF 1R21X21C21.... LineParameters(1:3,1:9)=[0.3465 0.0 0.0 1.0179 0.5017 0.4236 6.2998/(2*pi*60) _1.2958/(2*pi*60) 1.2595/(2*pi*60); 0.0 0.3375 0.0 0.5017 1.0478 0.3849 1.2958/(2*pi*60) 5.9597/(2*pi*60) 1.2595/(2*pi*60); 0.0 0.0 0.3414 0.4236 0.3849 1.0348 1.2595/(2*pi*60) 1.2595/(2*pi*60) 5.6386/(2*pi*60)] LineParameters(1:3,10:18) =[1.3238 0.0 0.0 1.3569 0.4591 0.4591 4.6658/(2*pi*60) 0.8999/(2*pi*60) 0.8999I(2*pi*60); 0.0 1.3294 0.0 0.4591 1.3471 0.4591 0.8999/(2*pi*60) 4.7097/(2*pi*60) 0.8999/(2*pi*60); 0.0 0.0 1.3238 0.4591 0.4591 1.3569 .0.8999/(2*pi*60) .0.8999/(2*pi*60) 4.6658/(2*pi*60)] LineParameters(1:3,19:27) =[1.3238 0.0 0.0 1.3569 0.4591 0.4591 4.66581(2*pi*60) 0.8999/(2*pi*60) 0.8999/(2*pi*60); 0.0 1.3294 0.0 0.4591 1.3471 0.4591 0.8999/(2*pi*60) 4.7097/(2*pi*60) 0.8999/(2*pi*60); 0.0 0.0 1.3238 0.4591 0.4591 1.3569 0.8999/(2*pi*60) 0.8999I(2*pi*60) 4.6658/(2*pi*60)] LineParameters(1:3,28:36) [0.7526 0.0 0.0 1.1814 0.4236 0.5017 5.6990/(2*pi*60) 1.0817/(2*pi*60) 1.6905/(2*pi*60); 0.0 0.7475 0.0 0.4236 1.1983 0.3849 1.0817/(2*pi*60) 5.1795/(2*pi*60) 0.6588/(2*pi*60); 0.0 0.0 0.7436 0.5017 0.3849 1.2112 1.6905I(2*pi*60) 0.6588/(2*pi*60) 5.4246/(2*pi*60)] LineParameters(1:3,37:45) =[O.3807 0.0000 0.0000 0.6922 0.0000 0.0000 0.0001 0.0000 0.0000; 0.0000 0.3807 0.0000 0.0000 0.6922 0.0000 0.0000 0.000 1 0.0000 0.0000 0.0000 0.3807 0.0000 0.0000 0.6922 0.0000 0.0000 0.0001 LineParameters(1:3,46:54) [ 0.3465 0.0 0.0 1.0179 0.5017 0.4236 6.2998/(2*pi*60) 1.2595/(2*pi*60) 1.2595I(2*pi*60); 0.0 0.3375 0.0 0.5017 1.0478 0.3849 1.2595I(2*pi*60) 5.9597/(2*pi*60) 1.2595/(2*pi*60); 0.0 0.0 0.3414 0.4236 0.3849 1.0348 1.2595/(2*pi*60) 1.2595/(2*pi*60) 5.6386/(2*pi*60)] LineParameters(1:3,55:63) =[1.3238 0.0 0.0 1.3569 0.4591 0.4591 4.6658/(2*pi*60) 0.8999/(2*pi*60) 0.8999/(2*pi*60); 0.0 1.3238 0.0 0.4591 0.4591 0.4591 0.8999/(2*pi*60) 4.70971(2*pi*60) 0.8999/(2*pi*60); 0.0 0.0 1.3294 0.4591 0.4591 1.3471 0.8999I(2*pi*60) 0.8999I(2*pi*60) 4.7097/(2*pi*60)] LineParameters(1 :3,64:72) =[1.3292 0.0000 0.0000 1.3475 0.0000 0.0000 4.5 193/(2*pi*60) 0.0000 0.0000; 0.0000 1.3292 0.0000 0.0000 1.3475 0.0000 0.0000 4.5193/(2*pi*60) 0.0000; 0.0000 0.0000 1.3292 0.0000 0.0000 1.3475 0.0000 0.0000 4.5193/(2*pi*60)] \ LineParameters(1:3,73:81)=[1.3425 0.0000 0.0000 0.5124 0.0000 0.0000 88.9912/(2*pi*60) 0.0000 0.0000; 0.0000 1.3425 0.0000 0.0000 0.5124 0.0000 0.0000 88.9912/(2*pi*60) 0.0000; 0.0000 0.0000 1.3425 0.0000 0.0000 0.5124 0.0000 0.0000 88.9912/(2*pi*60)] Lineparameters(1:3,82:90) =[ 0.3465 0.0 0.0 1.0179 0.5017 0.4236 6.2998/(2*pi*60) .1.9958/(2*pi*60) ..1.2595/(2*pi*60); 0.0 0.3375 0.0 0.5017 1.0478 0.3849 1.9958/(2*pi*60) 5.9597/(2*pi*60) 0.7417/(2*pi*60); 0.0 0.0 0.3414 0.4236 0.3849 1.0348 -1 .2595/(2*pi*60) .0.7417/(2*pi*60) 5.6386/(2*pi*60)] LineParameters(1:3,91:99) =[0.7982 0.0 0.0 0.4463 0.0328 -0.0143 96.8897/(2*pi*60) 0.0000 0.0000; 0.0 0.7891 0.0 0.0328 0.4041 0.0328 0.0000 96.88971(2*pi*60) 0.0000; 0.0 0.0 0.7982 -0.0143 0.0328 0.4463 0.0000 0.0000 96.8897/(2*pi*60)] % Load data nLoad=8;% Three phase load LoadBus(1)\u00E2\u0080\u00993; LoadBus(2)=4; LoadBus(3)5; LoadBus(4)=10; LoadBus(5)6; LoadBus(6)=7; LoadBus(7)1 1; LoadBus(8)=9; rLoad\u00E2\u0080\u0094zeros(8,3); 125 Appendix C Test Case Data rLoad(:,:)1e6; xLoad=zeros(8,3); xLoad(:,:)=1e6; rLoad(5,1)=(4. 16* 1000/sqrt(3))s2/(385* 1000); % phase A resistance xLoad(5,1)=(4. 16*1000/sqrt(3))2/(220* 1000); % phase A reactance rLoad(5,2)=(4. 16* l000/sqrt(3))/s2/(385* 1000); % phase B resistance xLoad(5,2)=(4.16*1000/sqrt(3))f2/(220*1000); % phase B reactance rLoad(5,3)=(4.16* l000/sqrt(3))f2/(385* 1000); % phase C resistance xLoad(5,3)=(4. 16*100O/sqrt(3))l2/(220*1000); % phase C reactance % Induction machine data IndMacBus=1 1; Ttorq=3; % torque applied I S I VI I P J rs I Lls I LM I iT I Llr IndMacParameters = [500 /(1.341e-3) 4600 4 11.06 0.262*4 1 .206*4/ws 54.02*4/ws 0.187*4 1 .206*4Iws]; % Voltage source data nVSourcel; VSourceBus(1)=12; Vabcs(1,1) = (4.16*1000/sqrt(3))*sqrt(2) ; % Vamp Vabcs(1,2) = (4.16*l000Isqrt(3))*sqrt(2) * exp(-jay * 2 * pj / 3); Vabcs(1,3) = (4.16*l000Isqrt(3))*sqrt(2) * exp(jay * 2 * pi /3); return; 126"@en .
"Thesis/Dissertation"@en .
"2009-05"@en .
"10.14288/1.0067234"@en .
"eng"@en .
"Electrical and Computer Engineering"@en .
"Vancouver : University of British Columbia Library"@en .
"University of British Columbia"@en .
"Attribution-NonCommercial-NoDerivatives 4.0 International"@en .
"http://creativecommons.org/licenses/by-nc-nd/4.0/"@en .
"Graduate"@en .
"Shifted frequency analysis for EMTP simulation of power system dynamics"@en .
"Text"@en .
"http://hdl.handle.net/2429/7727"@en .