UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Chaotic ferroresonance in power transformers Mozaffari, Said 1996

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

Item Metadata

Download

Media
831-ubc_1997-196283.pdf [ 5.06MB ]
Metadata
JSON: 831-1.0064846.json
JSON-LD: 831-1.0064846-ld.json
RDF/XML (Pretty): 831-1.0064846-rdf.xml
RDF/JSON: 831-1.0064846-rdf.json
Turtle: 831-1.0064846-turtle.txt
N-Triples: 831-1.0064846-rdf-ntriples.txt
Original Record: 831-1.0064846-source.json
Full Text
831-1.0064846-fulltext.txt
Citation
831-1.0064846.ris

Full Text

C H A O T I C F E R R O R E S O N A N C E IN P O W E R T R A N S F O R M E R S  By SAID MOZAFFARI B.Sc, Civil Eng., United States International University, 1987 B.A.Sc, Elec. Eng., University of British Columbia, 1990 M.A.Sc, Elec. Eng., University of British Columbia, 1993 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF T H E REQUIREMENTS FOR T H E DEGREE OF D O C T O R OF PHILOSOPHY in T H E F A C U L T Y OF G R A D U A T E STUDIES DEPARTMENT OF ELECTRICAL ENGINEERING We accept this thesis as conforming t(^he^eqyire^l standard  T H E UNIVERSITY OF BRITISH COLUMBIA August 1996 © Said Mozaffari, 1996  In presenting this thesis in partial fulfilment  of the requirements for an advanced  degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department  or  by his or  her  representatives.  It  is understood that  copying or  publication of this thesis for financial gain shall not be allowed without my written permission.  Department The University of British Columbia Vancouver, Canada Date  DE-6 (2/88)  ABSTRACT  Ferroresonance can be categorized as a nonlinear resonance which is unpredictable and whose occurrence can cause damage in power distribution and transmission systems. In most instances, ferroresonance occurs when one or two of the source phases are lost while the transformer is unloaded or lightly loaded. This can happen due to the clearing of single phase fusing, the operation of single phase reclosers or the performance of single phase switching procedures. In these cases, if one of the three switches was open, and only two phases of the transformer were energized, there would be a voltage induced i n the "open" phase. This induced voltage will "backfeed" the distribution line, back to the open switch. Ferroresonance will occur if the distribution line is highly capacitive. This ferroresonance will involve the nonlinear magnetizing reactance of the transformer's open phase and the shunt capacitance of the distribution line.  O n e of the fundamental properties of ferroresonance is the fact that several stable solutions can  exist under  steady-state conditions  for  a given circuit.  The  realization  that  ferroresonance is a nonlinear and sometimes chaotic process opens up many possibilities. The newly developed techniques for analysis of nonlinear dynamical systems and chaos should now be evaluated for use with ferroresonance. T h e method of Slowly Varying Amplitude is used to derive an analytical solution for the equivalent ferroresonant circuit. The solution of the nonlinear equation for a typical ferroresonant circuit containing a power transformer  is shown to  be dependent ii  on the  accurate description of  the  magnetization curve. A detailed analysis of many simulation results demonstrates that the probability of chaos increases as losses decrease and the nonlinearity of the transformer magnetization rises. T h e effect of varying the transformer core losses, the value of the source voltage and the length of the transmission line on the chaotic solution of the system has been studied. T h e concept of transient chaos as compared with steady-state chaos is also discussed. T h e solution of the ferroresonant circuit is shown to be dependent o n the value of its initial conditions. It is shown that a small change in the initial conditions leads to a large difference in long-term behavior of the system, and this makes the future of the system unpredictable. W i t h detailed analysis of many simulation results, the basins of attraction for different chaotic regions of the system have been obtained. It is shown that inclusion of series losses is not an important factor in studying ferroresonance as compared to core losses. T h e Electromagnetic Transient Program ( E M T P ) is used to simulate ferroresonance and to obtain Poincare maps and bifurcation diagrams. T h e appropriate representation of nonlinear elements in the E M T P  when chaotic systems are to be  simulated is shown to be very important. A l l the simulation results are tested using three different integration routines. In all the simulations the results of using different integration routines were the same. T h e result of the case study performed by using the E M T P showed the existence of chaos and bifurcation in a simple power system and the need to develop accurate ways to reliably identify chaotic behavior.  111  TABLE OF CONTENTS  Abstract Table of Contents List of Figures List of Tables Acknowledgments  .  1. I N T R O D U C T I O N 2. T H E O R E T I C A L CONSIDERATIONS 2.1 Introduction  ...  2.2 Ferroresonance in Power Systems 2.3 Brief History of Ferroresonance 2.4 A Brief History of Nonlinear Dynamical Systems and Chaos 2.5 Definition and Categorization Methods for Chaotic systems. 2.6 Nonlinear Dynamics and Chaos Applied to Ferroresonance..  3. SYSTEM DESCRIPTION A N D M O D E L I N G 3.1 Introduction 3.2 T h e Equivalent Circuit 3.3 Differential Equation of the System 3.4 Determination of Parameters and Per U n i t System  4. A N A L Y T I C A L SOLUTIONS 4.1 Introduction 4.2 Analysis of Nonlinear Equation  5. SIMULATION A N D C A T E G O R I Z A T I O N O F FERRORESONANCE 5.1 Introduction 5.2 Varying the Magnitude of the Source Voltage 5.3 T h e Effect of the Transformer Core Losses 5.4 Transient Chaos 5.5 Varying the Length of the Transmission Line 5.6 Dependence on the Description of the Magnetization C u r v e . . 5.7 T h e Effect of Initial Conditions 5.8 T h e Effect of Copper or Series Losses 5.9 T h e Effect of Higher Order Harmonics  iv  6. USING EMTP F O R SIMULATION O F C H A O S IN POWER SYSTEMS  101  6.1 Introduction  101  6.2 A n Overview of E M T P  102  6.3 Ferroresonance Test Case  107  6.4 Representation of the Transformer Nonlinear Characteristic with E M T P  Ill  6.5 Integration of the Trajectories  113  6.6 Case Study  116  7. CONCLUSIONS A N D RECOMMENDATIONS  127  7.1 Conclusions  127  7.2 Recommendations  133  REFERENCES  135  APPENDIX A  .....  v  138  LIST O F F I G U R E S  Figure 2.1:  Series L C circuit with sinusoidal voltage source  Figure 2.2:  X -1  Figure 2.3:  Frequency response of a linear L C circuit  Figure 2.4:  X-I  characteristic for a nonlinear inductance  7  Figure 2.5:  Frequency response of a nonlinear L C circuit  8  Figure 2.6:  Jump phenomena for variation of the amplitude of the excitation  9  Figure 2.7:  Ferroresonant circuit  9  characteristic for a linear inductor  6 6 7  Figure 2.8:  Graphical solution of the ferroresonance circuit  12  Figure 2.9a:  Phase plane trajectory of Duffing oscillator for period-1 motion  20  Figure 2.9b: Phase plane trajectory of Duffing oscillator for period-3 motion  20  Figure 2.9c:  21  Phase plane trajectory of Duffing oscillator for period-5 motion  Figure 2.9d: Phase plane trajectory of Duffing oscillator for chaotic motion  21  Figure 2.10:  22  Poincare mapping  Figure 2.11a: Poincare map of Duffing oscillator for period-1 motion  23  Figure 2.11b: Poincare map of Duffing oscillator for period-3 motion  24  Figure 2.11c: Poincare map of Duffing oscillator for period-5 motion  24  Figure 2 . l i d : Poincare map of Duffing oscillator for chaotic motion  25  Figure 3.1:  31  A typical distribution system supplying a three phase load  Figure 3.2:  Model for ferroresonant circuit including line capacitance  32  Figure 3.3:  T h e circuit that feeds the disconnected coil  32  Figure 3.4:  Preliminary equivalent circuit  34  Figure 3.5:  Reduced equivalent circuit  35  Figure 3.6:  Approximation of the saturation region of the test transformer  37  Figure 4.1:  Equivalent ferroresonant circuit  41  Figure 4.2:  Graph of P verses A  46  Figure 4.3:  Jump phenomenon  47  vi  Figure 4.4:  Effect of clamping on the frequency response curve  48  Figure 4.5:  T h e graph of the amplitude of the flux linkage vs. eigenvalues  50  Figure 4.6:  Amplitude of the response as a function of amplitude of the excitation  53  Figure 5.1:  Bifurcation diagram with varying E  56  Figure 5.2a: Phase plane diagram of period one region for £ = 0 . 1 5 to 1.76 pu., R = l%,Cfor  100 k m  57  Figure 5.2b: Poincare map of period-1 region for £ = 0 . 1 5 to 1.76 pu., R = 1%, Cfor 100 k m  58  Figure 5.2c: Frequency power spectrum for £ = 0 . 1 5 to 1.76 pu., R = 1%, Cfor 100 k m  58  Figure 5.3a: Poincare map of period-2 region for £ = 1.77 pu  59  Figure 5.3b: Phase plane diagram of period-2 region for E= 1.77 pu  59  Figure 5.4a: Poincare map of period-4 for E= 1.87 pu  60  Figure 5.4b: Phase plane diagram of period-4 region for E= 1.87 pu.  60  Figure 5.5:  Poincare map of chaotic solution for E= 1.9 pu  61  Figure 5.6:  Frequency power spectrum for E= 1.9 pu  62  Figure 5.7a: Poincare map of period-2 region for R = 0.005% Figure 5.7b: Phase plane diagram of period-2 region for R = 0.005%  64 - 65  Figure 5.8a: Poincare map of period-4 region for R = 0.001%  65  Figure 5.8b: Phase plane diagram of period-4 region for R = 0.001%  66  Figure 5.9:  67  Poincare map of chaotic solution for R = 96.8 MQ  Figure 5.10: Poincare map of chaotic solution for conservative system  70  Figure 5.11: Poincare map of the quasi-periodic solution for the conservative system  72  Figure 5.12: Poincare map of the chaotic solution for the dissipative system  72  Figure 5.13a: Poincare map after 3 sec. simulation time  74  Figure 5.13b: Poincare map after 6 sec. simulation time  75  Figure 5.13c: Poincare map after 12 sec. simulation time  75  Figure 5.13d: Poincare map after 30 sec. simulation time  76  Figure 5.14a: Poincare map of transient chaotic solution  76  vii  Figure 5.14b: Poincare map after the transient solution is removed Figure 5.15:  77  Bifurcation map of varying the length of the line in k m . , £ = 1.7 pu.,  R-1%  78  Figure 5.16: Bifurcation map of varying the length of the line in k m . , E= 1.85 pu.,  R-1%  79  Figure 5.17: Bifurcation map of varying the length of the line in k m . , £ = 2 . 5 pu., # = 1%  80  Figure 5.18: Bifurcation map of varying the length of the line in k m . , E= 1.0 pu.,  R-0.1%  81  Figure 5.19: Bifurcation map of varying the length of the line in k m . , E= 1.6 pu.,  R=0.1%  82  Figure 5.20: Bifurcation map of varying the length of the line i n k m . , £ = 2 . 0 pu.,  R-0.1%  82  Figure 5.21:  Current linkage with varying exponents  84  Figure 5.22:  Bifurcation diagram with varying E for n = 5  85  Figure 5.23:  Poincare map of transient chaotic solution for « = 5, i? = 96.8  Figure 5.24:  Poincare map after the transient removed  87  Figure 5.25:  Example of attractors and basin of attractions in state space  88  Figure 5.26:  Poincare map of periodic solution R = 0.0005%, £ = 0 . 1 5 pu., Cfor 100km.... 90  Figure 5.27:  Basin of attraction for periodic attractor R=1%,  £ = 0 . 1 5 pu., Cfor 100km.... 91  Figure 5.28:  Basin of attraction for periodic attractor R=1%,  £ = 1.9 pu., Cfor 100km  Figure 5.29:  Poincare map of chaotic solution for R = 1%, £ = 5 pu., C for 100km  92  Figure 5.30:  Basin of attraction for periodic attractor i?= 1%, £ = 5 pu., C for 100km  93  MQ.  86  92  Figure 5.31: Poincare map of chaotic solution for R = 1%, £ = 2 5 pu., i?= 1%,  94  Figure 5.32: Bifurcation diagram with varying £ and R =6%  97  Figure 5.33: Bifurcation diagram with added 50% third harmonic  99  i  Figure 6.1:  Piece-wise linear inductance representation  103  Figure 6.2:  Obtaining flux in the E M T P  105  Figure 6.3:  Poincare map of the E M T P simulations for a period-1 response  108  Figure 6.4:  Poincare map of the E M T P simulations for a period-2 response  109  Vlll  Figure 6.5: Poincare map of the EMTP simulations for the chaotic response with £=1.9 pu  109  Figure 6.6: Poincare map of the EMTP simulations for chaotic response for low losses (i?=96.8 MQ)  110  Figure 6.7: Poincare map of the EMTP simulations for conservative system  Ill  Figure 6.8: Discretization of the saturation curve with 11 and 21 points  112  Figure 6.9: Poincare map of the region with £ = 1.9 pu. for 11-point representation  113  Figure 6.10: Poincare map of the simulation with Gear algorithm £ = 1.9 pu., R=\%,  Cfor  100km  116  Figure 6.11: Power system under analysis  117  Figure 6.12: Voltage across the A l - B l for the first 0.1 seconds  118  Figure 6.13: Voltage waveform across the A l - B l for 1 seconds  119  Figure 6.14: Frequency power spectrum for the voltage A l - B l  119  Figure 6.15: Phase plot diagram of the flux vs. the voltage across the node A l - B l  121  Figure 6.16: Poincare map of the flux vs. the voltage across the node A l - B l  121  Figure 6.17: Bifurcation diagram for varying the source voltage  123  Figure 6.18: Poincare map for period-1 response V  =7.0 kV.  124  =9.5 kV  124  Figure 6.20: Poincare map for chaotic response F , =5.5 kV  125  source  Figure 6.19: Poincare map for period-2 response V  sourci!  sol ree  ix  LIST O F T A B L E S  Table 5.1: Varying the magnitude of the source voltage (R= 1%, C for 100km)  63  Table 5.2: Varying the value of the transformer core losses ( £ = 0 . 1 5 p u . Cfor 100km)  67  Table 5.3: Simultaneous variation of core losses and source voltage C for 100km  68  x  ACKNOWLEDGMENTS  I would like to thank m y wife Suzanne Belanger for her moral support and constant encouragement.  I would also like to thank m y parents for their love, nurturing and  support.  I am indebted to m y supervisor, D r .  A . C . Soudack, for his constant guidance and his  indispensable role in the successful completion of this thesis.  I would also like to thank m y good friends and fellow students in the Power G r o u p of the Electrical Engineering Department at the University of British Columbia, for their dedication to providing a friendly and cooperative environment for research.  I would like to thank M r . James A . Halladay of the Electrical Research Section of American Electric Power for supplying the case study discussed in Chapter 6 of this thesis.  xi  Chapter  1  INTRODUCTION  Ferroresonance is a type of nonlinear resonance characterized by overvoltages whose waveforms are highly distorted and which can cause equipment damage and power quality problems [1].  In simple terms, ferroresonance is an L C "resonance" involving a nonlinear  inductance and a capacitance. F o r ferroresonance to occur in power networks, the system must have a voltage supply (generally sinusoidal), capacitance (due to lines), a nonlinear transformer or an inductance coil (including saturable ferromagnetic materials), and low losses. Because of the nonlinearity in these systems, more than one mode of ferroresonance can  occur. Ferroresonant voltages depend on the magnitudes of the source voltage,  capacitance, and losses. A t the same time, the nonlinear inductance characteristic which is described by the transformer excitation curve can also effect the ferroresonant voltage. In most instances, ferroresonance occurs when one or two of the source phases are lost while the transformer is unloaded or lightly loaded. This can happen due to the clearing of single phase fusing, the operation of single phase reclosers or the performance of single phase switching procedures. O n e of the fundamental properties of ferroresonance is the fact that several stable solutions can exist under steady-state conditions for a given circuit. F o r example, remnant magnetization of the core, voltage at the time of switching and the amount of charge on the capacitance are all initial conditions which determine the steady state response. W i t h a small difference in these initial conditions, the ferroresonance  1  overvoltages can have very different waveforms.  The purpose of this Ph.D. thesis is to  investigate ferroresonance from the standpoint of nonlinear dynamics and chaotic systems.  In nonlinear dynamical systems, due to the nonlinearities inherent i n the system, the behavior may not be predictable. Depending on the magnitude of the forcing function and on the system parameters, the system's output waveform may be either periodic or nonperiodic (chaotic). This thesis uses the theories of nonlinear dynamics and chaos to analyze, understand, and control ferroresonance. W o r k in this area has been done by others to categorize the various modes of ferroresonance in commonly used sizes of five-legged core transformers [2]. T h e object of this work was to compare the qualitative characteristics of voltage waveforms in the time domain. T h e realization that ferroresonance is a nonlinear and sometimes chaotic process opens up many possibilities. It can be shown that the transitions between periodic and non-periodic modes commonly occur due to small changes  in circuit parameters or  initial conditions. Therefore, the newly-developed  techniques for analysis of nonlinear dynamical systems and chaos should now be evaluated for use with ferroresonance.  In general, power systems are characterized by nonlinear differential equations, and unusual and unexpected behavior has been observed in both simple and complex networks. T h e behavior of power systems are expected to become even more complicated because of such factors as greater interconnection and the implementation of high-speed electronic control devices.  Another  power  main objective of this thesis is to conduct an exploratory  2  research to determine whether the potential for chaos exists in power networks and, if it exists, how to avoid or control the chaotic behavior. This kind of research can increase our fundamental understanding of power system dynamics and reduce the potential for catastrophic  failure as well as allow secure network operation closer to  performance  margins.  In Chapter 2, the ferroresonance phenomenon in power systems is explained along with a short historical review of works in this area. investigate ferroresonance  The method of graphical solution used to  is presented. This chapter also provides an introduction to  nonlinear dynamical systems and the concept of chaos. A description of all the methods of observation, categorization and analysis of chaotic systems is presented. T h e importance of applying the methods of nonlinear dynamics and chaos to the problem of ferroresonance is discussed.  In Chapter 3, the equivalent circuit of a typical ferroresonant circuit in a power system along with the set of differential equations that describes it is developed. A l l the parameters of this model and their per unit values are presented.  Chapter 4, describes the  method of slowly varying amplitude [3] and its application in  deriving an analytical solution for the equivalent ferroresonant circuit. Using this method, it will be shown that the model of the ferroresonant  circuit exhibits a saddle-node  bifurcation leading to the occurrence of two stable attracting points and one region of an unstable saddle point which leads to a chaotic response.  3  Chapter 5, presents an in-depth analysis of the behavior of the ferroresonant circuit. A wide range of solutions for different modes of behavior for the ferroresonant circuit is presented. The critical values of the circuit parameters that can lead the ferroresonant circuit into chaos are shown. T h e concept of transient chaos when the transformer core losses are small is  discussed. T h e dependence  of the  solution on the  accurate  description of  the  magnetization curve is addressed. T h e basin of attraction for different chaotic regions of the system  is obtained. The effect  of series  losses  on the  modes  of behavior of  the  ferroresonance circuit is presented.  Chapter 6, describes how the E M T P (Electromagnetic Transient Program) can be used for simulations of chaos in power systems. It explains  a method of using E M T P to obtain  Poincare maps and to perform bifurcation simulations for the study of chaos.  EMTP  features used i n the modeling of the ferroresonant circuit are outlined. T h e influence of the nonlinear modeling i n the E M T P on the simulation results will also be discussed in this chapter. This chapter also presents the results of  a case study which was performed by  using the E M T P . This study is an example of one of the main motivations behind this work, which is to show a systematic way of analyzing system stability based on bifurcation theory and observed sustained oscillations and chaos via the period-doubling route. This case study validates the existence of chaos and bifurcations in power systems and the need to develop accurate ways to reliably identify their presence.  Chapter 7 summarizes the contributions and the conclusions of this work. Based on the results, recommendations and suggestions for further work are made.  4  Chapter  THEORETICAL  2  CONSIDERATIONS  2.1 Introduction This chapter introduces ferroresonance in power systems, and describes the methods that have been used to investigate this phenomena. A n overview of the related literature and contributions to the area of ferroresonance is also provided. A brief history of nonlinear dynamical systems and chaos is presented, and a description of some of the techniques that indicate the existence of chaos will be shown. T h e methods of observation, categorization and analysis of various types of chaotic behavior will be presented, and the use of these techniques in analyzing ferroresonance will be discussed.  2.2 Ferroresonance in Power Systems The trend toward using higher distribution voltages and underground feeders has increased the number of instances in which ferroresonance overvoltages have been reported [4]. T h e problem of ferroresonance can be categorized as a nonlinear resonance which can cause damage in power distribution and transmission systems. In simple terms, ferroresonance is an L C resonance involving a nonlinear inductance and a capacitance. Ferroresonance can be better understood if it is compared with linear L C resonance. A typical linear L C resonant circuit consists of an ideal inductor connected in parallel or series with an ideal capacitor. There is no damping in the circuit and the behavior of this L C combination is observed as  5  the frequency of an applied sinusoidal voltage or current is varied. Figure 2.1 shows an example of a series L C circuit with linear elements. In this circuit L is constant and independent of current, regardless of the flux linked X by the inductor. T h e relationship between the flux and the current is shown i n the Figure 2.2.  F i g u r e 2.1: Series L C C i r c u i t w i t h S i n u s o i d a l V o l t a g e Source  X  F i g u r e 2.2: X-l C h a r a c t e r i s t i c f o r a linear i n d u c t o r  In this circuit, the resonance occurs when the total impedance of the circuit jX, - jX  c  equal to zero. T h e frequency at which this happens is co,. =  ,  is  . Therefore as co  approaches co , current I approaches oo. Figure 2.3 shows the frequency response of a linear r  L C circuit.  /  F i g u r e 2.3: F r e q u e n c y Response o f a L i n e a r L C C i r c u i t  If damping is considered as a resistance in parallel with the inductor in the above case, the current will be limited to I< oo and the resonant frequency will be shifted to a new damped resonant frequency co = co^ . In the circuit of Figure 2.1, if the linear inductance is replaced by a nonlinear saturable iron core inductor then ferroresonance can occur. A typical A.-I characteristic of a such an inductor is shown in Figure 2.4,  which is a typical for a  transformer core.  F i g u r e 2.4: X-l C h a r a c t e r i s t i c f o r a n o n l i n e a r i n d u c t a n c e  The circuit with a nonlinear inductance will have a much different type of frequency response. In this case, there is no single resonant frequency, since the frequency response  7  characteristic will be multi-valued [5]. Figure 2.5  shows the frequency response of a  nonlinear L C circuit when there is no damping in the system.  co  CO  c  F i g u r e 2.5: F r e q u e n c y Response of a N o n l i n e a r L C C i r c u i t  F o r co > co , there may be two stable modes of "resonance" as well as a third unstable c  mode. F o r co <  co , there is only one possible mode of operation. A s co is decreased and c  its value passes co , operation can make a sudden change or jump from one stable operating c  mode to another. This is one of the reasons why this type of behavior is sometimes called jump resonance. Therefore, jump resonance refers to a condition in a sinusoidally excited system where an incremental change in the frequency of the input to the system  or an  incremental change in the driving voltage or in the magnitude of one of the parameters of the system causes a sudden jump in signal amplitude somewhere in the system. This jump can be one of voltage, current, flux linkage or all three. Figure 2.6 shows the jump phenomenon when the amplitude of the excitation is varied slowly. In this diagram the effect of the damping is considered. In this figure, starting from point 1, as V is increased, s  V, slowly increases through point 2 to point 3. A s V is increased further, a jump takes place s  8  from point 3 to point 4  with an accompanying increase in V ! , after which V[ increases  slowly with V . If the process is reversed, V) decreases slowly as V s  s  decreases from point 5  to point 6. A s V is decreased further, a jump from point 6 to point 2 takes place, with an s  accompanying Rudenberg [5]  decrease in V  (  , after which V , decreases slowly with decreasing V . s  gives a very clear explanation of this jump phenomenon based on a  graphical method.  F i g u r e 2.6: J u m p p h e n o m e n a f o r v a r i a t i o n of the a m p l i t u d e of the e x c i t a t i o n  Consider the circuit of Figure 2.7, in which the linear inductance has been replaced with a nonlinear inductor.  F i g u r e 2.7: F e r r o r e s o n a n t C i r c u i t  9  W h e n the series resistance is ignored, the sum of the voltages around the only mesh of the circuit can be written as:  ^(0-^ (0-^/(0 = o  (2.1)  e  The value of V can be replaced by its time-integral expression and V, as the total derivative c  of  L(i)i(t).  T h e n equation 2.1 can be written as:  V (t) ~ ^ \i{t)dt- j\[L(i)i(t)] = 0 s  Evaluating equation 2.2 and substituting  q(t) - ^i(t)dt  (2.2)  will result in  ^ ( 0 - ^ ( 0 - 1 ( 0 - ^ — ^ - ^ =0  (2.3)  F r o m equation 2.3 it is evident that finding a closed-form solution for this nonlinear circuit will be quite difficult. This would be made more evident by adding a source impedance and by providing a complete equivalent of the transformer.  Historically, methods of graphical solution represent one of the earliest attempts to explain ferroresonance. T h e graphical solution for the circuit of Figure 2.7, including the series resistance, can be obtained from two independent relationships for the voltage across the inductance and the capacitance, (assuming sinusoidal variation of current). T h e voltage  10  across the inductor is proportional to the frequency, and the voltage across the capacitance is proportional to the current and inversely proportional to the frequency and capacitance.  V, = co  /(/)  co C  The total magnitude of voltage for the circuit is:  V, = J(V,+V ) +(RI) C  2  (2.5)  2  F r o m equations 2.4 and 2.5, the voltage across the nonlinear inductor can be written as:  V ^V -(IR) I=  2  S  2  +  - ^ -  (2.6)  CO C  The first term in the right-hand side of Equation 2.6 {sjv  2  - (IR)  2  ) represents an ellipse  whose main axes have the magnitude of V and — r , and the second term is a straight line R s  having slope of — — .  A d d i n g these two quantities represents an oblique ellipse, whose  CO  intersection with the characteristic of V, presents the three possible states of the oscillation of the circuit. Figure 2.8 shows the graphical solution for the ferroresonance circuit of Figure 2.7.  11  /(/)  Current  F i g u r e 2.8: G r a p h i c a l s o l u t i o n of the ferroresonance c i r c u i t  Points 1 and 2 in the Figure 2.8 represent  the stable solutions, whereas point 3 represents  an unstable solution. T o show this, rewrite Equation 2.6 as:  (IR)  2  If the quantity  = V;  (2.7)  co C  increases in magnitude with an increase in current, then  co C  according to equation 2.7, (IR)  2  tends to decrease, and this suppresses any further increase  in current. Thus stability is achieved. However, if the quantity  12  coC  decreases in  magnitude, with an increase of current, the magnitude of (IR) tends to increase and under 2  this condition, the current continues to increase and the solution is unstable. The dashed area i n the Figure 2.8 shows the variation of the magnitude of  Figure 2.8 corresponds to an unstable solution since  /(/)  coC  co C  . Point 3 i n  decreases with an  increase of current. In the same figure, points 1 and 2 represent the stable solutions.  2.3 Brief History of Ferroresonance T h e phenomenon of ferroresonance was first observed i n a power system i n the 1930's. Interest was generated when it was recognized that the use of series capacitors for voltage regulation  caused  ferroresonance  i n distribution systems,  resulting  i n damaging  overvoltages [6]. Since then, companies such as Pacific Power and Light, and Idaho Power, with distribution voltages of 19.9/34.6 k V , have reported a number of failures of underground wye-delta  transformers when switching at the transformer pole, due to  ferroresonant voltages [1]. In these reports it was shown that even with no overhead line connecting the transformer and the single-phase switch, ferroresonance overvoltages can occur as a result of internal capacitances  of the transformers. Another very severe  occurrence of ferroresonance i n a 1000 M V A , 525/241.5 k V Y-conneeted autotransformer in Big E d d y (Oregon-USA) was reported by the authors in [7]. A very high increase of temperature  i n the bank was reported and a chemical analysis showed significant  decomposition of both oil and solid insulation.  13  O n e of the first analytical works i n the area of ferroresonance was done by Rudenberg i n the 1940's [5,8]. Since then, the research i n this area has been divided into two main categories; first, improving the models used to predict the behavior of the transformers and second, studying ferroresonance involving these transformers as installed i n an electrical power systems. T h e work by Swift [9] and Jiles [10] has shown insight into the core behavior of transformers and the separation of hysteresis and eddy current losses. Feldman [11],  and Frame  [12],  nonlinearities i n saturable  have  developed piecewise-linear  inductances.  methods  of modeling the  Hopkinson [13], performed system tests and  simulations o n the effect of different switching strategies on the initiation of ferroresonance in  three  phase  systems.  Smith  [2],  showed  and categorized  different  modes  of  ferroresonance i n one type of three- phase distribution transformer based o n the magnitude and appearance of the voltage waveforms. W o r k has been done by others to categorize the various modes of ferroresonance,  but these works generally compare  the qualitative  characteristics of the voltage waveforms i n the time domain. F o r the first time in 1991, Kieny [14] showed that simple-time domain simulation is not enough to understand the behavior of the ferroresonant circuit. Given the wide variety of possible solutions, he recognized that besides normal and subharmonic states, a pseudo-periodic character of this phenomena can be observed. H e suggested a method of drawing bifurcation diagrams and recognized two types of bifurcation: limit points, and H o p f bifurcation.  14  2.4 A Brief History of Nonlinear Dynamical Systems and Chaos Today, dynamics is an interdisciplinary subject, a subject that deals with change, with systems that evolve in time. Interest in the subject began in the mid-1600's, when N e w t o n invented differential equations, discovered his laws of motion and universal gravitation, and combined them to explain Kepler's law of planetary motion. F o r some, the study of dynamics began and ended with Newton's law, F = ma . It was believed that if the forces between particles and their initial positions and velocities were given, one could predict the motion or history of a system forever into the future.  T h e breakthrough came with the  work of H e n r i Poincare, a French mathematician, who united nonlinear dynamics with the geometric ideas of topology. H e introduced a new point of view that emphasized qualitative rather than quantitative questions. F o r example, instead of asking for the exact position of the planets at all times, he asked, "Is  the solar system stable forever, or will  some planets eventually fly off to infinity?" [15]. H e was also the first person to glimpse the possibility of chaos. In understanding the sensitivity of chaotic behavior to the initial conditions of the system, thereby  rendering long-term prediction impossible, he wrote  it may happen that small differences in initial conditions produce very great ones in the final phenomenon. A small error in the former will produce an enormous error in the latter. Prediction becomes impossible....". In the first half of this century, dynamics was largely concerned  with nonlinear oscillators  and their  applications i n physics and  engineering . Nonlinear oscillators played an important role in the development of such technologies as radio, radar, phase locked loops, and the laser. It was with the invention of high-speed computers in the 1950s that a breakthrough in the study of nonlinear dynamical 15  systems was achieved. T h e computer made it possible to experiment with equations in a way that was impossible before, and thereby to develop some intuition about nonlinear systems. In the early 1960's a prominent meteorologist, Edward Lorenz, who was engaged in computer  modeling, wondered whether  existing science  was capable  of weather  prediction. H e studied a simplified model of convection rolls in the atmosphere to gain insight into the  notorious  unpredictability of the  weather. Lorenz made  a  major  breakthrough in the understanding of modes of operation of a nonlinear system by noticing that small variations in initial conditions (one digit in the fourth decimal) caused the result of exact weather simulations to diverge rapidly. H e found that the solutions to his equations never settled down to equilibrium or to a periodic state. Instead they continued to oscillate in an irregular, aperiodic manner. Lorenz also showed that there was structure in this irregular motion and his solution, when plotted in three dimensions, fell onto a butterfly-shaped structure. Therefore small differences in input could quickly become large differences in the output and this phenomenon became k n o w n as "sensitive dependence on initial conditions". This effect became half-jokingly k n o w n as the Butterfly Effect, referring to how a butterfly stirring the air today in Peking can transform storm systems next month in N e w Y o r k [15].  It took more than a decade and a half for the theory of chaos to gain  enough recognition to be named, and it took even longer for the investigation of chaos to earn scientific respectability. Other major theoretical contributions came in the 1970's. Feigenbaum [16] discovered that there exist certain universal laws governing one of the routes of transition from regular to chaotic behavior, and his work established a link between chaos and phase transition. In the m i d 1970's, 16  Benoit Mandelbrot developed a  new non-Euclidean "fractal" researchers  geometry  and coined the w o r d "fractal"  [17].  Today  such as Gollub, Libchaber, Swinney, Linsay, M o o n , and Westervelt [18] are  testing the new ideas about chaos in experiments on fluids, chemical reactions, electronic circuits, mechanical oscillators, and semiconductors.  2.5 Definition and Categorization of Methods for Chaotic Systems The theory of chaos, in addition to periodic attraction, recognizes non-periodic or chaotic solutions in the behavior of nonlinear dynamical systems. These solutions, although nonperiodic, are not only deterministic, but also bounded, and represent the long-term behavior of the systems. T h e recognition that a deterministic dynamic system can exhibit pseudo-random behavior is one of the most important conclusions of the new theory. Therefore the term chaos can be defined as aperiodic long-term behavior in a deterministic system that exhibits sensitive dependence on initial conditions. In this definition there are three important ingredients:  1.  Aperiodic long-term behavior means there are trajectories which do not settle down to  fixed points, periodic orbits, or quasi-periodic orbits as time goes to infinity. 2.  Deterministic means that the system has no random or noisy inputs or parameters. T h e  irregular behavior arises from the system's nonlinearity rather than from noisy driving forces. 3.  Sensitive dependence on initial conditions means that two trajectories starting very  close together will rapidly diverge from each other, and therefore will have a different  17  future. T h e practical implication is that long-term prediction becomes impossible in a system like this, where small uncertainties are amplified very quickly. There are some standard tools to study dynamic systems [18].  T h e system must be  simulated to get the time series, which is simply the values of the variables as functions of the time. T h e time series of a periodic motion has the appearance of a uniform trace but it is not easy to detect the presence of more than one basic period by examining the time series. So, for motions that contain more than one basic period, examining the time series may not be sufficient to characterize the motion. Also, it is difficult to distinguish between three or higher-period quasi-periodic motions and chaotic motions by examining the time history.  A  better understanding can be obtained by getting the phase plane trajectory of a chaotic  solution [18]. Phase space methods view the trajectory of a waveform vs. its derivatives as time progresses. These diagrams serve as a unique form of identification for any kind of waveform. A phase plane trajectory will have a closed path if the waveform is periodic. These types of trajectory paths are called periodic attractors, because final steady results will be attracted to one of  these trajectories regardless of initial conditions. T h e phase plane  trajectory of a chaotic waveform is attracted to a constrained region but seems to wander randomly within that region and never closes back on itself. T h e pattern formed by this trajectory  as it passes randomly through this constrained region is called a "strange  attractor". Strange attractors can also be defined as attractors that exhibit sensitive dependence on initial conditions. Strange attractors were originally called strange because  18  they were associated with a new geometric object called a fractal set [17]. In a threedimensional phase space, the fractal set of a strange attractor looks like a collection of an infinite set of sheets or parallel surfaces, some of which are separated by distances which approach  the  infinitesimal.  The  above  definitions were  for  damped systems  but  conservative systems can also exhibit the same types of motion as lossy systems, including periodic, subharmonic, quasi-periodic, and chaotic motions.  However, one of the main  differences between vibrations in lossy and lossless systems is that chaotic orbits in lossy systems exhibit a fractal structure in the phase space whereas chaotic orbits in lossless systems do not.  A s an example, consider a nonlinear spring where the resorting force is not linearly proportional to the displacement. F o r the case of equal effects of compression and tension, the equation of motion takes the following form:  dx  dx  dt  dt  2  ,  —y + 2y —-+a x + P x = cosco t  Equation  (2.8)  2.8 is called the Duffing oscillator after the mathematician who studied it [19].  Figure 2.9 shows the final state-space trajectories of Duffing's equation for a number of different values of the parametery  and for co = 1 .  19  (a)  2  -3  .  DB0  .  898  Period-1 m o t i o n  (y  = 0.05,  (b) Period-3 m o t i o n ( y  20  (X =  = 0.1,  a  1,  P=  1)  = 1, P =  1)  2 . BOO  e.aaa  (c) P e r i o d - 5 m o t i o n (y  (c) C h a o t i c m o t i o n  = 0.15,  0C = 1, P =  1)  (y = 0.175, a = 1, P = 1)  F i g u r e 2.9: (a) Phase-plane trajectory o f D u f f i n g o s c i l l a t o r f o r p e r i o d - l m o t i o n . (b) P e r i o d 3 m o t i o n , (c) P e r i o d - 5 m o t i o n , (d) C h a o t i c m o t i o n .  Another tool for understanding chaotic waveforms is the plotting of just one point per cycle of the periodic forcing function on the phase plane.  The family of points "strobed"  in this manner is referred to as a Poincare section [20].  To  21  show  this  consider  an  n-dimensional autonomous system X = f(X)  with  X(t0)  - X0,  and let S be an n-1  dimensional surface of a section. It is required that S be transverse to the flow and that all trajectories starting o n S flow through it and not  parallel to it. T h e Poincare map P is a  mapping from S to itself, obtained by following trajectories from one intersection with S to the next. If X  n  X  = P(X ).  n+X  n  e S denotes the nth intersection, then the Poincare map is defined as  Figure 2.10 shows the periodic orbit of an autonomous system in the  corresponding state space. In this figure if X  p  at X  p  returns to X  system of X = f (X).  p  is a fixed point of P, then a trajectory starting  after some time T , and is therefore a closed orbit for the original  In nonautonomous systems, the period associated with a periodic  orbit is explicitly k n o w n and this period can be used to construct a Poincare section.  s  Figure 2.10: Poincare Mapping Therefore, a period-one solution of a continuous time system will correspond to a fixed point in the Poincare map. A n nth-order subharmonic will correspond closed orbit of the Poincare map and will therefore have V  22  to a period-n  fixed points. If the motion  CO  consists of two incommensurate frequencies, where  is an irrational number, the  co  Poincare map will become a continuous closed orbit. Finally, if the Poincare map does not consist of either a finite set of points or a closed orbit, the motion may be chaotic. If the system is undamped, the Poincare map of the chaotic motion will appear as a cloud of unorganized points in the phase plane. For damped systems, the Poincare map will appear as an infinite set of highly organized points. Figure 2.11 shows the Poincare map of period1 (a), period-2 (b), period-3 (c) and the chaotic solution (d) of the Duffing equations.  2 . eea  -i—i  1  1  -* xC 2  1—i  1  1—r  1  T—i  i  r  I  L  H  J  -2.aee -z.eea  J  -x.aaa  I  I  i_  J  _ J  a . eea xCll  (a) P o i n c a r e m a p o f period-1  23  J  I  L  2 . eae  - I — I — I — I — I — I — I —  1  - I — I — I — I — I — I — I  r  I  r  1.090  0.000  H  1  1  H  H  H  H  1  H  1  1  1  '  "  1"  XC21 -1.O00  _l -2.000 -2.000  I  1  J L -l.000  1  1  -|  L  I 1 1.00O  I  I  1_  '  xCll  (b) Poincare map of period-3  8.000  ~i—i—i—i—i—i—r  T — i — i — i — I — i — i — i — r  1.000  0 . 000  H  1  1  1  K-  H  1  J 0 000  L  1  1  1  h-  -1.000  I  I  I  K  i  l  l  (c) Poincare map of period-5  24  J  I l _  2 .  eea  -1 - 1  1  1  1  1  1  1 1  1—r-  "i • "  i  |  i  i  -s^;...  •  1  1  i  i  1  1  l.aaa  "  .  .  x\x  1  1  1  i  i  i  1  1  -  ;  1  1  1  1  |  1  1  I  1  1  -x.eea  1  1  1  1  a. aea  - 1 .  X C l l  (d) P o i n c a r e m a p o f the c h a o t i c s o l u t i o n  F i g u r e 2.11: (a) P o i n c a r e m a p o f D u f f i n g o s c i l l a t o r f o r p e r i o d - l m o t i o n . (b) P e r i o d - 3 m o t i o n . (c) P e r i o d - 5 m o t i o n , (d) C h a o t i c m o t i o n .  Another visible characteristic of a chaotic signal is its frequency power spectrum. T h e power spectrum of a chaotic system will contain not only discrete frequency components due to the forcing function and its harmonics, but also distributed frequency components. Therefore, the appearance of a broad spectrum of frequencies in the output when the input is a single frequency harmonic, is a clue to detecting chaotic vibration.  So far the tests for chaotic motion described above are mainly qualitative and involve some judgment and experience on the part of the investigator. There is, however, a widely used quantitative test for chaos, the Lyapunov exponent. The Lyapunov exponents associated with a trajectory are essentially a measure of the average rates of expansion and contraction of trajectories surrounding it. T h e y are asymptotic quantities, defined locally in state space, and describe the exponential rate at which a perturbation to a trajectory of a system grows or decays with time at a certain location in the state space. This means that if two  25  trajectories start close to one another in phase space, they will move exponentially away from each other for small times on the average, if the Lyapunov exponent Therefore, if d  Q  is >  0.  is a measure of the initial distance between the two starting points, the  distance a small time later will be:  d(t) = d 2 ' 0  The  X  symbol A, is defined as the Lyapunov exponent. Because in a bounded system  d(t) cannot go to infinity, the divergence of chaotic orbits can only be locally exponential. Therefore, to define a measure of this divergence of orbits, we must average the exponential growth at many points along a trajectory. T o start, one begins with a reference trajectory  and a point on a nearby trajectory and measures — .  W h e n d(t) becomes too large and  "o  the growth departs from exponential behavior, one looks for a new nearby trajectory and defines a new d0(t).  T h e n the Lyapunov exponent can be defined by the following  expression [18].  r ^ ^ N  l  l  'O  o  g  2  7 7 T ^ n "OVA  k=\  l  )  The dependency of a chaotic system on initial conditions means that if two trajectories start close to each other in phase space, after a small time on the average, they will move exponentially away from each other. Since a positive Lyapunov exponent indicates expansion, what distinguishes chaotic motion from non-chaotic motion is the existence of  26  at least one positive Lyapunov exponent. T h e magnitude of the exponent reflects the time scale for which the system dynamics become unpredictable.  The calculation of the Lyapunov dimension can also be used to classify chaotic attractors from non-chaotic attractors. F o r non-chaotic attractors the dimension is an integer value, and all chaotic attractors have non-integer dimensions [19].  F o r example, if the attractor is  a line or a simple closed curve, we say that the dimension is equal to 1 because a line or a curve is a one-dimensional object. Similarly, a surface has a dimension of 2 and a  solid  volume a dimension of 3. Geometric objects with a dimension that is not an integer play an important role in the dynamics of chaotic systems. These geometric objects have been named 'fractals" [Mandelbrot, 1982] because their dimension is not an integer.  Kaplan and  Yorke [21] have suggested that one can calculate a dimension for a fractal attractor based on the Lyapunov exponents. The Lyapunov dimension d, is defined to be:  j  IX d  ' = J  +  J  =  r -  The modeling of nonlinear systems will always lead to a system of nonlinear differential equations depending upon various physical parameters. T h e type of solution (periodic, quasi-periodic, nonperiodic, etc.) for a such a system will depend on the values of these parameters. In order to predict which modes of ferroresonance may occur for a wide range of various parameters, a bifurcation diagram will be used [14]. Bifurcation, a French word introduced into nonlinear dynamics by Poincare, is used to indicate a qualitative change in  27  the features of a system, such as the number and type of solutions, under the variation of one or more parameters on which the considered system depends. T h e bifurcation diagram is a plot of the magnitude taken from a family of Poincare sections vs. the system parameters that are being varied.  2.6 Nonlinear Dynamics and Chaos Applied to Ferroresonance T h e problem of ferroresonance can be categorized as a nonlinear resonance which causes damage  in power  distribution and transmission  systems.  One  of the  fundamental  properties of ferroresonance is the fact that several stable solutions can exist under steadystate conditions for a given circuit. F o r example, remnant magnetization of the cores, the voltage at the time of switching and the amount of charge on the capacitance are all initial conditions which determine the steady-state response. W i t h a small difference in these initial conditions, the ferroresonance overvoltages can have very different waveforms. The purpose of this work is to investigate  ferroresonance from the standpoint of nonlinear  dynamics and chaotic systems. This thesis uses the theories of nonlinear dynamics and chaos to analyze, understand, and control ferroresonance. W o r k in this area has been done by others to categorize the various modes of ferroresonance in commonly used sizes of fivelegged core transformers [2], but this work simply compares the qualitative characteristics of voltage waveforms in the time domain.  The realization that ferroresonance  is a  nonlinear and sometimes chaotic process opens up many possibilities. It can be shown that the transition between periodic and nonperiodic modes commonly occurs because of small changes in circuit parameters or in initial conditions. Therefore, the newly-developed  28  techniques for analysis of nonlinear dynamical systems and chaos should now be evaluated for use with ferroresonance.  29  3  Chapter  SYSTEM DESCRIPTION A N D  3.1  MODELING  Introduction  In this chapter the equivalent circuit of a typical ferroresonant circuit i n a power system is developed. T h e set of differential equations that describes the model of the system is derived, and the parameters of the model and their per unit values are calculated.  3.2 The Equivalent Circuit Most distribution systems make wide use of grounded-wye to grounded-wye transformers to serve three phase loads. In modeling these systems, the distribution line is represented by its R L C ^-equivalent, and the three-phase circuit breakers and gang-operated switches are used at the substation where distribution lines originate. T h e distribution lines connecting transformers to the system source may be either overhead lines or underground cables. However, cables have a relatively large shunt capacitance compared to overhead lines, and for this reason ferroresonance most often involves underground cables. Threephase or single-phase transformers can appear at the end of a distribution line or at any point along the line. Figure 3.1 shows a simplified schematic of such a system.  30  X  - —^N\r-^ffi^  ^  1  T  X  X  X  X  F i g u r e 3.1: A t y p i c a l d i s t r i b u t i o n system s u p p l y i n g a three phase l o a d through a grounded-wye to grounded-wye transformer  In most instances, ferroresonance occurs when one or two of the source phases are lost while the transformer is unloaded or lightly loaded. This can happen due to the clearing of single-phase fusing, operation of single-phase reclosers or performance of single-phase switching procedures. In this case if one of the three switches was open, and only two phases of the transformer were energized there would be a voltage induced i n the "open" phase. This induced voltage "backfeeds" the distribution line, back to the open switch. Ferroresonance will occur if the distribution line is highly capacitive. This ferroresonance involves the nonlinear magnetizing reactance of the transformer's open phase and the shunt capacitance of the distribution line.  In recent years, there has been some work done on chaos in power systems. Five recent papers [22,23,24,25,26] addressed the solution of the nonlinear equation for a typical ferroresonant circuit in power systems. T h e system studied here consists of a source feeding  31  an unloaded power transformer with one of the supply conductors being interrupted (Figure 3.2).  F i g u r e 3.2: M o d e l f o r ferroresonance c i r c u i t i n c l u d i n g l i n e capacitance  The  transformer is then energized through the capacitive coupling with the other phases.  Figure 3.3  shows the circuits that feed the disconnected coil through the  capacitive  coupling.  F i g u r e 3.3: T h e c i r c u i t that feeds the d i s c o n n e c t e d c o i l  In order to derive the mathematical description of the above circuit, Thevenin's theorem is used to obtain an equivalent circuit. The equivalent capacitance can be found by shorting the two remaining voltage sources of phases 1 and 2. In doing so, both corresponding ground capacitance and transformer windings are shorted and can be omitted. W i t h phases 32  1 and 2 having the same potential, there is no current flow through the mutual capacitance connecting the nodes 1 and 2, and it can be omitted as well. Therefore, the remaining equivalent circuit consists of two mutual and one ground capacitances which all connect node 3 to ground potential. Thus, the equivalent capacitances can be grouped as:  C = C +2C g  (3.1)  m  Next, the equivalent source voltage is derived. Assuming steady-state conditions before interrupting phase 3, the current I  3  is supplied from the voltage source i n the third phase.  T h e amount of current which arrives at the transformer winding to contribute to the flux can be found by subtracting the fraction of the current that is  lost due to the ground  capacitance and the mutual capacitance, such that:  i =j<oC (V -V ) 2  m  3  ( i  2  h = h -i, - i " , ~h = h -ja> C V - / c o C (2V -V g  3  m  3  2  - >  -V ) x  In equations 3.2 i, stands for the current i n phase 3, arriving at the transformer end of the line. Since the voltages V  lt  V  2  and F  3  are of identical magnitude, each delayed by a 1 2 0 °  phase shift, the phasor addition of V , V x  and V is zero. Therefore, equation 3.2 can be  2  3  written as:  i =I -ja>(C +3C )V L  3  K  m  33  3  (3.3)  The voltages V , V x  2  and V differ only in their phase shift and they can be rotated by 120° 3  without loss of generality. After the switch opening in phase 3, the source current goes to zero and equation 3.3 can be written as:  i =-ja>(C +3C )V t  L  m  (3.4)  l  This current has to be equal to the current in the equivalent circuit of Figure 3.4.  c  F i g u r e 3.4: P r e l i m i n a r y equivalent c i r c u i t  i =ja>C(E-V ) L  Where C = C + g  (3.5)  L  2C„,, and V, = V . x  Relating equations 3.4 and 3.5, the equivalent source E can be calculated as:  £=  " y +2C  C  C  (3.6)  Under loaded operating conditions, the flux induced in the primary transformer winding is immediately compensated by the current flow in the secondary winding so that the core losses in the transformer can be neglected. However, in the case of an unloaded (or very lightly loaded) power transformer, current can develop  34  in the secondary side and this  causes the flux to emanate from the transformer leg and flow through the iron core. This in turn increases the transformer core losses, which therefore can no longer be neglected. Figure 3.5 shows the final reduced equivalent circuit with the  inclusion of R  which  represents the transformer core losses. C  F i g u r e 3.5: R e d u c e d equivalent c i r c u i t  3.3 Differential Equation of the System In  the peak current range for steady-state  operation, the  flux-current  linkage can be  approximated by a linear characteristic such as:  / = a where the coefficient a' 1  4>  (3.7)  of the linear term corresponds closely to the reciprocal of the  inductance (a « —). However, for very high currents the iron core might be driven into  saturation where the flux-current characteristic becomes highly nonlinear. In this study thecj) - / characteristic of the transformer is modeled by an eleventh-order polynomial [22].  35  i = aisf +bty  (3.8)  11  where a = 2.8xl(T ; 3  2> = 7 . 2 x l ( T ; 3  / is the current i n pu value; and <|) is the flux in the transformer core in pu value. The polynomial of order eleven and the coefficient ' b' of equation 3.8 along with only two terms of the polynomial are chosen for the best fit of the saturation region. Figure 3.6 shows the comparison between different approximations of the saturation region against the true magnetization characteristic that was obtained from field measurements by D i c k and Watson [27].  Using a polynomial of order less than eleven to  represent  the  magnetization curve might be adequate for small capacity transformers, but it does not bend sharply enough at the knee point to satisfy the magnetization characteristic of modern high-capacity transformers. Also as is shown in [24],  during ferroresonance  with the  operating point of the transformer located in the saturated region, and for the core materials used in modern high-voltage power transformers, significant,  and a single-valued curve  can be used for the  magnetization curve.  36  hysteresis  loops are  representation  of  not the  i=* xirr 7  2  /  i = £>x10-2  ia  (0.28*.0.724")x10-  2  • transformer  Figure 3.6: Approximation of the saturation region of the test transformer  In chapter 5, the effect of changing the exponent of the nonlinear term o n the modes of behavior of the ferroresonance circuit will be investigated.  The voltage across the transformer winding of Figure 3.5 can be defined as:  (3.9)  This is the same voltage across the resistance R, and the current through it can be written as:  h  ~  (3.10)  R~ R  The current through the equivalent capacitance C is  i =C  , dt  L  = C E - C 6 Y  37  (3.11)  Knowing  all the branch  currents  from  equations  3.9  to  3.11,  along  with the  $ - / characteristic of the transformer, one can write the differential equation for the circuit of Figure 3.5 as:  d6  1 di*  2  1  ;(a*+&* RC dt+ 7C  dt  1  H  , ,  ) = ©,£cos(a>,0  (3.12)  where C0j is the frequency, l.Opu; and E is the peak value of the voltage applied to the equivalent circuit of Figure 3.5 in pu.  3.4 Determination of the Parameters and Per Unit System The transmission line that feeds the transformer in Figure 3.5, has a line-to-ground capacitance,  C =5A\ F/j g  and a  n  (m  line-to-line capacitance,  C„, = 1.18  n  •  The  transformer that is used in this study is a 3-phase 2 5 M V A , 110/44/4 k V autotransformer with a very sharp bending <j) - / characteristic obtained from field measurements by D i c k and Watson [25]. It is assumed that the transformer is connected i n wye o n the line side, and the delta winding o n the load side is open. The transformer is assumed to be unloaded or very lightly loaded. The voltages applied to the line are assumed to be equal to the rated voltage of the transformer V = V = V = V x  2  3  raled  . The current / = 1.0 pu corresponds to the  transformer's rated current of 131A and flux (j) = 1.0 pu corresponds to the transformer's phase-to-neutral rated voltage of V = 63.5 kV.  T h e base impedance is calculated as  Kase  Zbase ~ base  =  484.0Q.  T h e transformer  core losses  38  are taken as 1% of the rated  transformer R  base  capacity.  A  100%  transformer  loss  corresponds  to  a  resistance  of  = 484.0Q. However , the loss is inversely proportional to R , so that 1% loss equals a  value of 100 x R  iase  = 48.4 kQ.. W i t h  the above assumptions, the parameters of the  equivalent circuit of Figure 3.5 are calculated as:  C = C  2C =7.77  t +  E = (-£)V  a  l  =0.159F,  n F  /  h n  (3.13)  R (1 % loss) = 48.4KJ  The equivalent circuit of a typical ferroresonant circuit was presented i n this chapter. T h e set of differential equations describing the circuit was developed and the parameters of the model and their per unit values were calculated. T h e model presented in this chapter will be used in Chapter 4 to obtain the analytical solution for the ferroresonance circuit.  39  4  Chapter  ANALYTICAL SOLUTION  4.1  Introduction  In this chapter, the method of Slowly Varying Amplitude or S V A [3] is used to derive an analytical solution for the equivalent ferroresonant circuit of Figure 3.5. Using this method, it will be shown that the system exhibits a saddle-node bifurcation leading to the occurrence of two stable attracting points and one region of an unstable saddle point which leads to a chaotic response.  4.2  Analysis of N o n l i n e a r E q u a t i o n  T w o approaches to the existence and stability of periodic solutions are credited to V a n der P o l and K r y l o v and Bogoliubov [28]. These two methods are similar in many ways, and are mostly confined to a first approximation. The method of slowly varying amplitude (SVA), devised by V a n der Pol, is an approximate method that can be applied to a class of nonlinear oscillation systems. Details of this method have been discussed by Nayfeh and M o o k [3]. T h e S V A method will be used to derive an analytical solution for the equivalent ferroresonance circuit of Figure 4.1.  40  -K  ^  v,  F i g u r e 4.1: E q u i v a l e n t ferroresonance c i r c u i t  W i t h respect to Figure 4.1, the following equations can be written:  (4.1)  h ~ h ~ h = °  V =-£  = $ =  r  R (i -i ) 2  ]  (4.2)  2  V =e(t)-RJ,-R (i -i ) c  1  2  C  x  (4.3)  2  dt de(t) dt  1  di,  dV,  dt  dt  J_d§_ R  2  (4.5)  dt  i = a§ + b§ "  (4.6)  2  Ce- R C^jx  -  C<j>'-4- -  <|)  (4 +  a<j) +  J/,  iTc '~dT c + R  +  41  (4.4)  H ") = 0  "  = e  (4.7)  (4.8)  di, di. If R, « R? and — »> — then: dt dt 1  2  di  di d§  dt  dt dfo  x  2  (4.9)  di  di d§  x  2  dt ~ disf dt 1  * lTcv x +  *  +[Rl(cf+nb  2  .  •  m  +  a<b + bi>"  c  =  e  '  (4 10)  The above equation is a variant of the Duffing equation [18], with the magnetizing curve having an index of n - 11, instead of n = 3. F o r e = E sinco t, equation 4.10 can be written as:  * + /,(<t> )4> + fi(4>) = co Ecos(co  /)  (4.11)  wnere  /.(•) = ^ /(<l>) =  4  2  + [*,(fl + »&T" )] 1  (4.12)  + &<)>"  C  Equation 4.11 can be time scaled by setting x = co t.  d$ 2  M)  —~2~ +  dx  co  4 / «o E 2  or  —+  co  j— = — C O S T  co  Equation 4.13 can be amplitude scaled by letting cp = — .  42  (4-13)  dx  su*  di  co  sco  T h e n o n l i n e a r restoring c o m p o n e n t o f equation 4.14 is:  fiJW  sco  )  _£_  =  _A_  Ceo  2  2 9  »  (p+TT^^'q)" Ceo  (4-15)  2  T h e n o n l i n e a r d a m p i n g c o m p o n e n t o f t h e e q u a t i o n 4.14 is  •"  V  V  co  ' — ( a * , + — ) + -* RC  co  2  —Z—  (4.16)  co  W i t h the f o l l o w i n g definitions:  - ( « * « , ' = n b R co RC P  l  s ' - \  8=—,  2  as  s=  '-'p-^ V b  e q u a t i o n 4 . 1 6 c a n b e w r i t t e n as:  ^  ^ co  a  + pep"-'  (4.17)  F o r a s m a l l d i f f e r e n c e b e t w e e n t h e n a t u r a l f r e q u e n c y a n d t h e e x c i t a t i o n f r e q u e n c y , w h i c h is defined b y  a d e t u n i n g f a c t o r p,  where  p  is t h e d e v i a t i o n o f t h e o s c i l l a t i o n f r o m  frequency of excitation, w e can write:  43  the  (4.18)  whe re: a  2  CO  =  Ceo  :  W i t h the above definition of p, equation 4.14 can be written as:  ^ - y + (<x +p<p"" )^7 + (l-/?)q> + <p" =8cosx  (4.19)  d <f> d(D „ „_, d<p „ —Y + cp = - a ^ - - P c p + /?cp - cp + 8 COST dx dx dx'  (4.20)  1  dx  dx  or 2  0  T h e flux linkage is now considered to vary slowly both in amplitude and in phase, and can be considered to be of the form:  (4.21)  cp = A(x )cos(0 (x ))  Using the method of K r y l o v and Bogoliubov from [28] and Migulin [29,30], and taking the first and second derivatives of the equation 4.21 and substituting it into equation 4.20, it can be shown that:  • 8 .„ aA .4 =-sin6 2  2  „ P H  1  „-z  2"  2  1  The average rate of change in frequency is calculated as:  44  A".  (4.22)  .  -pA  AQ =—^2  A" 5 + — --cos0 2" 2  (4.23)  T o determine the character of the solution, the singular points must be found and the motion i n their neighborhood must be examined. T h e nature of trajectories  i n the  neighborhoods of singular points shows whether a small perturbation i n the steady-state motion decays or grows; that is, they illustrate the stability of the steady-state motion. T h e  steady-state motions occur when A = 0 = 0 , which corresponds to the singular points of equations 4.22 and 4.23. Therefore, the steady-state amplitude and the frequency of oscillation of the system are obtained by equating the right-hand sides of equations 4.22 and 4.23 to zero.  pA A" 5 — — = — — - — COS0 2 2" 2 rt  6 a - s i n 0 =~A 2 2  + Q K  1 ,71-2  (4.24)  1 2"  A"  (4.25)  Eliminating 0 from both equations 4.24 and 4.25, we calculate the value of p as:  ,n-l  n-l P  ~  2""  1  \A  ±  2  (4.26)  where p is the measure of detuning of the system from resonance. T h e plot of A as a function of p for a given E and system parameters R and C is called a frequency-response  45  curve [3]. Each point on this curve corresponds to a singular point in a different state plane; there is one state plane for each combination of parameters. For the case of no series resistance (./?, = 0.0) and n = 11, equation 4.26 can be written as:  A  xo  I E  :  1  From equation 4.27 it can be seen that p is a function of excitation E and the system parameters i? andC. Figure 4.2 shows the graph of A versus p, for the case of 2  R = 1%(4.84&Q), C = 111nF for \00km and E = 0.15pu. 2  Detuning Factor  Figure 4.2: Graph of A as a function of  46  p  Figure 4.2 shows that the effect of the nonlinearity is to bend the amplitude curve which causes multivalued regions. This multivaluedness of the  response curve will lead to jump  phenomena. T o explain this phenomenon, in Figure 4.3, the amplitude of the excitation is fixed, and the frequency of the excitation is varied slowly up and down through the linear natural frequency, and the amplitude of the harmonic response is observed.  F i g u r e 4.3: J u m p p h e n o m e n o n  Starting from point 1 in Figure 4.3, as the frequency is reduced,  p decreases and A slowly  increases through point 2 until point 3 is reached. A s p is decreased further, a jump from point 3 to point 4 takes place with an increase in A, after which A decreases slowly with decreasing p.  However, if we start at point 5 and p  is increased,  A increases slowly  through point 4 until point 6 is reached. A s p is increased further, a jump from point 6 to point 2 takes place with an accompanying decrease in A, after which A decreases slowly with increasing p.  T h e maximum amplitude corresponding to point 6 is attainable only  when approached from a lower frequency. The portion of the response curve between point 3 and 6 is unstable. W h e n two stable steady-state solutions exist, the initial conditions 47.  determine which of these represents the actual response of the system. Figure 4.4 shows the influence of the damping coefficient R  on the response curve. In the absence of the  damping the peak amplitude is infinite, and the frequency-response curve consist of two branches. However, in the presence of damping, the peak amplitude is finite.  A  P  F i g u r e 4.4: Effect of d a m p i n g o n the frequency-response curve.  F r o m Figure 4.2 it can be seen that beyond a critical value of the amplitude of the flux linkage, hysteresis appears, leading to occurrence of two stable attracting points and one region of unstable saddle point, which could possibly lead to chaotic response. W h i c h of the two stable solutions is picked will depend on the initial conditions. T h e stable attracting point exhibits both fundamental as well as harmonic and subharmonic regions  of  operation. T h e stability of the different portions of the response curves can be determined by investigating the nature of the singular points of equation 4.24  and 4.25. This is done  by eigenvalue analysis of the Jacobian matrix of equations 4.24 and 4.25 [3].  48  F o r i?, — 0.0 and n = 11, matrix of 4.23 can be written as:  pA J(A,B) =  2 p  2  + 5.3711 x 1 0 ~ ^ 3  + 4.882 x l O  - 4  ^  1 1  (4.29)  a  10  2  The eigenvalues of the above matrix are calculated by setting X I - J(A,Q )| = 0 .  X + A. (-0.5a +0.5a A)- — ( - a 2  2  -p ) 2  + 2.929x IO' pA 2  u  -2.622x 1 0  - 6  ^  2 1  =0  Solving the above quadratic equations will result i n two solutions. T h e real part of the eigenvalues is then plotted against the amplitude of the flux linkage at a fixed value of the amplitude of the driven function. Figure 4.5 shows the graph of the two eigenvalues vs. A . The existence of a saddle for flux linkage magnitude beyond 2pu is seen from the graph. U p to the value of 2 pu, the two eigenvalues have the same negative sign and hence the steady-state motions are stable. However, the eigenvalues of the characteristic equation after the flux linkage of 2 pu are real and opposite i n sign which corresponds to a saddle point. T h e graph clearly demonstrates the existence of a saddle point or, i n other words, the existence of an intersection of a stable and an unstable manifold in the phase plane.  49  It can be said that one of the keys to chaos is the search for the existence of horseshoe-like maps [18]. In many systems, as with ours, strange attractors are organized around saddletype fixed points. These are the points at which the orbit comes into a small region of phase space along a stable manifold and goes out along a path close to an unstable manifold. These saddle-type fixed points create the contraction and stretching mechanisms for a horseshoe map. Repeated iteration of the stretching, contracting, and folding of regions of the phase space leads to unpredictability and, in the case of dissipative systems, to the fractal nature of the dynamics in the map.  0.5  1  Amplitude  1.5  F i g u r e 4.5: T h e graph of the a m p l i t u d e of the f l u x l i n k a g e vs. the eigenvalues  50  2.5  Therefore, the calculation of fixed points and determination of the nature of their stability are not just academic exercises; the results can be used along with computational tools to give a clue to the nature of the chaotic attractor.  In the previous procedure, the concept of the jump resonance phenomenon, as the driving frequency was varied, was considered. In the graphical technique discussed in Chapter 2, it was shown that jumps in the amplitude could be expected as the driving source amplitude was varied and the driving frequency was held constant. T o gain further insight into this aspect of the ferroresonance , the Equivalent Ritz Method [31] will be used to obtain the graph of the variation of the source voltage vs. the amplitude of the flux linkage. The damped ferroresonant equation from the Figure 3.5 can be written as:  In equation 4.30,  we have lumped the phase angle due to the damping in the system into  the driver; this will enable us to assume  a solution with no phase angle present. T h e  solution is i n the form of:  <|) (7) = ,4 cosco t  (4.31)  Substituting equation 4.31 into equation 4.30 and using the expansions:  cos" co t = ^ y [ 2 c o s ( l l c o 0 + 22cos(9co 0 + H0cos(7co 0 + 330cos(5co /) + 660 cos(3co 0 + 924 cosco t] cos(co t + 9) = cosco t cosG - sinco t sinG  51  we obtain: i  a  (— -co  11 1  0  )Y4COSCO  C  t+—A C  1  —f7-[2cos(llco /)+  924cosco  t\--—AGi  2  sinco  t  i?C  = £ ( c o s c o rcosG - sinco f sinG) Neglecting the terms cos3co t and the higher orders, by the Equivalent Ritz Method, we obtain the coefficient of sinco t and cosco / as follows: coscor:  (-^-co )A + ~A 2  C  V  U  C  J  -TT924 = £ c o s 0 2"  ^4 co - — — = -E sinG  sinco t:  RC  T o eliminate the term in G , we square each of the above equations and add them to obtain:  a C  , , ,  (--co ) ^ 2  2  2  Aa b + RC C 2  2  2  2  2  2  „ (924)  \,  2  2 2  2  22  1848 b a +^rr-(--co  2" C C  ,  )^  „  = £  2  (4.32)  In equation 4.32, the values of the parameters a, b, R, C and co can be fixed and the value of the excitation can be varied slowly against the amplitude of the flux linkage. Figure 4.6 shows the amplitude of the response (A) as a function of amplitude of the excitation (E) for the following parameter values:  R=0A%, C for  100 k m . ,  a= 2.8 x 10"  3  ,  b = 7.2x  10" and 3  co = 1.0 p u .  Figure 4.6 shows the multivaluedness of the response curve due to the nonlinearity and the concept of the jump phenomenon.  Beyond a critical value of the source voltage, a jump  52  occurs, leading to the occurrence of two stable attracting points and one region of unstable saddle points. This is known as the cyclic-fold bifurcation and is associated with jump phenomenon or hysteresis where, for a limited parameter range, two stable limit cycles coexist separated by an unstable limit cycle.  2  R  1.8J-  A m p l i t u d e of the s o u r c e v o l t a g e E (pu)  F i g u r e 4.6: A m p l i t u d e of the response as a f u n c t i o n o f a m p l i t u d e of the e x c i t a t i o n  53  5  Chapter  SIMULATION A N DCATEGORIZATION OF  5.1  FERRORESONANCE  Introduction  This chapter provides an in-depth analysis of the behavior of the ferroresonant circuit of Figure 3.5.  In order to find all modes of behavior for the ferroresonant circuit of Figure  3.5, a wide range of solutions for equation 3.12 must be calculated by varying the values of E, R and C . Finding the critical values of these parameters will make it possible to predict which parameter values are dangerous in a given configuration and to allow for a safety margin with respect to these values in order to operate a network without risk. T h e fourthorder Runge-Kutta method of integration will be used as the primary tool to integrate the system's differential equation. F o r solving the nonlinear equation, the I N S I T E package was used [32]. I N S I T E consists of interactive and graphically based programs for the simulation of nonlinear dynamical systems. T h e concept of transient chaos when the transformer core losses are small is described. T h e dependency of the model represented by equation 3.12 for describing  electromagnetic  phenomena  in ferromagnetic  materials,  on  an  accurate  description of the magnetization curve, is presented. T h e solution of the nonlinear equation is shown to be dependent on the value of its initial conditions. W i t h detailed analysis of many simulation results, the basins of attraction for different chaotic regions of the system will be obtained.  54  5.2 Varying the Magnitude of the Source Voltage F o r the first set of simulations, the effect of varying the magnitude of the source voltage on the behavior of equation 3.12 was considered. T h e length of the transmission line was fixed at  100km, the transformer core losses were set at 1% of the rated transformer capacity and  the initial conditions were kept constant at < | > (0) = 0.0  and < | > (0) = 4l  pu. This initial  condition was chosen (maximum voltage and zero flux) because it is the one that most frequently occurs  in a circuit breaker when the current extinguishes. If the current goes  through zero (flux zero) and the circuit load is inductive, the voltage is at a maximum, and for a  1.0 pu. voltage, the peak value is at V2pu. In order to predict which modes of  ferroresonance may occur for a wide range of the magnitude of the source voltage, a bifurcation diagram is used. The bifurcation diagram is a plot of the magnitude taken from a family of Poincare sections  vs. the system parameter that is being varied. Figure 5.1  shows the bifurcation diagram when the varying parameter is the magnitude of the source voltage.  55  4. e a a  2 .  ~l  1  r  1  T  r  1  1  T  r  r  sea  1 . 0 0 0  xCl  J  H  1  1 K  I  I  I  Ii  H  1  1-  I  I  -a.300  J -3.000 0-000  I  J 1.87S  I  L  J  3.730 bifuroatlon  j > A m H « tmr-1  L J 3.62S  L.  E  F i g u r e 5.1: B i f u r c a t i o n d i a g r a m w i t h v a r y i n g £  The above bifurcation diagram shows what modes of ferroresonance are possible for a given value of the source voltage. The single-valued area of the bifurcation map indicates period one, the dual-valued area indicates period two, etc.  By changing one of the  parameters, one looks for a pattern of periodic responses. One characteristic precursor to chaotic motion is the appearance of subharmonics. In this type of bifurcation, a stable limit cycle loses its stability, while another closed orbit is born whose period is twice the period of the original cycle. In many physical examples a sequence of n period-doubling bifurcations is observed, in which a stable limit cycle with period 2"T is finally obtained. If the sequence is infinite, with a finite accumulation point, the resulting motion beyond this point is chaotic. The period-doubling route to chaos has been observed experimentally in all branches of physics, chemistry and biology [18]. In Figure 5.1, the region with many scattered points that looks blurred indicates chaotic ferroresonance. Figure 5.1 shows several regions of chaos interconnected by single curves that bifurcate before entering the  56  next chaotic section. The starting point of the bifurcation diagram is characterized by the flux under normal operating conditions, since the magnitude of the source voltage is about one per unit. With increasing forcing voltage, the flux rises and is driven into the saturation region. Further increments in the forcing voltage lead to a decaying bifurcation curve that reaches the first chaotic region at E = 1.9 pu. In order to find the critical values of E where the type of solution changes, E was varied from 0.15 pu. to 3.5 pu. Up to a voltage of 1.76 pu., the response was a period-one motion at the fundamental frequency. Figure 5.2 shows the phase plane diagram (a), Poincare map (b), and the frequency-power spectrum (c) for this period-one region. In these diagrams x[l] denotes the flux and x[2] indicates the voltage. For the frequency-power spectrum, with co =1.0pu., the fundamental frequency is /  fun  = — = Q.\59 pu. In the frequency-power spectrum, there is also small amounts of the  third (f  lrd  = QAllpu.) and fifth (f  5lh  - 2 . 0 0 0  -1.080  = 0.795/?w.) harmonics.  0.000  1.000  x t l l  (a) Phase-plane diagram for period-one motion  57  2.000  T—i  2. oeo  0. 000  1 *  H  I  I  I  I  I  I  -i  I  I I  1 1  H  1 h*  I  J  0.  1—i—i—i—r  I  I  1  1  I-  I  I  L  000  X[ll  (b) P o i n c a r e m a p f o r p e r i o d - o n e m o t i o n  1  0.  1  1  1  1  1  1  1  _l  I  I LX  1  |~  1  1  1  1  I  L  p  I  i  r  X88  0.  000  i  -f  0.2S0  freq  L  0.500  J  (Hz), 1 B 2 4 - P O i n t  0.750  1.000  FFI  (c) F r e q u e n c y p o w e r s p e c t r u m f o r p e r i o d one m o t i o n F i g u r e 5.2: P e r i o d - o n e r e g i o n of E =0;15  t o 1.76 p u . , R = 1%, C f o r 100 k m  a) Phase-plane d i a g r a m b) P o i n c a r e m a p c) F r e q u e n c y p o w e r s p e c t r u m  A t the point when E was increased to 1.77 pu. a bifurcation occurs and a period-two waveform was observed. Figure 5.3 shows both the Poincare map and the phase plane trajectory of this period doubling, i n which the two points i n the Poincare map represent two distinct frequencies in the output waveform. 58  ~ i — i  x.  r  i  1  1  r  1  1  1  I-  I  I  I  T  i  i  1  1  1  1  1  1  1  1-  I  I  I  1 — r  aee  H  1  1  1  J  1  1  1  -+  xC21  a  1  -l.eee  I  I  e.aaa  I  I  I  I  i.eaa  I  1  ;  xCl j  (a) Poincare map for period-two motion  (b) Phase-plane diagram for period-two motion Figure 5.3: Period-two region of E = 1.77 pu. a) Poincare map b) Phase plane diagram  A t the point when E was increased to 1.87 pu. a second bifurcation occurs and a period four waveform was observed. Figure 5.4 shows both the Poincare map (a), and the phase plane trajectory (b), of this period four waveform.  59  1.  0 0 0  a.eee  x [ 2 ]  -1.aea  - 2 . e e a  - l . a a a  a.eae  l . o e e  2 . 0 0 0  x t l l  (b) Phase plane diagram for period-four motion Figure 5.4: Period-four region of E = 1.87 pu.  a) Poincare map b) Phase-plane diagram  Further small increments in the forcing voltage lead to a decaying bifurcation curve that reaches the first chaotic region at E =1.9 pu. Figure 5.5 shows the Poincare map of this chaotic region.  60  3 . 000 1  1  1  1  1  1  i  i  i  |  i  i  i  i  |  i  i—i  i -  ^ :  XC21 - 1 . SC •• ^ r i  i  i  i  1  i  1  1  1  4  ^ ( V-  1  1  1  1  1  1  1  1  1  1  0. S 0 0 XCU  Figure 5.5: Poincare map of chaotic solution for £ = 1 . 9 p u .  The Lyapunov exponents X and the dimension d, for this chaotic region were calculated as: X, =0.090063,  X =-0.160491 2  d,= 2.562  From chapter 2, we know that a Lyapunov exponent is a diagnostic tool for deciding whether or not a system is chaotic. The dependency of a chaotic system on initial conditions means that if two trajectories start close to each other in phase space, after a small time, on the average, they will move exponentially away from each other. Since a positive Lyapunov exponent indicates expansion, what distinguishes chaotic motion from non-chaotic motion is the existence of at least one positive Lyapunov exponent. The magnitude of the exponent reflects the time scale on which system dynamics become unpredictable. The calculation of the Lyapunov dimension d, can also be used to classify strange attractors from nonstrange attractors. For nonstrange attractors d, is an integer, 61  and all chaotic attractors have noninteger dimensions [18]. Note that for this example, X, > 0 and d, is non-integer.  Another visible characteristic of a chaotic signal is its power spectrum. Figure 5.6 shows the frequency power spectrum for the above chaotic  region. The appearance of a broad  spectrum of frequencies in the output when the input is a single-frequency harmonic is a clue to detecting chaotic motion.  n 1 — r  "I  0 . 2 5 0  i  i  i  i  0 . 5 0 0 f r e q  ( H z ) ,  1 0 2 4 - p o i n t  1  1  1  1  r  0 . 7 5 0 F F T  Figure 5.6: Frequency power spectrum for E = 1.9 pu.  Table 5.1 shows the complete response of the system behavior as the magnitude of the source voltage is increased.  62  Table 5.1: Varying the magnitude of the source voltage (R = l%,  CforlOOkm)  Response (steady state)  E(pu)  0.15 1.77 1.87  -> ->  1.875  1.76  Period one  1.86  Period two (bifurcation)  1.874  Period four  1.89  Period eight  1.90  ->  1.91  Chaotic  1.92  ->  1.93  Period five  ->  1.949  Period three  1.95  ->  2.14  Period one  2.15  ->  2.16  Period two  2.165  ->  2.169  Period four  2.17  ->  2.30  Chaotic  2.90  Period one  3.10  Period two  up  Chaotic  1.94  2.31 2.91 3.15  ->  5.3 The Effect of the Transformer Core Losses F o r the second set of simulations, the effect of varying the transformer core losses R on the behavior of equation 3.12 was studied. In this case, the transmission line was taken as 100 63  km and the equivalent source  voltage was fixed at  E = 0.15pu. (corresponding to  V = \.0pu at the supply voltage). The initial conditions were kept constant at < | ) (0) = 0.0 and (j>(0) = V2. pu. With the transformer core losses at 1% (R = 4SAkQ), the ferroresonant steady-state response was of period one. As the transformer losses were decreased the ferroresonant responses of period two and four were observed. The first period-doubling bifurcation occurred at R = 0.005%. Figure 5.7 shows the Poincare map and the phase-plane diagram of this period-two response when the transformer losses are at 0.005% (R = 9.68MQ).  T—i—i—i—|—i—i—i—i—T—i—i—i—i—|—n—i—i—i— •  x t 2  _o  H  1  1  1  1  1  1  1  1  J  I  I  I  I  I  I  I  I  —  1  1  1  1  1  1  1  1  1  I  I  I  I  I  I  I  I  1  J  n  a  n  - 2 . e e a  - 1 . e e a  I  a . e e a x t l ]  (a) Poincare map for period-two motion  64  1 . a e e  2 . e e a  -rri i | i i i i  -i i i i | i i i i  0 . 0 0 0  i  i  i  1  i  i  i  i  •  i  i  i  1  i  i  i  i  i  0 . 0 0 0 x C l  J  (b) Phase-plane d i a g r a m f o r p e r i o d - t w o m o t i o n  F i g u r e 5.7: P e r i o d - t w o r e g i o n R = 0.005% a) Poincar6 map b) Phase-plane diagram  A further reduction of core losses to a value of 0.001%, resulted i n a period-four response. Figure 5.8 shows the Poincare map and the phase-plane diagram of this period-four region.  xr.2  T—i  1  1  1  -i—i  1  1  —i  1  1-  i  i  i  i  i  1  1  1—r  "i  1  -H  1  1  1—i  1  1  1  1  1  I  I  J  J  i  i  i  J 0 .  I  0 0 0  XC11  (a) P o i n c a r e m a p f o r p e r i o d - f o u r m o t i o n  65  H  1—r  1  1-  - 2 . 0 0 0  - X . 0 0 0  0 . 0 0 0  X . 0 0 0  2 . 0 0 0  x t l l  (b) Phase plane d i a g r a m f o r p e r i o d - f o u r m o t i o n  F i g u r e 5.8: P e r i o d - f o u r r e g i o n R = 0 . 0 0 1 % a) Poincar6 map b) Phase-plane diagram  Reducing the transformer core losses below 0.0006% (R = 80.67MQ) introduced chaotic behavior. Any further decrement kept the response in the stable chaotic region. Figure 5.9 shows the Poincare map of the chaotic region when the losses are at 0.0005% and the value of the driving voltage E = 0.15pu. The Lyapunov exponents k and the dimension d, of the chaotic region were calculated as:  X = 0.0831815 x  X = -0.0831817 2  d, = 2.99  We recall that the Lyapunov exponents measure the rate at which trajectories on the attractor diverge from one another and trajectories off the attractor converge toward the attractor. Therefore, for a chaotic system with low losses, where the divergence occurs very slowly, the calculation of a dimension for a fractal attractor based oh the Lyapunov exponents will result in a value for d, which is almost an integer. In fact, for a system  66  without any losses (conservative system), the sum of the Lyapunov exponents must be zero [33]. In a three-dimensional chaotic system, with Lyapunov exponents A., > 0 > A , , we 2  know that with even a small amount of damping A,, + A, < 0 , from which it follows that 2  2 < d, < 3.  1 .  0 0 0  — i — i — i — i — ( — i — i — i — i  ;—1—1—1—I—1—1—1—1  -_  •  - ' *'  " " *  " '  " '  '  I'.  1  '  "  ••... 1  '  ^ • . . • r ' . O "\  '  .  '.  1  V  •'- -Ii-:,  1  1  '  1  ' .' I  .~  x  1 • '•  \  * •  1  *  •  "i  h  - <  1  x I 2 1  ,  - 2 . 0 0 0  1  1  1  1  1  1  1  1 0.  - 1 . 0 0 0  '  1  1  1  - 1 — 1 — 1 —  0 0 0  xCll  Figure 5.9: Poincare map of the chaotic solution for R = 96.8  MO  Table 5.2 shows the response of the system behavior as the value of the transformer core loss R is decreased.  Table 5.2: Varying the value of the transformer core loss  (E = 0.\5pu., C for 100km) Response (steady state)  R in % 1  ->  0.006  Period one  67  0.005 0.001 0.0005  -> -> ->  0.002  Period two (bifurcation)  0.0006  Period four  below  Chaotic  A s the value of the transformer core losses was decreased, a smaller source voltage was needed to drive the system into the chaotic region. F o r example as shown i n section 5.1, i n order for the system to reach the first chaotic region for a 1% core loss,  E has to be  increased to a value of 1.9 pu. However, if the losses are now reduced to 0.1%, the first chaotic region appears sooner, at E = 1.65 pu. A general study was done i n which both parameters R and E were varied simultaneously. Table 5.3 shows the value of the source voltage for which the first period-doubling bifurcation and the first chaotic region appears, as the value of the transformer core loss is decreased.  Table 5.3: Simultaneous variation of core losses and source voltage  C for \QQkm First bifurcation appears at  First chaotic region appears at  R in %  E(pu.)  E{pu.)  1.0  1.77  1.90  0.5  1.75  1.86  0.1  1.51  1.65  0..05  1.46  1.6  0.01  0.45  0.6  68  0.005  0.15  0.25  0.001  0.15  0.18  0.0005  -  0.15  The above result indicates the importance of the inclusion of the transformer core losses i n ferroresonance studies. Even small losses can make the difference between the chaotic behavior and ordinary ferroresonance or periodic operation. Therefore, the recent trend of decreasing the losses of power transformers increases the possibility of chaotic behavior.  In the ferroresonant circuit of Figure 3.5, if the transformer losses are considered small enough (R = 1 x 1 0 Q ) so it can be considered virtually lossless, the system can then no 20  longer be considered a dissipative system. Consequently, the ferroresonant circuit with lossless power transformers should be considered as a conservative system. In practice, many real systems are nearly conservative. T h e most famous conservative system is the solar system, and it was the consideration of the dynamics of the solar system that led Poincare to introduce many of the methods already described i n Chapter 2 for dealing with nonlinear dynamics.  What distinguishes a conservative system from a dissipative  system is that there are no transient motions i n the dynamics which are followed by limiting motions. Therefore, there are no attractors i n conservative systems. Each initial condition results i n a unique orbit, which may be periodic, quasi-periodic or chaotic. In conservative systems, the chaotic motion does not have the kind of fractal structure that one finds in dissipative systems. A chaotic orbit i n conservative systems will visit all parts  69  of a subspace of the phase space uniformly; that is, the orbit exhibits a uniform probability density over restricted regions in the phase space. This lack of attractors in conservative systems means that trajectories starting with different initial conditions may behave quite differently as time goes on; there is no common attractor onto which they settle. The Lyapunov exponent test can be used for both dissipative and non-dissipative systems, whereas the fractal dimension measured by Lyapunov dimension test only makes sense for dissipative systems. Figure 5.10 shows the Poincare map of the chaotic region for the ferroresonant system of Figure 3.5, when the transformer core losses are considered negligible. In this case R = 1 x 10 Q, 20  C for 100km, E = 3.0pu and the initial conditions  are § (0) = 0.0 and <J) (0) = V2 pu. The same Poincare map will be obtained if the core losses are ignored ("rr; = 0.0).  RC i  -i—i—i—r  -1—I  1  H  I  l_  1 1—r  ~~i  i  i  r  I  I  x C 2 J  -2.  _1 s e a  I  '  i  I  - 1 . 2 5 0  1_ .  0 0 0  I 1 .  I  2 5 0  Figure 5.10: Poincare map of the chaotic solution for conservative system  70  1 1  L  As  can be seen the above Poincare map, the chaotic motion for the lossless system does  not have any fractal structure, there are no attracting regions in phase space and the points visit all parts of the subspace of the phase space uniformly. F o r this reason there is no basin of attraction for a conservative system. T h e set of initial conditions for which the motion tends toward a given attractor is called a basin of attraction [33]. each set of initial conditions results in a unique orbit. This  In a conservative system dependency is shown by  keeping the same parameters used for the chaotic region of Figure 5.10 constant, and  changing the initial conditions to <j) (0) = —1.5 and <j) (0) = 2.5 pu. This time, instead of a chaotic solution, the response is quasi-periodic. A quasi-periodic solution is a dynamic solution characterized by two or more incommensurate frequencies. In fact, if the simulation is started with any initial conditions chosen from outside the bounded region of the Poincare map of Figure 5.10, the response of the system will be quasi-periodic. Figure 5.11 shows the Poincare map of the quasi-periodic response of the conservative system.  71  1  1  1  1  1  1  1  1  1  i i  i  i  i i  i i  i i  * i  i  i  i  ,.|  «..._ i  i  i  « —  f  j  i i  i  _  i  i' \  \ \ \  —  i — i — j — i  — i — i  i — i — i — i  V ,  "* — i  .i...  i  i  _  1  i  i  i  i a.  eea  x t l  J  F i g u r e 5.11: P o i n c a r e m a p of the quasi-periodic s o l u t i o n f o r the conservative system  For  the same system parameters of Figure 5.10, if the clamping factor is increased, the  Poincare map in the chaotic region exhibits a much different shape with a fractal structure. Figure 5.12 shows the Poincare map of the chaotic region when R = 0.001% and E  = 3.0 pu. i  2.  r  i  i  i  r  ~1  1  1  1  J  I  I  I  1  1  J—  s e a  m  J  I  l_  J  I  I  L e.  _l  eea  K C l ]  F i g u r e 5.12: P o i n c a r e m a p o f the chaotic s o l u t i o n f o r the dissipative system  72  I  I  1_  5.4 Transient Chaos In some situations, when the losses in the system are small, chaotic behavior may appear for some parameter changes but after some time, the solution settles into a periodic or quasiperiodic mode. Such transient chaos might be a consequence of a sudden disappearance of sustained chaotic dynamic. Therefore, numerical simulations should be allowed to run for some time after one thinks the system is in chaos, even if the Poincare map seems to be mapping out the fractal structure characteristic of a strange attractor. This was evident as the losses of the system in this study were reduced below 0.01%. In many instances, after the transient solution was taken out, it was evident that the system was moving toward a period-one solution. F o r the same reason, obtaining a bifurcation diagram by varying parameter R , when R is very small, will become practically impossible. A n example of the transient chaos phenomena  E = 033pu., the  is shown in Figure 5.13. F o r  R-  0.001%,  C for \00km  and  Poincare map starts exhibiting a chaotic region which is indistinguishable  from the Poincare map of the steady-state chaotic behavior of region shown in Figure 5.9. This behavior may possibly extend for a long time, depending on the initial conditions, but the orbit will eventually move away from the region of the chaotic attractor and diverge to a period five attractor. Therefore, judging the outcome of the system behavior too early based on the Poincare map, may result in a wrong prediction. A s to how long a system must be simulated in order to be sure the response is in a steady state is a difficult question to answer. Different dynamical systems seem to behave quite differently depending on the values of both the parameters and the initial conditions. This complication arises because,  73  in general, a dynamical system may have several attractors that coexist for a given range of parameter values. The system can change its behavior because the attractors have different characters or the basins of attraction can interact in such a way as to give rise to a periodone, period-five or a chaotic response. For example in Figure 5.14 the Poincare map initially indicates a chaotic response; however if the simulation is run for a long time, the system will settle  into a period-one response. This phenomenon, where the initial  dynamics of the system response appears to be chaotic, is called transient chaos [18]. The time length of a chaotic transient depends sensitively on the initial conditions.  1  T  |  1  1  1  1  1  •'  1  1  • LL~I! ' • • '  1  1  1  • - • -  <••'••  —• >'  ± .  0 0 0  1  I  u-:—i  TTT  1  | 1 1  1 1 | 1 M 1  1  TI  e . e e e  1  1 1  1  0 .  0 0 0  x C l l  (a) Poincare map after 3 sec. simulation time  74  2 .  0 0 0  1.250  0 . 0 0 0  1 1 1 1 1  —i  •i  i  i  |  i  i  i  i  1  see  n i i | i i M• 1 1 1  2.  II : 3 -.8  -2.S00 -2.000  i  i  i  i  1  i  1  1  1  1  1  1  I—i—  i  1  1  I  0. 000  -1.000  xtl J  (b) Poincare map after 6 sec. simulation time  2.S00  1  1 1 1 1  1  1  1  1  i  1  11  1  1  1 1  ~=M —-jtrr  1  1  1  1  1  111  1  i  i  i  0. 000 K i l l  (c) Poincare map after 12 sec. simulation time  75  i  1  i  i  i  1  —  :  i  r  1  1—i  e.eee  J  I  i_  i  1—r  -i  1  1  H  _i  I  I  1_  i—i  _1  1 r  -H  1  I  i_  "i—i—r  1  1  1  J  e.eeo  1  »  H  1  I  l_  (d) P o i n c a r e m a p after 30 sec. s i m u l a t i o n t i m e  F i g u r e 5.13: P o i n c a r e m a p o f the transient chaos settling i n t o a p e r i o d - f i v e response R = 0 . 0 0 1 % , C for \00km,E = 0.33pu.  T—i—i—i—I—i—\—i—i—I—i  i  i  r  I  1  L  0.008 xt21  J -4.eea -e.see  I  L  I  I  I  0.125  I  I  I  I  I  I  0.750 xCl i  (a) P o i n c a r e m a p o f transient chaotic s o l u t i o n  76  1  I  4.  080  i  r  H  1  1-  J  1  L  I  i—i—i—r  "i—i—i—r  ~i—r  2 . 0 0 0  0 . 0 0 0  H  1  1  (-  I  L  H  1  h-  I  1_  xC2]  - 2 . 0 0 0  - 4 . 0 0 0  —  -0.500  J 0.12S  I  I  L 0 . 7 5 0  J  1.375  J  1  2 . 0 0 0  x t l ]  (b) Poincare map after the transient solution was removed  Figure 5.14: Poincare map of transient chaos, settling into a period-one response R = 0.1 %,C for \00km,E = 2.3\pu.  For the simulations of Figure 5.13 and 5.14 in order to make sure a type of behavior called intermittency [33] dose not occur, the simulations were run for a long time even after the transient solution was removed. For example for the Figure 5.13, the simulation was run to 90 seconds and the result was still a period five motion. Intermittency occurs whenever the behavior of a system seems to switch back and forth intermittently between apparently periodic and chaotic behavior. In these type of systems the behavior is predominantly periodic for some control parameter value with occasional bursts of chaotic behavior. As the value of the control parameter is changed, the time being spent chaotic increases and the time being spent periodic decreases until, eventually the behavior is chaotic all the time.  77  5.5 Varying the Length of the Transmission Line It has been shown i n previous sections that the variation i n applied voltage and core losses can result i n different modes of ferroresonance and different peak voltages. B y increasing the applied voltage or by decreasing the core losses, the ferroresonant circuit of Figure 3.5 w i l l go through a period-doubling bifurcation route to chaos. In this section the effect of varying the series capacitance C, which represents the length of the transmission line on the behavior of the equation 3.12, was studied. In this study the length of the line was varied f r o m 10 k m to 2000 k m . The value of C not only affects the damping factor,  but it has also a direct effect o n the nonlinear term, ~^X ^ b§ a  +  11  YRC'  ) • F o r this reason the  bifurcation map is used i n order to find the critical length of transmission line for which the system exhibits chaotic behavior, for different sets of values for E and R .  F o r the transformer core losses set at 1% and the magnitude of the source voltage fixed at any value i n the range f r o m 0.15 p u . to 1.7 p u . , varying the length of the line f r o m 10 k m to 2000 k m has no effect on the response of the system, and the system response remains i n a periodic motion w i t h a fundamental frequency. Figure 5.15 shows the bifurcation map when the varying parameter is the length of the line i n k m , and R = 1%, E - 1.7pu. The single-value area of the bifurcation map indicates a period-one response for all values of line length.  78  2.  see  T—i—i—i—|—i—i—i—i—|  r  1  1  I  1  1.25B  0.000  ^—i  I  1—i—1—H  • • •  -1.250  -2.500  e.eae  I  i  i  i  i  I—i—i—i—i—L  0.500 1.000 1.500 bifurcation parameter-: C <1E3)  2.000  F i g u r e 5.15: B i f u r c a t i o n m a p f o r a v a r y i n g l e n g t h o f l i n e i n k m . ,  E = X.lpu., R = 1%.  Increasing the value of the forcing function makes the length of the line an important factor. This is evident from the bifurcation diagram of Figure 5.16, in which E is increased to 1.85pu. The critical values of the line length for which the response of the system is chaotic are between 40 to 100 km, and at 1090 km. As figure 5.16 shows, there are regions of periodic behavior between regions of chaotic behavior. Such periodic windows occur frequently in chaotic systems. In Figure 5.16, the response start in a period one region and as the length of the line is increased the first period-doubling bifurcation occurs at 20 km. The response goes into chaotic region at 40 km, and comes out of it into a period-two response at 110 km. At the length of 120 km, the response will again become a period-one response. A further increase in the magnitude of the source voltage results in more values of line length for which the system response goes into the chaotic region. This is shown in the  79  b i f u r c a t i o n m a p o f F i g u r e 5.17 i n w h i c h t h e source  E = 25pu.,  voltage is i n c r e a s e d t o  R = \%.  2.500  1  1  1  1  1  1  1  1  1  1  1  I  1  I  |  1  1  1  R  1.250  0.000  -IH—i—i—i—|—i—i—i—i—|—H—i—i—i—|—•—'—'—>-  xlll -1.2S0  i -2.500 0.000  i  i I I i i I I I I I 1 1 1 0.500 1.000 1.500 bifurcation parameter: C <1E3>  I  F i g u r e 5.16: B i f u r c a t i o n m a p f o r a v a r y i n g l e n g t h o f l i n e i n k m . , E =  3.000 r - — | — i  1  1  1  1  1  1  1  1  1  I  I  I  |  I  1 1  L 2.000  \.S5pu.  R = 1%  I  1  I ^  ... . ::'-"'!;i'.l:i !!•:*'• :ii  1.500  0.000  i  , i.  ;  " H — |  • i—i—\.\; ;fV ;  ;  i—i—i—i—I—i—i—'—  Kill -1.500  -3.000  j B0  i  I I I I I I I I i I — i L 0.500 1.000 1.500 bifurcation parameter: C <1E3>  i — i — i — L .  F i g u r e 5.17: B i f u r c a t i o n m a p o f v a r y i n g the l e n g t h o f t h e line i n k m .  80  2.000  With transformer core losses reduced to 0.1%, the value of E had only to be increased to 1.0 pu., in order for the first chaotic region to appear at a line length of 400 km. Increasing the line length from 400 km to 640 km kept the response in this chaotic region. For the line length of 650 km to 1200 km the system response was of period two. Increasing the length of line further brings the system back to a period-one response for which continues up to a length of 2000 km. Figure 5.18 shows the bifurcation map for varying length of line with R = 0.1% and E = l.Opu.  2 •  5B8  I  1.258  —  8. 888  '  xin  "  -1.258  —  -2.5B0  '  e.eea  1  1  1  1  1  1  1  1  1  1  1  I  I  I  |  I  I  I  I  —  1  1  1  1  '['• | ' V ; — I  1  1  1  1  •  1  1  1  1  1  1  1  I--:!.!:  —  I  I  I  I  I  I  I  I  a.see  I  I 1.000  1  1  1  1  1  1.500  1  1  1  1 2.000  bifurcation parameter: C <1E3>  Figure 5.18: Bifurcation map of varying the length of the line in km., E = l.Opu. R = 0.1%  Increasing the magnitude of the source voltage above 1.0 pu. causes the system response to go into the chaotic region at lower lengths of the transmission line. Figure 5.19 shows the bifurcation map of varying length of the line when E = 1.6pu., in which the first chaotic region appears at 40 km. Further increases in the source voltage cause the response to go  81  chaotic at any length of line above 500 k m . This is shown i n bifurcation map of Figure 5.20.  T — i — i — i — | — i — i — i — i — | — i — i — i — i — | — i — i — i — r .i'  i". V i;  1.250  e.eee  1  1  #••>:•'*) •. ti  1  1  1  1  1  1  1  1  1  1  1  1 h  1  xCll -1.250  J B.BB8  I  I  1  I  I  I  I  0 . 5 0 0  I  I  I  I  I  1 . 0 0 0  bifurcation parameter: C  I  I  I  I  I.SBB  e.saa  i.eee  I  J 2 . 0 0 0  (1E3)  F i g u r e 5. 19: B i f u r c a t i o n m a p f o r v a r y i n g l e n g t h o f l i n e i n k m . , E = \.6pu.  e.eee  I  R = 0.1%  l.see  2 . 0 0 0  bifurcation parameter: C <1E3>  F i g u r e 5.20: B i f u r c a t i o n m a p f o r a v a r y i n g l e n g t h o f the l i n e i n k m . , E = 2.0pu.  82  R - 0.1 %  The above variation of the parameter  C shows that the length of the transmission line  becomes an important factor i n avoiding the chaotic regions only when the source voltages are high or when the transformer core losses are low. These results indicate that, when the length of  line increases, the damping factor decreases and the possibility of chaos i n  ferroresonance rapidly increases. However, at the same time, by increasing the value of  the effect of the nonlinear term ~^( § a  +b§  ")  l s  C,  reduced. A n d that is w h y in the range of  lower voltages, for higher lengths of the transmission line, the system response is periodic.  5.6 Dependence on the Description of the Magnetization Curve T h e performance of the mathematical model represented by equation 3.12 for describing electromagnetic phenomena in ferromagnetic materials depends on the accurate description of the magnetization curve. In the peak current range for steady-state operation, the fluxcurrent linkage can be approximated by the first term of the equation 3.8 which is:  / = a<j) + b § "  where the coefficient a of the linear term corresponds closely to the reciprocal of the inductance (a « YL). saturation  where  the  However, for very high currents the iron core might be driven into flux-current  characteristic  becomes  highly  nonlinear. Using  a  polynomial of order less than eleven (n < 11) to represent the magnetization curve might be adequate for small-capacity transformers, but it does not bend sharply enough at the knee  point  to  describe  the  magnetization 83  characteristic  of  modern  high-capacity  transformers. Figure 5.21 shows the details of the saturation region for different exponents of the nonlinear term in equation 3.8. In this section the results of an investigation into the effect of changing the exponents of the nonlinear term on the modes of behavior of the ferroresonant circuit are presented.  Current / (p.u.) F i g u r e 5.21: C u r r e n t l i n k a g e w i t h v a r y i n g exponent  In section 5.1, Figure 5.1 shows the bifurcation diagram when the order n of  the  nonlinearity is taken as eleven and the varying parameter is the magnitude of the source voltage. The diagram shows what modes of ferroresonance are possible for a given value of the source voltage. In that diagram, an increase in the forcing function leads to a decaying bifurcation curve that reaches the first chaotic region at E = 1.9pu. Figure 5.22 is the same bifurcation diagram as in Figure 5.1, except that the order n of the nonlinearity is taken as five instead of eleven.  84  Comparing the two bifurcation diagrams it becomes apparent that the first chaotic region begins much sooner for the higher exponent n . For n = 11, the first chaotic region starts at E = 1.9pu. whereas for n =5 chaos first appears at E = 4.0pu. It was found that the higher the order n, the lower the magnitude of the source voltage at which chaos first occurs. However, when the order of the nonlinearity is above nine the first chaotic region always starts at about E = 1.9pu. This fact implies an important conclusion: the sharper the bending of the transformer magnetization characteristic, the higher the likelihood that chaotic ferroresonance will appear.  Figure 5.22: Bifurcation diagram with varying E for n =5  For improved power transmission, utilities have a high interest in developing transformers that are characterized by sharper bending around the knee point. However, this also increase the potential for realistic network conditions to approach chaotic regions.  85  In section 5.2 we have shown that with n = 11, reducing the transformer core losses below 0.0005% introduced chaotic behavior. And any further decrement kept the response in the stable chaotic region. When the order n of the nonlinearity was chosen as five instead of eleven, reducing the losses did not drive the system into the stable chaotic region. Instead, long- lasting transient chaos was observed. Figure 5.23 shows the Poincare map for the case when n = 5 and losses are 0.0005%. In this case, if the simulation is run for a long time, the system will start settles into a period-one response.  3  • QOO  1.300  1  1  1  1  —  |  1  1 1 1  1 1  1  '  ' '  1 | 1  1  1 1  '£§  -  " '"MvlJ  .. : * .  1 -2.000  1  1  1  1 1 1  -B.75B  1  1  1  1  "Jj  1 1 1 1 1 1  1  0.S00  1.7S0  1  1  1 1 3.000  xCl)  Figure 5.23: Poincare map of transient chaotic solution for n = 5,  R = 96.8 Mil  Figure 5.24 shows that as the transient solution is taken out, the system moves toward a period one solution. In fact, if the system is allowed to run for a longer time the Poincare map will show only one dot which indicates a period one response. Therefore for cases of n = 5, and n = 7, regardless of how small the losses get, chaos does not appear in the  86  solution. This result shows the importance of appropriate choice of n in representing the magnetization curve of the transformer.  T — i — r  ~i—i—i—i—I—i—i—i—r  T — r  l.saa  '  1  0.800  I '  x(21  -l.saa  -3.00a  J  I  I  L  e.saa  -a.  j  i_  1.730  xCll  Figure 5.24: Poincare map after the transient is removed  5.7 The Effect of Initial Conditions One of the most important signatures of chaotic systems is their extreme sensitivity to initial conditions, and the divergence of nearby trajectories. If a system displays divergence of nearby trajectories for some values of its parameters, then the behavior of that system becomes unpredictable. However, the system is still deterministic, and the future behavior of a trajectory could be predicted if we knew the exact value of the initial conditions for a given trajectory. But if we make the smallest change in those initial conditions, the trajectory will follow a completely different path.  87  Therefore, a small change in initial conditions leads to a large difference in long-term behavior of the system. The set of initial conditions that generate trajectories which approach a given attractor is called the basin of attraction. If there are two or more attractors, there are some initial conditions that lie on the border between the two basins of attractions. These attractors are all associated with geometric objects in phase space. The equilibrium state is associated with a point, the periodic motion is associated with a closed curve, the quasi-periodic motion is associated with a surface in a three dimensional phasespace, and the strange attractor is associated with a fractal set [34]. Figure 5.25 shows an example of these attractors and the basins of attractions. The trajectories that start inside the dotted basin settle on the attractor inside the dotted region. Those trajectories that start outside of the dotted basin head for other attractor.  88  In the previous sections the initial conditions <|) (0) and < | > (0) were kept constant and the critical values of parameters E,  R and C were determined in order for the system to  operate outside the chaotic region. In this section the results of a study into the effect of varying the initial conditions and drawing the basins of attraction for the ferroresonance circuit of Figure 3.5 are shown. T h e system described by equation 3.12  is studied for  different sets of initial conditions, and the sensitivity of the system to these changes is investigated.  F o r the first set of simulations, the effect of varying the initial conditions when the transformer core losses are small is studied. It was shown i n section 5.2 the transformer E = 0.15  core losses to 0.0005% and keeping the equivalent source voltage at  pu. with cj) (0) =  Figure 5.9 in section 5.2 Figure 5.9,  that by reducing  0,  <)) (0) = V2 ,  the response enters chaotic  ferroresonance.  shows the Poincare map of this region. F o r the chaotic region of  if all the parameters are kept the same and the initial conditions are changed to  (() (0) = 0 and <j> (0) = 1.0, the response, instead of being chaotic, settles into a periodic motion, after the long-lasting transient chaos dies out. Figure 5.26 shows the Poincare map of this periodic attractor, in which the single-valued area of the Poincare map indicates periodic-one.  89  1  1  1  1  1  1  1  1  H  1  1  1  •  1  "  '  i  i  I  i  i  i  i  T  1  1  1  '  '  '  i  i  i  1  1  1  1  I  I  '  I  ~  i . eee  1  1  *  ' -_  xC2 J  i  -i.eee  T  e.eea xtl  J  F i g u r e 5.26: P o i n c a r e m a p of the p e r i o d i c s o l u t i o n  R  = 0.0005%,  i  I  i.eee  i — i — i — i —  E = O.lSpu., C for  a.eae  lOOfcm  The above result indicates that the system has more than one attractor for a given set of parameters and different initial conditions. Depending on the choice of initial conditions, the system may settle on either a periodic or a chaotic attractor. The basin of attraction for the periodic attractor in Figure 5.27 illustrates that any combination of initial conditions inside the shaded area gives rise to a periodic attractor. All other sets of initial conditions lead to the chaotic attractor of Figure 5.9. Comparing Figures 5.27 and 5.9, it is evident that when the system of equation 3.12 is lightly damped, the set of initial conditions which give rise to a chaotic attractor are all the points that the trajectory visits in the Poincare map, plus the points that lie outside of the shaded region of the map. Another region of chaos for equation 3.12 was observed when the magnitude of the source voltage was increased to E = 1.9 pu. Figure 5.5 shows the Poincare map of this region  90  when the transformer losses are increased to 1% (R = 48.4 kQ), the length of the line is 100 km., and the initial conditions are taken as < | > (0) = 0 and < > f (0) = V2 .  <D(0)  -1.55  Figure 5.27: Shaded area is the basin of attraction for the periodic attractor R = 1%, E = 0.15p«., C for 100£m  In this simulation, if the initial conditions are changed to <(> (0) = 0 and < | ) (0) = 0.6, the system response settles on a periodic motion, similar to Figure 5.26.  Therefore depending on the set of initial conditions, the response can either be periodic or chaotic, and this shows the coexistence of a chaotic attractor and a stable periodic orbit. Figure 5.28 shows the basin of attraction for the periodic attractor, in which any initial condition chosen from inside the shaded area leads the system into a periodic motion, and all other points chosen outside of the shaded area results in a chaotic response. In this case it  91  can be seen that the basins of attraction are somewhat more complex in shape, and some large-amplitude starts lead to small-amplitude motions, and vice versa. The third region of chaos occurs when the excitation voltage is increased to a higher value of E = 5 pu. The Poincare map of this chaotic region is shown in Figure 5.29 whencb (0) =0 and <|>*(0) = -J2 .  Figure 5.29: Poincare map of the chaotic solution for R = 1%, C for 100km, E = S.Opu  92  The set of initial conditions that give rise to a periodic attractor is shown in the shaded area of Figure 5.30. All other sets of initial conditions lead to the chaotic attractor of Figure 5.29.  i(0)  «*(0)  Figure 5.30: The shaded area is the basin of attraction to the periodic attractor for R = 1%, C for 100km, E = 5.0pu  Comparing Figures 5.28 and 5.30, it is evident that the area of the basin of attraction for the periodic attractor is decreased as the forcing function is increased.  To further show this point, the excitation voltage is increased to a high value of 25 pu. The Poincare map of this chaotic region is shown in Figure 5.31, when tj)(0)=0 and < | > (0) = V2 . For this simulation, the basin of attraction for the chaotic attractor is the entire phase space (state space), and the basin of attraction for the periodic attractor has  93  totally disappeared. Therefore the area of the basin of attraction for the chaotic attractor is increased as the forcing function is increased.  For a low loss system, the basin of attraction which results in a chaotic trajectory consists of all the points that the trajectory visits in the Poincare map, plus the entire space outside the attractor. However, for normal losses and higher values of source voltage, the basin of attraction for a chaotic trajectory is somewhat more complex in shape, and has to be calculated on a point by point basis. It is evident from this study that as the forcing function is increased the basin of attraction for the chaotic attractor is increased, and at a high value of the forcing function, the basin of attraction for the chaotic attractor includes the entire phase plane.  12 .  1  1  1  1  1  1  1  1  —I  1  1  T"  1  1  i  i  i  i  |  1  1  1  1  i  i  i  i  -V*  -  ,  ,  ,  / ^ p p p  i  i  i  i  1  i  i  i  i  1  1  -1.5BB  Figure 5.31: Poincare map of the chaotic solution for R = 1%, C for 100km, E = 25.0p«  94  1  5.8 The Effect of Copper or Series Losses In a well-designed transformer, for maximum efficiency the losses that vary as the square of the current must be equal to the constant losses [35] . In a transformer, the core loss is defined as the sum of hysteresis and eddy-current losses, and its value depends o n the peak flux density produced in the core. T h e peak flux density depends upon the terminal potential difference, and as long as the voltages at the transformer primary and secondary remain within a few percent of rated value, the core loss can then be considered nearly constant. T h e sum of the primary and secondary losses are called the copper or series losses. The series losses can be determined for any condition of operation if the winding currents and their resistance are known. T h e ratio of the core losses to copper losses at rated output may be different in transformers designed for different types of service. However, if the transformer is designed for the maximum efficiency at rated current and voltage, the fullload copper losses will be equal to the core losses. During  ferroresonance the R M S current  flowing through the capacitance is a small percentage of the system-rated current, and the voltage across the transformer coil is higher than the coil's rated voltage. Therefore, the series losses are much smaller than the shunt losses and can be neglected [22]. In order to verify the above statement and to find the effect of adding the copper losses o n different modes of behavior for the ferroresonant circuit of Figure 3.5, a value of series with the capacitance C .  T h e value of R to represent the copper losses is calculated as: {  95  was added in  Pcu — R\l  — Pcore  P  no  min al  no  min al  For 1% of the rated transformer capacity, the value of R is calculated as: x  core = 0.01 x 25 MVA = 0.25 MVA  p  0.25 x l O ,  '=^W  6  ( l l 0 x l o  »=  4 l 4 f l  W i t h the copper losses set at 1% of the rated transformer capacity, the length of the line fixed at lOOkm, and the core looses set at 1%, the bifurcation diagram was obtained when the varying parameter was the source voltage. T h e bifurcation diagram obtained was exactly the same as the one shown i n Figure 5.1. This indicated that the addition of 1% of copper losses had no effect on the modes of behavior of the ferroresonant circuit. In fact, the copper losses had to be increased to 3% of the rated transformer capacity i n order for the first chaotic region to appear at E = l.95pu., instead of at E - 1.9pu. with no copper losses added. Increasing the copper losses to 6% caused the first chaotic region to appear at a value of  E - 2.0pu. Figure 5.32 shows the bifurcation diagram when the varying  parameter  is the  magnitude  of the source  Cfor \00km.  96  voltage  with  R = 6%, R = 1%, and x  -0.500  -2 . 000 0  L_I  1  1  1 1 1 1 1 1 | I I 1 I I I 2.000 4.000 6.000 bifurcation parameter: E  I  I  I  8.000  F i g u r e 5.32: B i f u r c a t i o n d i a g r a m w i t h v a r y i n g E a n d R = 6 % x  The above result shows that including the series losses is not an important factor in studying ferroresonance as compared to including of core losses. This is true until the value of copper loss approaches a very high value of about 6% of rated transformer capacity, which is not generally the case in practical distribution systems.  5.9 T h e Effect of H i g h e r - O r d e r H a r m o n i c s  It was shown in previous sections that the appearance of subharmonics plays an important role in the transition of ferroresonant circuit behavior from period-one response en route to chaos. For a sinusoidal function, there is a single peak of the spectrum at its frequency /,.  But in general, for a well-behaved periodic function, there are peaks of rapidly  diminishing amplitude at / , (fundamental frequency), 2/, (second harmonic), 3 / , (third harmonic), etc. In an ideal model with a continuum of precise measurement over an infinite interval, these peaks have the zero width and infinite height of a delta function, but noise in  97  the system, in the measurement, and in processing of the signal gives narrow peaks of finite height in practice. In order to find the influence of adding these harmonics on the chaotic motion of the ferroresonant circuit,  a third harmonic with a magnitude of 10% of the  fundamental was added to the source: d fo "df 2  +  1 di> RC dt  +  1 C * (fl  ,, +  6<t>  }  =  n ^  0  0  ^ ^  +  c  o  ,0-l£c°s(3co t) s  From the bifurcation diagram, when the varying parameter has the magnitude of the source voltage with the addition of 10% of the first harmonic, the first chaotic region appears at E = 1.8/?w. But increasing the magnitude of the third harmonics to a higher value has a very limited effect on the appearance of chaos. In fact with the addition of a third harmonic with a magnitude of 50% of that of the fundamental, the first chaotic region appears at E = \.l\pu. Figure 5.33 shows the bifurcation diagram with varying E , when a 50% third harmonic was added to the source voltage. Comparing this bifurcation diagram with Figure 5.1 where there is no 3 harmonic added, one can see that the first chaotic region appears rd  sooner.  In this chapter the results in investigating the impact of different sets of circuit parameters and different magnetization characteristics of the transformer on the behavior of  the  ferroresonant circuit have been presented. It was seen that variations in the exponent of the nonlinearity, losses, length of the transmission line, or applied voltage can result in different modes of behavior, which may be periodic or chaotic.  For example, when the applied  voltage is gradually increased, the responses may first be period one, then period two, then  98  chaotic. A more accurate representation of the < | > - / characteristic of the transformer results in chaos occurring at a much lower magnitude of the source voltage.  ~i—i—i—r  i  i i  r  T—i—i—r  >'itt V  1.G00  xr.il  -I  1 H  (.: . I  H  1  1  1 H  I  I  -0.500  0.000  J  I  I L 1.800 3.600 5.400 bifurcation parameter: E  L.  7.200  Figure 5.33: Bifurcation diagram with added 50% third harmonic  This chapter also presented results in investigating the impact of different sets of initial conditions on the modes of behavior of a ferroresonant circuit. It was shown that the ferroresonant circuit can have more than one attractor for a given set of parameters and different initial conditions. Depending on the choice of initial conditions, the system may settle on either a periodic or chaotic attractor. For a low-loss system, the basin of attraction for a chaotic trajectory consists of all the points that the trajectory visits in the Poincare map plus the entire space outside the attractor. However, for normal losses and higher values of the source voltage, the basin of attraction for a chaotic trajectory is somewhat more complex in shape and has to be calculated on a point-by-point basis. It is evident from our study that as the amplitude of the forcing function is increased, the basin of attraction  99  for the chaotic attractor is increased, and eventually at a high value of the amplitude of the forcing function, the basin of attraction for the chaotic attractor will include the entire phase plane.  In all the cases presented in this chapter, the numerical simulations were allowed to run for a long time so that transient chaos is not misinterpreted as steady-state chaos.  100  Chapter  USING EMTP T O SIMULATE  6.1  6  C H A O S IN P O W E R SYSTEMS  Introduction  A n important practical problem is the interpretation of signals from a nonlinear system. A nonlinear system can be in the form of a natural phenomenon, a laboratory experiment, or a numerical experiment. Most often the description of a system as a set of differential equations is not easily obtainable. However if we do not know the equations which govern a system, then we may use the measurements of its output to find the states of the system, their stability and bifurcations as the control parameters vary. A n advantage of the timeseries analysis is that, since the system is referred to as a "black box", it is not necessary to have the knowledge of the system's differential equations in order to gain insight into it or simply to study its behavior. In previous sections the complexity of a real ferroresonant circuit in a power system was reduced to a single-phase equivalent circuit i n order to obtain the differential equation of the circuit and to facilitate the analysis and simulation. In order to avoid this kind of simplification, this chapter describes the implementation of a ferroresonant circuit using the Electromagnetic Transient can be used to model the individual components  Program ( E M T P ) . T h e E M T P  of power systems containing any  nonlinear elements. A m o n g the various important features presented by the E M T P lies its ability to model nonlinear elements. Therefore, it also can be used to study nonlinear dynamic systems. This chapter describes a method of using E M T P to obtain Poincare maps  101  and to perform bifurcation simulations for the study of chaos. T h e influence of the nonlinear modeling i n the E M T P  on the simulation results will also be discussed i n this  chapter.  6.2  A n O v e r v i e w of E M T P  Different methods have been historically used to analyze the transient behavior of power systems. These studies were traditionally performed with an analog computer k n o w n as a Transient N e t w o r k Analyzer ( T N A ) . However, in the late 1960's, the Transient  Program  ( E M T P ) for digital computers,  developed at  Electromagnetic  Bonneville  Power  Administration by D r . H . W . D o m m e l [36], replaced T N A . The E M T P is currently one of the most widely used general purpose programs for transient analysis. E M T P has  a wide  range of circuit-modeling capabilities, from linear and nonlinear R L C circuit elements to circuit  breakers,  lightening  diodes,  thyristors,  transformers,  arrestors, transmission lines, ...etc.  dynamic  EMTP  synchronous  machines,  uses the trapezoidal rule of  integration, and in the M i c r o T r a n ® version of E M T P used i n this work, it sometimes (when discontinuities occur) uses the backward Euler rule of integration to damp numerical oscillations [37].  T h e most important element in studying the ferroresonant conditions  with E M T P is the modeling of the saturable inductor representing the core magnetization inductance. There are two ways of simulating transformer nonlinearities with the E M T P . One  way is to use a piece-wise linear representation  of the transformer's  nonlinear  characteristic. T h e other possibility is to apply the compensation method. T h e piece-wise linear method is used for nonlinear inductances with a characteristic of two slopes Z, and  102  L  2  as shown in Figure 6.1. M o r e than two slopes can be used as long as the parallel  inductances are connected from node to ground.  Figure 6.1: Piece-wise linear inductance representation  In the M i c r o T r a n ® version of the E M T P , the piece-wise linear representation forces the use of the critical damping adjustment ( C D A ) to damp the over-shoot at corner points. T h e CDA  scheme developed by J.R. Marti and J. L i n [38]  uses the trapezoidal rule of  integration for the normal part of the simulation, and temporarily changes to the backward Euler rule to go over discontinuities such as switching operation, change of slope in nonlinear elements, step inputs, etc. W i t h this scheme the network is solved in each time step with the trapezoidal rule, and as soon as a discontinuity is encountered, the network is  solved twice with the backward Euler method with a step size of  . It then resumes the  simulation with the trapezoidal rule. The two half-step backward Euler solutions will critically damp the oscillations created by the discontinuity. This two half-step backward Euler solution might also damp out the chaotic behavior.  T o have a purely trapezoidal solution, the compensation method was chosen to model the transformer nonlinearity. In the compensation method, the nonlinear inductance L(i) can  103  be  specified point  by  point  as  a piecewise  linear  characteristic  § = /(/).  interpolation was used between points, and the points that specify the (<|) ,/)  Linear  characteristic  were monotonically increasing. T h e compensation method used to represent the nonlinear inductance is also called a 'true type' representation, because the operation is always on the proper (j) —i segment. In this method the factorized admittance matrix, which includes only the linear elements, is formed once at the beginning of the simulation. A N e w t o n iteration, involving the nonlinear elements and a  reduced Thevenin equivalent of the  system, is performed at each time step to establish the correct segment for each element. More information on these models and E M T P usage in general can be found in the E M T P Theory Book [36].  F o r representation of the core losses, a separate resistance in parallel  with the nonlinear inductance is used.  A s it was shown in Chapter 2, a Poincare map can be obtained by sampling any system quantity and its derivative once each cycle of the forcing function. F o r the  ferroresonant  circuit of Figure 3.12, the voltage across the nonlinear inductance can be directly outputted by the E M T P . However, in order to obtain the derivative of the voltage, namely the flux across the nonlinear inductance, a large linear inductor must be placed in parallel with the nonlinear inductance as shown in Figure 6.2.  104  L=\O H n  Nonlinear L  Figure 6.2: Obtaining flux in the E M T P  F r o m Figure 6.2, the following equations can be written:  V  dt  dt  <t>, =<t> = LJe e  Therefore i  e  can be measured by E M T P and then multiplied by the constant L , which e  gives the flux through the nonlinear inductor.  The other practical problem is the sampling of the voltage and its derivative i n the E M T P . In order to do this, the E M T P integration time step and output frequency had to be adjusted so that their product equaled one 6 0 H Z period. T h e time step ( A r ) selected for the simulation of the circuit is 18.5185p, s, that is 900 points per 6 0 H Z cycle (16.6667 ms).  60  y  900  = 18.5185ps  105  The parameter I O U T in the data description card of E M T P must be set at 900.  This  parameter is used in the E M T P to control the amount of printed output. T h e number of data fields allocated for the input of the I O U T  in the E M T P is three; therefore  maximum number that can be used for this value is 999.  the  T h e choice of 900 points per cycle  is completely arbitrary; however the A t of around 20u. s is adequate for this study.  F o r the Poincare map, the results of the simulation were printed every nine hundredth point (1 point per cycle, or equivalently, with a sampling frequency of 60HZ) for a total of 4000 points. This number of points for drawing the Poincare map is necessary in the case of chaotic trajectories in order to capture the  fractal structure. Therefore, the maximum  simulation time for the E M T P files is calculated as:  ^max = 4000(900 x 18.5185p. s) = 66.667 s  T o obtain the data to construct the bifurcation diagram from the E M T P , the procedure is very similar to the one used for the Poincare maps, that is, points were obtained with a sampling frequency of 6 0 H Z , as one of the parameters was being varied. F o r example, a ferroresonant simulation must be run for a long time, with the value of the source voltage being gradually increased or decreased. T h e flux values are then sampled at the same rate as for Poincare sections, that is once each cycle of the sinusoidal forcing function. However, for the bifurcation diagram, the number of points ( I O U T ) was reduced to 300, since many simulations were necessary in order to obtain a meaningful bifurcation diagram and since there is no fractal structure associated with the bifurcation diagram.  106  6.3 Ferroresonant Test Case In this section the ferroresonant circuit of Figure 3.5  will be simulated using the E M T P .  The E M T P results will then be compared to the results presented in chapter 5 which were obtained by solving the differential equation of the system by fourth-order  Runge-Kutta method  of integration.  Three  INSITE,  regions  of  which uses the behavior  considered: period one, period two and chaos. A very detailed representation transformer nonlinearity was used. The (j) - i characteristic of the transformer  were of the  which was  modeled by the eleventh-order polynomial is discretized using the compensation method with 21 points. It must be pointed out that using so many points to represent the saturation curve is very unusual i n power-system simulations. The c o m m o n practice is to use a small number of points (5 to 9 points). The effect of this representation is discussed in the following sections. T h e complete E M T P input file is given in Appendix A . Figure 6.3 shows the Poincare map obtained from the E M T P simulation of the system operating in a period-one region. This is the same region as shown in Figure 5.2, in which R = 1% , C  for \00km  and  E = 0.15pu.  The two simulations give the same result.  107  -  2  -  1  0  1  2  Rux(pu) Figure 6.3: Poincare map of the E M T P simulation for a period one response  As was shown in Figure 5.3a of section 5.2, by increasing the source voltage to E = \.Spu., the response bifurcates into a period-two waveform. This result again is confirmed by simulating the case with EMTP. Figure 6.4 shows the Poincare map of the EMTP simulation when R = 1% , C for 100km  and E = \.8pu. which, as expected, shows a  period-two response.  Increasing the source voltage to E = 1.9pu. pushes the response into the chaotic region and the Poincare map of the EMTP simulation matches the one obtained from the INSITE program shown in Figure 5.5. Figure 6.5 shows the Poincare map of the EMTP simulation for the chaotic region of E = \.9pu.  108  2  3  &  8, CO  •3  >  -2  -2  0  -1  Flux(pu) Figure 6.4: Poincare map of the E M T P simulation for a period two response  u  CD  I a  Figure 6.5: Poincare map of the E M T P simulation for the chaotic response with E = \.9pu.  In section 5.3, it was shown from simulation results of Figure 5.9 that the ferroresonant circuit responses were driven into the chaotic region when the core losses were reduced below 0.0005% of the rated transformer capacity. The same simulation was done using the EMTP, and again both methods agree on the type of behavior the system presents. Figure  109  6.6 shows the Poincare map of the chaotic response for R = 96.8 MQ. when the simulation is done with the EMTP.  2  2  Figure 6.6: Poincare map of the E M T P simulation for the chaotic response for l o w losses (R = 96.8 M f l )  In Figure 5.10, it was shown that, when the transformer core losses were considered negligible, the Poincare map of the chaotic region does not have any fractal structure. The EMTP simulation of the same conditions shown in Figure 6.7 produces the same Poincare diagram, and the result of two different integration routines agree on the type of the behavior the system exhibits.  110  5  "-2.5  0  -1.25  1.25  2.5  Flux(pu) Figure 6.7: Poincare map of the E M T P simulations for the conservative system  6.4 Representation of the Transformer Nonlinear Characteristic with EMTP In section 5.6, it was shown that the performance of the mathematical model for describing electromagnetic phenomena in ferromagnetic materials depends on the accurate description of the magnetization curve. In Chapter 5, the differential equation of the system was solved using the analytical expression of the nonlinearity. The EMTP requires a discretization of the transformer nonlinear characteristic which can influence the way the system operates. For this reason when modeling the nonlinear inductance in the EMTP with the compensation method, thecp - / characteristic of the transformer which was modeled by an eleventh-order polynomial is discretized using 21 points. In order to show the importance of the number of points used in the discretization of the saturation curve, a series of simulations was done with only 11 points. Figure 6.8 shows the two saturation curves used in the EMTP , one with 11 points and the other with 21 points. The main difference between these two representations is the smoothness of the curve when bending from the  111  linear region into the saturation region. In the 11-point representation, there is an abrupt transition from the linear region into the saturation region.  400 350 300 250 1200 150 100 50 0 0  2000  4000  6000 Current  8000  10000  12000  Figure 6.8: Discretization of the saturation curve with 11 and 21 points  Figure 6.9 shows the EMTP simulation of the ferroresonant circuit in which the nonlinearity is represented by the 11 points discretization. The parameters used are the same as in Figure 6.5. In Figure 6.5 the response was in the chaotic region, but in this case it is apparent that that the chaotic behavior disappeared completely, having been filtered out by the inappropriate representation of the nonlinearity.  112  1.5 3  "ST o>  a  o > -1.5  J -0.25  I  0.5  1  1.25  1  2  Flux(pu) Figure 6.9: Poincare map of the region with E = 1.9pu. for 11-point representation  This result shows the importance of using an appropriate representation of nonlinear elements in the EMTP when chaotic systems are to be simulated. The important point to consider is the number of the points used in the bending region of the saturation curve.  6.5  Integration of the Trajectories  The purpose of an integration routine is to approximate the solution of continuous-time systems on a digital computer. These computers operate in discrete time; therefore, an integration routine models continuous-time systems with discrete-time systems. Different integration routines map the same differential equation into different discrete-time systems. Therefore, the simulation results of chaotic systems should always be verified by using two or more different integration routines. This is always a good check on the validity of the integration results. In the case of chaotic systems, the sensitive dependence on initial conditions  implies that an arbitrarily small error eventually affects the macroscopic  behavior of the system. Therefore, simulations of chaotic systems require careful 113  interpretation  and  should always  be  verified. Using  different  integration  routines  introduces different errors and the errors propagate differently. T h e validity of integration results are therefore obtained by using different integration routines which introduce different types of error. There are two different types of integration error: the local error, which is the error introduced by a single step of the integration algorithm, and the global error, which is the overall error caused by repeated application of the integration formula. The local error is divided into two categories: round-off error and truncation error. T h e round-off error depends on the number and types of arithmetic operations per step and is, therefore, independent of the integration step size. T h e only way to reduce the round-off error is to increase the floating-point representation.  This means using either a single-  precision (32 bits) or a double-precision (64 bits) representation. T h e truncation error get its name from the algorithm based on the Taylor series. If an infinite number of terms could be used i n the Taylor series, the truncation error  would be zero, but due to the  truncation of the series, an error is introduced. Therefore, the truncation error depends on the choice of the algorithm and is independent of the computer used.  F o r a given  algorithm of a given order, the local truncation error decreases with the step size. Also for a given family of algorithms and a given step size, the higher-order algorithms are more accurate. T h e global error is mainly the accumulation of the local round-off and local truncation errors. Therefore, this kind of error accumulates with each step and is inversely proportional to the step size.  114  A s we can see from the above definitions, a useful integration routine must automatically adjust the step size and the order of the integration routine in an attempt to use the largest step size that satisfies a user-supplied error tolerance. In these routines the user supplies the error tolerance, and the integration routine chooses a step size. After calculating the next step, it then estimates the error; if the error is larger than the specified error, the step size is reduced and the next step is calculated again. However, if the error is much less than the specified error, the step size is increased for the next step. In this study three different algorithms were used in order to check all the simulation results. T h e first was a fourthorder single-step algorithm with error control that is based o n the Runge-Kutta-Fehlberg algorithm [32] . This routine has a built-in error control which adjusts the step size so that the routine's idea of the local error does not exceed a user-specified tolerance. T h e second was a variable-order multi-step algorithm with error control based on the Gear algorithm [32].  T h e motivation behind the development of the Gear algorithm was the integration of  stiff systems. T h e third routine used i n this study was the trapezoidal algorithm. In section 6.2,  it was shown that the trapezoidal integration routine gave the same result as the fourth  order Runge-Kutta algorithm. Figure 6.10 shows the Poincare map of the simulation with Gear's algorithm for the same conditions shown in Figure 6.5 and Figure 5.9. different integration routines agree o n the type of  A l l three  behavior the system exhibits and  produce the same Poincare diagrams. A l l the simulation results that were presented in this study were tested using three different integration routines. In all the simulations the results obtained by using different integration routines were the same.  115  3.008  T—i—i—|—i—p—i—i—|—i—i  i  i  |  ^ V " ^ ^ ^ ^ ' ' ' . . . . .  i  -1.080I  I  I  I  I  I  -0.250  T l  I  1  I  I  0.S00 xCll  i  i  i  r  ~U : x X  I  I  I  I  1  1  L  1.250  Figure 6.10: Poincare map of the simulation with Gear algorithm E = \.9pu. R = 1%, C for lOOfcm.  A technique that was used in this study to judge the accuracy of the output of all the integration routines was to re-integrate the system using a tighter error tolerance. A similar technique was used in the case of choosing the time step. Results were then compared in order to verify that the choice of A / had no effect on the simulation results.  6.5 Case Study  In January 1995, the American Electric Power Service Corporation reported the occurrence of "questionable voltage spikes", in a very simple circuit, which did not stabilize after a relatively long simulation run time of one second. These spikes occurred in the low voltage side of a 480/2400 volt, transformer. Figure 6.11 shows the one-line diagram of the transformer with the load connected to the low-voltage terminal. These voltage spikes occurred when the transformer was operating in the nonlinear region of the saturation curve due to the shorting to the ground of switch "s".  116  BB1 LB1 _/480/2400F  19.761m//  S  \Fault  69.7Q  Load  445.6mH  12.63LIF 8.7Q  1  Figure 6.11: Power system under analysis  The source strength in this system is very weak and under steady-state conditions with switch "s" open, the load current supplied by the source is 55.0 A-peak through the low voltage winding of the transformer. W i t h the turn ratio of 1:5 , the corresponding current through the high voltage winding is 11.0 A-peak. The transformer saturation characteristic is connected i n the low voltage winding, and the current i n the nonlinear branch is 0.73 A peak for these steady-state conditions. The problem of "voltage spikes" occurs at the transformer  windings ( A l - B l ) and ( X l - Y l ) when the transformer  is operating in the  nonlinear region of its saturation curve due to the shorting to the ground of switch "s". In order to study the cause of these voltage spikes and identify this nonperiodic motion, the following checklist must be followed:  1. Observe the time history of the measured signal. 2. Examine the Fourier spectrum of the signal. 117  3. Look at the phase-plane history. 4. Take the Poincare map of the signal. 5. Obtain the Lyapunov exponent of the signal. 6. Vary the system parameters and look for the bifurcation and routes to chaos. The simulations for the study of this circuit were done with the EMTP. The complete EMTP input data file is given in Appendix A. Figure 6.12 shows the voltage vs. time across the nonlinear inductance after the switch "s" is closed; the simulation is for the first 0.1 second. The voltage diagram shows a high value of the fundamental component along with the series of integer and non-integer harmonic components that are imposed on it. In order to make sure the steady-state solution had been reached, the simulation was run for 1 second. Figure 6.13 shows the voltage diagram across the nonlinear inductance for 1 second of simulation time. 4000|  1  1  i  '  1  '  '  '  Time (second)  Figure 6.12: Voltage across the A l - B l for the first 0.1 second  118  1  1  4000r 3000 •  Time (second) F i g u r e 6.13: V o l t a g e w a v e f o r m across the A l - B l b r a n c h f o r 1 s e c o n d  The  next step is to obtain the frequency power spectrum of the above voltage diagram.  Figure 6.14 shows the frequency spectrum of the voltage across the nonlinear inductance.  200  300 Frquency (Hz)  400  F i g u r e 6.14: F r e q u e n c y p o w e r s p e c t r u m f o r the voltage A l - B l  119  600  F r o m the power spectrum diagram the presence of a dominant fundamental frequency, along with integer and non-integer higher-order harmonics, is confirmed.  T h e next step is to draw the phase-space plot of the flux vs. the voltage. A s we know, the phase plot of the chaotic motions has orbits which never close or repeat and tends to fill up a section of the phase space. Figure 6.15 shows the phase-plane diagram of the flux across the nonlinear inductance vs. its derivative, the voltage.  Again the phase-plot diagram  indicates a possible chaotic motion.  Next the Poincare map of the waveform is constructed by choosing a sampling frequency of 60 H z . This allows one to distinguish between the periodic motions and nonperiodic motions. F o r example, if the sampled motion were a subharmonic of period 3, the Poincare map would consist of a set of three points. But from Figure 6.16 it can be seen that the Poincare map of the waveform consists  of neither a finite set of points nor a  closed orbit. Here the Poincare map appears as an infinite set of highly organized points arranged i n what appear to be parallel lines. The appearance of fractal-like patterns in the Poincare map of a waveform is a strong indicator of chaotic motion.  120  4000 3000r 2000 h 1000  -1000 -2000 -3000 h -4000  Figure 6.15: Phase plot diagram of the flux vs. the voltage across the node A l - B l  Figure 6.16: Poincare map of the flux vs. the voltage across the node A l - B l  The next test is the calculation of the Lyapunov exponents. T h e Lyapunov exponent test measures the sensitivity of the system to changes i n initial conditions. Since a positive Lyapunov exponent implies the existence of chaotic dynamics, one needs only to calculate  121  the largest Lyapunov exponent, which tells whether nearby trajectories diverge (A. > 0) or converge (X < 0) on the average. F o r periodic orbits, all Lyapunov exponents are negative. The Lyapunov exponent calculated for the waveform of the Figure 6.12 is X = 0.064. The positive value of the Lyapunov exponent implies chaotic dynamics. T h e unit of X is in units of bits per data sample. Therefore, a value of +1 means that the separation of nearby orbits doubles on the average in the time interval between data samples.  The next technique for examining the pre-chaotic or post-chaotic changes in a dynamical system under parameter variations is the bifurcation diagram. B y drawing the bifurcation diagram one can see whether the system has periodic behavior for some range of the parameter space. It provides knowledge of the modes of the ferroresonance that are possible for any given value of the source voltage. T h e bifurcation diagram presented here is derived from the Poincare sampling of the voltage across the nonlinear inductance with the use of the E M T P . F o r a periodic waveform, there is one point per driving parameter and for chaos, there are many points. Figure 6.17 shows the bifurcation diagram for the waveform of Figure 6.12, when the varying parameter is the source voltage.  In Figure 6.17  the voltage across the nonlinear inductance (transformer voltage) is plotted  as the source voltage is increased from transformer  1.5 k V to 12 k V in steps of 100 V . Below 2200 V the  voltage moves linearly along the curve and is virtually sinusoidal. A t 2.1 k V  there is a discontinuity in the line of plotted points, which indicates the state of the system  122  has changed and the waveform is distorted due to the saturation of the transformer. This is one of the signatures of ferroresonant circuits.  T  1  1  1  1  • ll|  ii |  "  , u u u  0  1000  2000  3000  • 4000 5000 6000 7000 Source Voltage (V)  ' 8000  l. 1 ...J 9000 10000  1  Figure 6.17: Bifurcation diagram for varying the source voltage  At 2.2 kV a transition from period one to period two occurs. In the period-doubling route to chaos, the system starts with a fundamental periodic motion and as the source voltage is increased, the motion undergoes a bifurcation or change to a periodic motion with twice the period of the original oscillation. At 2.3 kV chaos occurs for this system; this is shown by a multiple of points per source voltage which indicates an aperiodic and distorted transformer voltage. As the voltage is increased, different windows of chaotic response, along with the periodic response, can be seen from the bifurcation diagram. This diagram can be used as a road map to find the safe mode of operation. Figures 6.18 to 6.20 show the Poincare maps for three different regions of period-one, period-two and chaotic  123  responses of the ferroresonance circuit. The Period-one region is for the source voltage of 7.0 kV, the period-two region is for the 9.5 kV, and the chaotic region is at lower source voltage of 5.5 kV. 2500  2000  _1500  1000  500  1  0.5  Flux[V.s]  2.5  1.5  V  F i g u r e 6.18: P o i n c a r e m a p f o r p e r i o d one response  200  -i  r-  2.05  2.1  source  = 7.0kV  180 160 140 _120 »  1100 5  BO 60 40 20j 02  2.15  2.2  _i  2.25 2.3 Flux [V.s]  F i g u r e 6.19: P o i n c a r e m a p f o r p e r i o d t w o response  124  1  i  2.35  V  2.4  source  i  2.45  = 9.5kV  2.5  2500  2000  1500  E §MOOO o  >  w  4  .1.5  -1  -0.5  0 0.5 Flux tV.s)  1  1.5  Figure 6.20: Poincare map for chaotic response period V  source  2  Z.5  = 55kV  In this chapter it was shown how the EMTP can be used for simulations of chaos in power systems. Simulations of ferroresonance, bifurcations and Poincare maps have been successfully produced with the use of the EMTP. The bifurcation diagram obtained from the EMTP proved valuable in showing many possible modes of ferroresonance that can exist over a wide range of source voltages. It was shown that, provided the ferroresonant circuit is appropriately modeled, the EMTP can be useful in analyzing the occurrence of chaos in power systems. The simulation results for chaotic systems should always be verified by using two or more different integration routines. For this reason, all the simulations presented in this chapter were run by three different integration routines: a fourth-order Runge-Kutta algorithm, a variable-order multistep Gear algorithm, and a trapezoidal algorithm. In all the simulations the results of using different integration routines were the same. This chapter also presented the result of a case study which was 125  performed by using the E M T P . This case study signified the existence of the chaos and bifurcations i n power systems and the need to develop accurate ways to reliably identify their presence.  126  Chapter  CONCLUSIONS A N D  In this chapter  7  RECOMMENDATIONS  the conclusions of this research are presented. T h e contributions of this  work in the areas of ferroresonance, nonlinear dynamics, chaos, and power systems are outlined. Recommendations and suggestions for further research in this area are made.  7.1 Conclusions The  trend toward using higher voltages for transmission of bulk power, along with  improvements in the characteristics of transformer steel, has increased the possibility of system disruption due to the occurrence of ferroresonance. Ferroresonance generally refers to  nonlinear oscillations in power systems  involving the excitation  of a nonlinear  inductance with a linear capacitance in series. These oscillations could cause damage to equipment such as transformers. Because of the nonlinearity in these systems, more than one mode of ferroresonance can occur. Ferroresonant voltages depend on the magnitude of the source voltage, length of the transmission line, losses, initial conditions and nonlinear inductance characteristics.  The aim of this thesis is to investigate ferroresonance from the standpoint of nonlinear dynamics and chaotic systems. T h e newly-developed techniques for analysis of nonlinear dynamical systems and chaos are evaluated for as a means for studying ferroresonance. In  127  general, we are to conduct exploratory research to determine whether  the potential for  chaos exists i n power networks and, if so, how to avoid or control such behavior.  In Chapter 2, the problem of ferroresonance in power systems is discussed and the method of graphical solutions which represents  one of the earliest attempts to explain this  phenomenon presented. A l l the methods of observation, categorization and analysis of various types of chaotic behavior are  presented, along with their use to  analyze,  understand, and control ferroresonance.  In Chapter 3, the equivalent circuit of a typical ferroresonant circuit in a power system presented along with the set of differential equations that describes it. T h e parameters of the ferroresonant model and their per unit values all calculated.  In Chapter 4, the method of Slowly Varying Amplitude used to derive an analytical solution for the equivalent ferroresonant circuit of Figure 4.1. A study of the frequencyresponse curve obtained from the solution shows that beyond a critical value of the amplitude of the flux linkage, hysteresis appears, leading to the occurrence of two stable attracting points and one region of unstable saddle points, which could lead to a chaotic response. T h e stability of the different portions of the frequency-response curve is determined by investigating the nature of the singular points. It shown that the eigenvalues of the characteristic equation after the flux linkage of 2 pu. when R=1%,  C for 100 k m , and  £•=0.15 pu., are real and opposite in sign, which corresponds to a saddle point. These saddle-type fixed points create the contraction and stretching mechanisms for horseshoe-  128  like maps. Multiple iterations, such as horseshoe maps, lead to complex orbits in phase space and start of chaotic behavior. Also in this chapter the Equivalent R i t z method used to obtain the graph of the variation of the source voltage vs. the amplitude of the flux linkage. This graph shows the multivaluedness of the response curve and the concept of the jump phenomenon. It is shown that beyond a critical value of the source voltage, a jump occurs leading to the occurrence of two stable attracting points and one unstable saddle point.  In Chapter 5, different modes of behavior for the ferroresonant circuit of Figure 3.5 were found by varying the values of the circuit parameters E, R, and C . T h e critical values of these parameters which drive the ferroresonant circuit into the chaotic region are found. Phase-plane  diagrams,  ferroresonant  behavior.  Poincare  maps,  and bifurcation diagrams  Variations in the  applied voltage,  losses,  used to or  categorize  length of  the  transmission line can result in different modes of behavior, which may be either periodic or chaotic. Some of the major findings of this chapter follow:  •  A s the applied voltage is gradually increased, the response goes into a chaotic region  through a series of period-doubling bifurcations. Lyapunov exponents and power spectra are calculated to confirm the existence of chaos beyond a critical value of the source voltage. A t the same time, as the value of the transformer core loss is decreased, a smaller source voltage is needed to drive the system into the chaotic region. This shows the importance of including the transformer core loss in ferroresonance studies. T h e inclusion of core losses, even small ones, ferroresonance  can  make the  difference between  chaotic  behavior  and ordinary  or periodic behavior of the solution. Therefore, the recent trend of  129  decreasing the losses of power transformers increases the possibility of chaotic behavior. T h e chaotic motion for lossless systems does not have any fractal structure; there are no attracting regions in phase space, and the points in the Poincare map hit all parts of the subspace of the phase space uniformly. F o r this reason, there is no basin of attraction for a conservative system. W i t h low losses in the systems, all numerical simulations must be allowed to run for a long time so that transient chaos is not misinterpreted as steady-state chaos and the intermittent chaos does not occur.  •  F r o m the variation of the parameter C , it is evident that the length of the transmission  line is an important factor in avoiding the chaotic region only when the source voltage is high or when the transformer core loss is low. Generally, with increased line length, the damping factor decreases and this rapidly increases the possibility of chaos. However, at the same time, by increasing the value of C , the effect of the nonlinear term is reduced. This is the reason why, i n the range of lower voltages, for greater lengths of the transmission line, the system response becomes periodic.  •  A more accurate representation of the § - i characteristic of the transformer results in a  prediction of chaos occurring at a much lower magnitude of the source voltage. The sharper the bending of the curve representing the transformer magnetization characteristic, the higher the likelihood that chaotic ferroresonance will appear.  •  Depending on the choice of initial conditions, the response of the ferroresonant circuit  may settle o n either a periodic or chaotic attractor. F o r a low-loss system, the basin of  130  attraction which results in a chaotic trajectory consists of all the points that the trajectory visits in the Poincare map plus the entire space outside the attractor. However, for normal losses and high values of source voltage, the basin of attraction for a chaotic trajectory is more complex, and has to be calculated o n a point-by-point basis. A s the forcing function is increased the basin of attraction for the chaotic attractor is increased, and at some point with a high value of forcing function, the basin of attraction for the chaotic attractor will include the entire phase plane.  •  T h e inclusion of series losses is not an important factor in studying ferroresonance as  compared to core losses. This is true until the value of series (copper) losses approaches a very high value of about 6% of rated transformer capacity, which is not generally the case in practical distribution systems.  •  T h e addition of a third harmonic with a magnitude of 10% of the fundamental to the  source voltage has a very limited effect on appearance of chaos in the solution.  Chapter 6 described the implementation of the E M T P to study ferroresonant circuits. It is shown that the E M T P can be used to model the individual components of power systems containing any nonlinear elements. T h e method of using E M T P to obtain the Poincare map and the bifurcation diagram for study of chaos outlined. In addition the influence of nonlinear modeling in the E M T P o n the results of the simulation is discussed. Some of the major findings of this chapter follow:  131  •  In order to avoid the use of the C D A scheme and to have a purely trapezoidal solution,  the compensation method is chosen to model the transformer nonlinearity.  •  T h e appropriate representation of nonlinear elements in the E M T P  when chaotic  systems are to be simulated is very important. T h e main point to consider is the number of the points in the bending region of the saturation curve.  •  T o obtain the  Poincare map in the E M T P , the integration time step and output  frequency must be adjusted so that their product equaled a 60 H z period.  •  T h e results obtained from the E M T P simulations for the ferroresonant circuit of Figure  3.5 completely match with those obtained by solving the differential equation of the system by I N S I T E .  •  The simulation results of chaotic systems should always be verified by using two or  more different integration routines. F o r this reason, all the simulation results that are presented in this study are obtained by three different integration routines: a fourth-order Runge-Kutta algorithm, a variable-order multistep Gear algorithm, and a trapezoidal algorithm. In all the simulations the results of using different integration routines are the same.  •  The results of the case study performed by using the E M T P shows the existence of chaos  and bifurcations in a sample power system and the need to develop accurate ways to reliably identify chaotic behavior.  132  7.2 Recommendations Chaos theory and the techniques of nonlinear dynamics can provide new approaches for understanding transient and long-term behavior of complex power systems,  and for  ensuring network stability. In recent exploratory studies at the University of California at Berkeley, and at Cornell University, researchers concluded that events such as voltage collapse, low-frequency electromechanical oscillations and transient stability may be linked to chaotic behavior [39].  These studies along with the results of this thesis all recommend  that there should be more research done in the development of more complete power system simulation tools in order to facilitate the efforts to determine if the various chaotic behaviors observed in simplified nonlinear models occur near the ordinary operating conditions of real networks or only under extreme conditions that are of little practical concern. Some of the suggestions for continued research in this area follow:  •  M o r e research  needs to  be done in determining whether  the  large  number  of  interlocking subsystems that compose an actual network allow chaotic fluctuations to develop.  •  Modeling results are dependent on proper representation of the core loss. A better way  of modeling transformer-core loss is needed. Characteristics of new low-loss core materials should be investigated and modeled. W i t h recent trends of reducing system losses, the potential for problems with transients is increased.  133  •  T h e results presented in Chapters 5 and 6 should be used to investigate the real-time  detection and suppression of ferroresonance.  134  REFERENCES  [I] .  H o p k i n s o n , R . H , "Ferroresonance During Single-Phase Switching of a 3-Phase Distribution Transformer Banks," I E E E Trans, o n Power A p p . and Systems, A p r i l 1965, pp. 289-293.  [2].  Smith, D . R . , Swanson, S.R. and Borst, J . D . , "Overvoltages with RemotelySwitched Cable-Fed Grounded Wye-Wye Transformers,"  I E E E Trans, on  Power Apparatus and Systems, V o l . 94, N o . 5, Sept./Oct. 1975 pp.1843-1853. [3]. [4].  Nayfeh, A . and M o o k , D . , Nonlinear  Oscillation, John  Wiley & Sons, 1979.  Schultz, R . A . , "Ferroresonance in Distribution Transformer Banks on 8/34.5 k V Systems," Presented at 1964 Spring Conference, R o c k y Mountain Elec. League.  [5].  Rudenberg, R.,  Transient Performance of Electrical Power Systems.  McGraw-Hill  Book Company, N e w Y o r k , N Y , Copyright © 1949. [6].  Butler, J . W . and Concordia, C , "Analysis of Series Capacitor Application Problems," A I E E Trans. V o l . 56, August 1937, pp. 975-988.  [7],  Dolane, E . , Gillies, D . and Kimbark, E . "Ferroresonance in a Transformer Switched W i t h an E H V Line," I E E E Summer Meeting, Portland, July 18-23, 1971, pp. 289-293.  [8].  Rudenberg, R., "Nonharmonic Oscillations as Caused b y Magnetic Saturation," Trans. A I E E , V o l . 68 I, 1949, pp.676-685.  [9].  Swift,  G . W . , "Power  Transformer  Core  Behavior  Under  Transient  Conditions," I E E E Trans, on Power Apparatus and Systems, V o l . 90, N o . 5, Sept./Oct. 1971, pp. 2206-2209 [10].  Jiles, D . C . and Atherton, D . L . , "Ferromagnetic Hysteresis," I E E E Trans. Magnetics, V o l . 19, N o . 5, Sept. 1983, pp. 2183-2185.  [II] .  Feldman, J . M . and Cappabianca, A . L . , " O n the Accuracy and Utility of Piecewise-Linear Models of Ferroresonance," I E E E Trans, on Power Apparatus and Systems, V o l . 97, N o . 2, M a r c h / A p r i l 1978, pp. 469-477.  [12].  Frame, J . G . , M o h a n , N . , and L i u , T . , "Hysteresis Modeling i n an ElectroMagnetic Transient Program," I E E E Trans, on Power Apparatus and Systems, V o l . 101, N o . 9, Sept. 1982, pp. 3403-3411. 135  [13].  Hopkinson, R . H . , "Ferroresonant Overvoltages C o n t r o l Based o n T N A Tests on  Three-Phase  Wye-Delta Transformer  Banks," I E E E Trans,  on Power  Apparatus and Systems, V o l . 87, N o . 2, Feb. 1968, pp. 352-361. [14].  Kieny,  C.  "Application of the  Bifurcation  Theory  i n Studying and  Understanding the Global Behavior of a Ferroresonant...," I E E E Trans. Power Delivery, V o l . 6, N o . 2, pp. 866-872, A p r i l 1991. [15].  Gleick, J . D . ,  Chaos: Making a New Science.  Viking, N e w Y o r k , N Y , Copyright  © 1987. [16].  Dynamics and Chaos Geometrical Methods for Engineers and Scientists, John Wiley & Sons, N e w  Thompson, J . M . T .  and Stewart, H . B . , Nonlinear  Y o r k , N Y , Copyright © 1986. [17].  Mandelbrot, B.B.,  The Fractal Geometry of Nature,  W . H . Freeman and  Company, San Francisco, C A , Copyright © 1983. [18].  Moon, F.C.,  Chaotic Vibration - An Introduction for Scientists and Engineers,  John Wiley & Sons, N e w Y o r k , N Y , Copyright © 1987. [19].  Nayfeh, A . H . and Balachandran, B.,  Applied Nonlinear Dynamics , John  Wiley  & Sons, N e w Y o r k , N Y , Copyright © 1995. [20].  Strogatz, S . H . , Nonlinear Dynamics  and Chaos,  Addison-Wesley, Reading, M A ,  Copyright © 1994. [21].  Yorke, J . A . , "Period Three Implies Chaos", American Mathematical Monthly, V o l . 82, pp. 985-992, 1975.  [22].  Marti,  J.R.,  and Soudack,  A . C . , "Ferroresonance  i n Power  Systems:  Fundamental Solutions", I E E Proc. C , pp. 321-329, July 1991. [23].  Araujo, A . E . A . , Soudack, A . C . and Marti, J.R., "Ferroresonance i n Power Systems: Chaotic behavior", I E E Proc.-Gener. Trans. Distrib., V o l . 140, N o . 3, pp. 237-240, M a y 1993.  [24].  Araujo, A . E . A . , Mozaffari, S., Soudack, A . C . and Marti, J.R., "Chaos in Power Systems:  EMTP  Simulations",  Proc.  11  th  Power  Systems  Computation  Conference, pp. 671-677, Sept. 1993. [25].  Mozaffari, S., Henschel, S. and Soudack, A . C , "Chaotic Ferroresonance i n Power Transformers", I E E Proc.-Gener. Trans. Distrib., V o l . 142, N o . 3, pp. 247-250 M a y 1995.  136  [26].  Mozaffari, S., Sameti, M . and Soudack, A . C , "The Effect of Initial Conditions on Chaotic Ferroresonance in Power Transformers", Under Review.  [27].  Dick, E.P., and Watson, W., " Transformer Models for Transient Studies Based on Field Measurements", IEEE Trans., pp. 409-417, 1981.  [28].  Kantorovich, L.V. and Krylov, V.I., Approximate Methods of Higher Analysis, Interscience Publishers, Inc., Copyright © 1958.  [29].  Migulin, A., Basic Theory of Oscillations, Mir Publishers, Moscow, 1983.  [30].  Chakravarthy, S.K. and Nayar, C.V., "Series Ferroresonance in Power Systems", Journal of Electrical Power & Energy Systems, Vol. 17, No. 4 , pp. 267-274, 1995.  [31].  Soudack, A . C , Electrical 557 Lecture Notes, Department of Electrical Engineering, The University of British Columbia, Vancouver, Canada  [32].  Parker, T.S. and Chua, L.O., Practical Numerical Algorithms for Chaotic Systems, Springer-Verlag, Copyright © 1989.  [33].  Hilborn, R . C , Chaos and Nonlinear Dynamics, Oxford University Press, Copyright © 1994.  [34].  Strogatz, S.H., Nonlinear Ehnamics and Chaos, Addison-Wesley Publishing Company, Copyright © 1994.  [35].  McPherson, G., An Introduction to Electrical Machines and Transformers, John Wiley & Sons, New York, N Y , Copyright © 1981.  [36].  Dommel, H.W., Electromagnetic Transients Program Reference Manual (EMTP Theory Book), Department of Electrical Engineering, The University of British Columbia, Vancouver, Canada, 1986.  [37].  Mozaffari, S., Dommel, H.W. and Marti, J.R., "EMTP Studies Without Numerical Oscillations in Networks With Convertor Circuits", Canadian Electrical Association, Engineering and Operating Division Transactions, April 2, 1992.  [38].  Marti, J.R. and Lin, J., "Suppression of Numerical Oscillations in the EMTP", IEEE Trans, on Power System, Vol. 5, No. 2, May 1990, pp. 394402.  [39].  Douglas, J., "Seeking Order in Chaos", Electric Power Research Institute Journal, June 1992.  137  APPENDIX A  This A p p e n d i x shows the E M T P input listings  for the simulation of the ferroresonant  circuit of Figure 3.5, and the test case study of Figure 6.11.  The first E M T P listing is for the ferroresonance circuit of Figure 3.5.  * Case identification card C H A O T I C BEHAVIOUR I N POWER SYSTEM * Ferroresonance circuit of Figure 3.5 * W i t h 21 points representation of the saturation curve * E = 1.89 p u - R = 1.0% -L = 777 nF(lOOkm) - Laux = 1.0*10**12 H * T i m e card 1.85185E-566.6666667900 *  *  1  L u m p e d R L C branch NODE1NODE2 .777 NODE2 48.4e3 NODE2 1.E15  * 93  1  Nonlinear inductance (card 1) NODE2  * -12973.6 -400.55 -7270.93 -380.02 -3945.99 -359.47 -2066.00-338.93 -1038.93-318.39 -499.19-297.85 -227.76-277.30 -97.99 -256.77 -39.44 -236.22 -14.78 -215.69 -5.1876-195.14 -1.78-174.60 -0.69 -154.06 -0.37 -133.52  Nonlinear inductance (card 2)  138  -0.26-112.97 -0.20  -92.44  -0.16  -71.89  -0.12  -51.36  -0.07  -30.81  -0.03  -10.27  0.000  0.000  0.03  10.27  0.07  30.81  0.12  51.36  0.16  71.89  0.20  92.44  0.26  112.97  0.37  133.52  0.69  154.06  1.78  174.6  5.1876  195.14  14.78  215.69  39.44  236.22  97.99  256.77  227.76  277.30  499.19  297.85  1038.93  318.39  2066.00  338.93  3945.995  359.47  7270.93  380.02  12973.62  400.55  9999999. $  =  $ =  = =  E n d o f l e v e l 1: L i n e a r a n d n o n l i n e a r e l e m e n t s == = =  =  =  E n d o f l e v e l 2: S w i t c h e s a n d P i e c e w i s e l i n e a r e l e m e n t s  *  Voltage or current sources  14NODE1 $ =  =  =  120031.121  60.0  E n d of level 3: Sources  *  -90.0 =  =  =  =  =  I n i t i a l c o n d i t i o n s (voltages) 2 N O D E 2  *  89814.0 I n i t i a l c o n d i t i o n s (currents)  3 N O D E 1 N O D E 2 3  -89814.0  N O D E 2  * *  Voltage-output  nodes 139  =  =  =  =  =  = =  NODE2 $ = = = E n d of level 4: User-defined voltage output $ = = =  Level 5: E n d of data case = = = = = = =  The next E M T P listing is for the American Electric Power case, i n which the occurrence of questionable voltage spikes was reported. T h e circuit is shown in Figure 6.11.  *  Case identification card  Case study from the A E P as shown in Figure 6.11 *  -1  T i m e card  .20008e-4 1000.0e-3 1 . *  1  Lumped R L C branch  * Equivalent source Z's at 12.47 K v L - L R . M . S . Generator nodelnodelA  8.70019.761  * Measuring inductance across the primary for flux nodeAlnodeBl  l.e23  *  1  * 2.1 K v winding nodeXlnodeYY *  12.631  . 480V./2400V.  .  . Single-phase transformer (simplified)  51nodeAlnodeBlINVERSE  8.0000000000E-02 1.8848720172E+00  52nodeXlnodeYl  0.0000000000E+00-3.7697440344E-01 2.0500000000E+00 7.5394880688E-02  *  Nonlinear inductance (card 1)  93nodeAlnodeBl *  1.24921.9806 Nonlinear inductance (card 2)  -905000. -42.58 -5000.0  -2.58  -500.0  -2.38  -80.0  -2.29  40.0 -2.284 4.39880-2.28124 -3.83015-2.28057 -1.92095-2.18060 140  -1.43780-2.08060 -1.24920-1.98060 -1.13430-1.88060 0.0  0.0  1.13430 1.88060 1.24920 1.98060 1.43780 2.08060 1.92095 2.18060 3.83015 2.28057 4.39880 2.28124 40.0  2.284  80.0  2.290  500.0  2.380  5000.0 2.580 905000.0 42.58 9999999. *  Lumped R L C branch nodeAlnodeBl  100.0  * Load impedance nLBl  69.72445.62  $ = = E n d of level 1: Linear and nonlinear elements *  Time-controlled switch nodelAnodeAl  -1.0  1196.00  5.00  nodeYlnodeYY  -1.0  1196.00  5.00  nodeBlnodeBB  -1.0  1196.00  5.00  1.0e-5  nodeLBnodeBB  -1.0  1196.00  5.00  1.0e-5  *  Time-controlled switch nodeBB  -1.0 996.000  5.00  $ = = = E n d of level 2: Switches and Piecewise linear elements *  Voltage or current sources  14nodel  10180.0  60.0  0.0  $ = = = E n d of level 3: Sources = = = = = = = = = = = =  *  Voltage-output nodes nodel  $ = = = E n d of level 4: User-defined voltage output $ = = = Level 5: E n d of data case = = = = = = = 141  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

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

Comment

Related Items