S I M U L A T I O N O F LIGHTNING SURGES O N TRANSMISSION LINES By H U Y E N V. N G U Y E N B.Sc, Elec. Engr., Gannon University, 1987 M.Engr., Electric Power Engr., Rensselaer Polytechnic Institute, 1988 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF T H E REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE 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 to the required standard THE UNIVERSITY OF BRITISH COLUMBIA February 1996 © Huyen V. Nguyen, 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 of ZieCTXlCAL €M£IlV€£fZlh)fr-The University of British Columbia Vancouver, Canada Date h^L re, w c , DE-6 (2/88) ii ABSTRACT Computer simulations with programs such as the EMTP are widely accepted for light-ning transient studies, to design the system in such a way that equipment failures and trans-mission line outages are prevented. In these studies, an accurate representation of the transmission lines and towers is required. Two new and improved time-domain models for the representation of these components are presented in this thesis. They are the frequen-cy-dependent and the exponential transmission line models. The former is based on syn-thesizing the line functions directly in phase co-ordinates, and is used to represent multiphase transmission lines, while the latter is a single-phase model, which also syn-thesizes the line functions, to be used for tower representations. These models are valid for the required frequency range of 1 Hz to 1 MHz. They can be used for lightning studies as well as for many other types of transient studies. Extensive time-domain simulations have been performed to compare the results with those produced by existing line models. The re-sults indicate that the new models are more accurate than some of the existing models, and numerically stable. Simulations of typical lightning surges have also been performed to compare the sol-utions produced by the developed models with actual measurements. With these new mo-dels, simulations match field recorded transient waveforms very well. In addition, a simple, fast and accurate method based on the Newton-Raphson iteration technique is introduced for finding transmission line modal parameters. It overcomes the ei-genvector and eigenvalue switching problems experienced by most other methods. With a proper constraint equation for scaling the eigenvectors, this proposed method produces modal parameters which can reasonably be fitted with rational functions of minimum-phase type. With this, the modal transmission line modelling can be improved. The routine is very stable and can be used to diagonalize the Y Z product for uncommon lines with strong asymmetry. It has also been used to obtain the phase-domain parameters for the new fre-quency-dependent line model. iii Table of contents Abstract ii Table of contents iii List of Tables — . : •. — v List of Figures v Acknowledgement ix CHAPTER 1: Introduction 1 1.1) Background 1 1.2) Types of lightning transients on overhead lines 2 1.3) Lightning impacts on high voltage power systems 2 1.4) Thesis outline 4 1.5) Thesis contributions 5 CHAPTER 2: Literature Review 12 2.1) Lightning waveform measurements in power systems 12 2.2) Lightning stroke representation in engineering calculations 13 2.3) Tower surge impedance and its response to lightning strokes 15 2.4) Transmission line characteristics for lightning studies 17 CHAPTER 3: Modelling of Single-Phase Nonuniform Transmission Lines for Representation of Transmission Towers 20 3.1) Overview : 20 3.2) Exponential transmission line 21 3.3) Synthesis of exponential line functions 25 3.4) Reflection factors — 31 3.5) EMTP simulations 35 CHAPTER 4: Frequency-Dependent Transformation Matrices for Untransposed Transmission Lines using Newton-Raphson Method 37 4.1) Overview 37 4.2) Equations for multiphase transmission lines — 38 4.3) Eigenproblem solution with Newton-Raphson (NR) method .. .42 4.4) Modal parameters for multiphase transmission lines from Newton-Raphson method 45 4.5) Open/Short -circuit tests in the frequency domain : — 54 4.6) Validation of the Newton-Raphson algorithm with the power method 60 CHAPTER 5: Frequency-Dependent Overhead Transmission Line Model in the Phase Domain 62 5.1) Overview 62 5.2) Phase-domain transmission line modelling 62 5.3) Synthesis of the propagation and characteristic matrices 65 5.4) Synthesis of Y c and A for asymmetrical transmission line configurations 73 5.5) Open/Short -circuit tests in the frequency domain 77 5.6) Time-domain simulation tests 87 CHAPTER 6: Simulation of Lightning Surges and Field Test Comparisons 91 6.1) Overview 91 6.2) Tower/ground wire surges 91 6.3) Simulation of a recorded shielding failure event in a 765 kV system 96 6.4) Shielding-failure flashover 103 6.5) Back flashover 104 CHAPTER 7: Future Research in Transmission Line and Transmission Tower Modelling I l l CHAPTER 8: Conclusions 113 APPENDIX A: Asymptotic Behaviour of Modal Transformation Matrices at Very Low and High Frequencies 115 APPENDIX B: Recursive Convolutions ..118 APPENDIX C: Interfacing User-defined Models with the EMTP version MicroTran . .121 APPENDIX D: Estimation of Maximum Lightning Shielding Failure Stroke on a Typical 765 kV Line Based on Electrogeometric Theory 124 APPENDIX E: Representation of a Long Line with Linear A C Steady-State Initial Conditions 128 REFERENCES 130 V List of Tables Table 5.1: Average number of poles allocated for synthesized functions 74 Table 5.2: Steady-state induced voltages (RMS) 90 List of Figures Figure 1.1- Multiple shielding failures waveforms recorded on March 26, 1991 6 Figure 1.2 - Expanded view of stroke #1 of Figure 1.1 7 Figure 1.3 - Expanded view of stroke #2 of Figure 1.1 8 Figure 1.4 -Expanded view of stroke #3 of Figure 1.1 9 Figure 1.5 - Shielding failure waveforms recorded on April 12, 1991 10 Figure 1.6- Expanded view of strike location of Figure 1.5 11 Figure 3.1 - Single-phase exponential line 21 Figure 3.2 - Equivalent circuit for the exponential line 25 Figure 3.3 - Magnitude of Z c k in Q 27 Figure 3.4 - Phase of Z c k in degrees — 27 Figure 3.5 - Magnitude of A 28 Figure 3.6 - Phase of A in radians 28 Figure 3.7 - Propagation function a(t) 29 Figure 3.8 - Magnitude of Z ' c m in Q 30 Figure 3.9 - Phase of Z ' c m in degrees 30 Figure 3.10- Single-phase exponential line terminated at node m with Z L 31 Figure 3.11 - Magnitude of Z i n in Q from Eq. (3.24) 33 Figure 3.12- Magnitude of Reflection factor for Z L = 150 Q 34 Figure 3.13- Magnitude of Reflection factor for Z L = 100 Q 34 Figure 3.14 - Receiving end open-circuit voltage (p.u.) 36 Figure 3.15- Receiving end short-circuit current (p.u.) 36 Figure 4.1 - Multiphase transmission line system 39 Figure 4.2 - Vertical double-circuit line 47 VI Figure 4.3 - Elements of column 1 of the matrix T, 50 Figure 4.4 - Elements of the matrix Y C M 51 Figure 4.5 - Elements of the matrix A M 52 Figure 4.6 - Phase angle of element A M . , calculated from Eq. (4.21) 53 Figure 4.7 - Phase angle of element AM_, calculated with atan function 53 Figure 4.8 - Receiving end open-circuit voltages (p.u. magnitude) 56 Figure 4.9 - Receiving end short-circuit currents (p.u. magnitude) 57 Figure 4.10 - Receiving end open-circuit voltages (p.u. magnitude) 58 Figure 4.11 - Receiving end short-circuit currents (p.u. magnitude) 59 Figure 4.12- Locations of largest elements of the correlation matrix 61 Figure 5.1 - Multiphase transmission line equivalent circuit in the time domain 65 Figure 5.2 - Elements of column 1 of matrix Y c for a double-circuit line 69 Figure 5.3 - Elements of column 1 of matrix P 70 Figure 5.4a - Phase angle (degrees) of elements P 2 1 and P 3 I before multiplying with -1.0 . 71 Figure 5.4b - Phase angle (degrees) of elements A n and A 2 1 from Eq. (5.12) 71 Figure 5.5 - Time-domain propagation functions: au(t), -a21(t), and -a31(t) 72 Figure 5.6 - Current source with matching admittance Yc(co) connected to node m 73 Figure 5.7 - Railroad near a three-phase 138 kV line 74 Figure 5.8 - Elements of column 1 of matrix Y c for an asymmetrical line 75 Figure 5.9 - Elements of column 1 of matrix A for an asymmetrical line 76 Figure 5.10 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 4.2 79 Figure 5.11 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 4.2 80 Figure 5.12 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 5.7 81 Figure 5.13 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 5.7 82 Figure 5.14- Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 4.2 83 Figure 5.15 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 4.2 84 Figure 5.16 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 5.7 85 Figure 5.17 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 5.7 86 Figure 5.18 - Receiving end phase 1 open-circuit voltage (p.u.) 87 V l l Figure 5.19 - Magnetically induced voltages on rails 89 Figure 5.20 - Magnetically induced voltage on rail #1 89 Figure 6.1 - Circuit for EMTP simulations 93 Figure 6.2 - Conductor geometry at average height 94 Figure 6.3 - Tower representation -94 Figure 6.4 - Simulation of Ishii et al experiment and measurement from [18] 95 Figure 6.5 - Simplified diagram of 765 kV system surrounding Rockport 99 Figure 6.6 - Lightning current source (A) 99 Figure 6.7 - Lightning shielding failure waveforms recorded on Oct. 7, 1990 100 Figure 6.8 - Expanded view of recorded shielding failure waveforms 101 Figure 6.9 - EMTP simulations of the shielding failure event in Figures 6.7 & 6.8 102 Figure 6.10 - EMTP simulations of a shielding-failure flashover event 107 Figure 6.11 - EMTP circuit for simulation of a back flashover event 108 Figure 6.12 - Lightning current source (A) 108 Figure 6.13 - Simulated voltages at tower #7 with a 15 Q footing resistance 109 Figure 6.14 - Simulated voltages at tower #7 with a 180 Q footing resistance 110 Figure C l - Equivalent circuit seen by the transmission line terminals 121 Figure D l - Shielding angle of the right conductor bundle 124 Figure D2 - Incomplete Shielding, Width X s is uncovered 127 Figure D3 - Effective shielding, Unprotected width X s is reduced to zero 127 Figure E l - Single transmission line 128 Figure E2 - Equivalent circuit of a long line without initial conditions 129 Figure E3 - Equivalent circuit of a long line with linear A C initial conditions 129 viii In memories of Reverend Monsignor Wilfrid J. Nash and my two dear sisters Dung and Xuan ix Acknowledgement I would like to express my deepest appreciation to all the people who contributed directly or indirectly to this thesis project. The following thanks are expressed specifically to those who have and had enabled me in realizing my full potential: To my thesis supervisors, Drs. Hermann Dommel and Jose Marti, whom I respect and ad-mire tremendously. Not enough can be said about my appreciation for all they have helped in giving guidance, encouragement, as well as arranging much-needed financial support. Their unique, patient, and reasonable ways of getting ideas across had allowed me to learn and per-form so freely. To Dr. Martin Wedepohl for his patience and effort in introducing the Newton-Raphson method of Chapter 4. To Dr. Luis Marti for his help and support as Industrial Technical Advisor. To Drs. A.J.F Keri, J.M. Schneider and C.H. Shih for providing field recordings and sys-tem parameters for the shielding failure and railroad cases, and also for training and guiding me in conducting assigned projects during my stay at the American Electric Power Service Corporation, Columbus, Ohio. It was during this period that I have become so interested in what is called "power system transients" To Dr. Jerry Selvaggi for passing his love of engineering on to me. To all my colleagues, especially to Fernando Castellanos and Luis Linares, in the power group at UBC for creating such a warm and friendly research environment. The financial assistance of the Natural Sciences and Engineering Research Council of Ca-nada through a grant to Dr. Hermann Dommel, and the graduate fellowships awarded by the University of British Columbia are gratefully acknowledged. I am also deeply grateful to my dearest uncle Van for his encouragement and help in liber-ating me to America. Thanks are addressed to all my brothers Doan, Tuan , Huy, Van, and to my oldest sister Hang, for their help and support. Thanks are also given to all my friends at the Vietnamese Catholic Church choir (1992-1996) in Vancouver, B.C., for their kindness in many ways. It is through this organization that I have been blessed to have known my wonderful wife Hien. Thanks are deeply expressed for her love, patience and support during these years. And finally, I would like to dedicate all my achievements to my dear parents. Thanks for their love and support through the years. 1 C H A P T E R 1 Introduction •i , . . • . 1.1) Background Analyses of field recording events over many years [2,8,9] have suggested that high overvoltage transients caused by lightning are considered to be the chief source of disturb-ances in high voltage transmission systems. Lightning phenomena have been studied throughout this century, although the process that initiates lightning is still not completely understood. This has led to the development of various theories in attempting to describe the so-called lightning-mechanism. There is, however, a consensus that lightning starts from the charge separation process (positive and negative), which is due to the transportation of lightweight particles to higher regions by the rapid updrafts of moist air, usually in hot humid areas. This vertical charge separation is known as the vertical thunderstorm dipole. It can be formed within the cloud or between the cloud and the earth, which creates electric fields that eventually bring about the electric breakdown known as lightning [3,4]. It is not within the scope of this thesis to fully explain the details of how the lightning phenomenon develops. What concerns the power system engineer is the severity of overvoltages caused by lightning in the transmission grid. The overvoltages introduced by lightning have tradi-tionally been estimated using conventional and simplified methods [3,5,7]. More involved calculations have become possible with digital computer programs such as the 2 Electromagnetic Transients Program (EMTP) [1]. In such programs, each power system component can be modelled in great detail. 1.2) Types of lightning transients on overhead lines The characteristics of lightning surges on overhead transmission lines, which result from lightning strokes, depend on how they are caused. They can be broadly divided into four types (ignoring the induction case): • a) Tower/ground wire surges - The stroke terminates on the tower structure/ground wires without any back flashover to the phase conductors. • b) Shielding failure - The stroke passes through the protective zone of the ground wires and terminates on the phase conductor(s). • c) Back flashover - The same as (a), but followed by a flashover to the phase conduc-tors). This type of flashover is called back flashover. • d) Shielding-failure flashover - The same as (b), but followed by a forward flashover to the ground/ground wires or tower. 1.3) Lightning impacts on high voltage power systems Lightning overvoltages in power systems are either caused by a lightning stroke to the line conductors, the tower structure and the ground wires, or they can also be created by an induction process from a nearby stroke to ground. The latter causes insignificant overvoltages in high voltage systems, but the former type may result in system faults from either forward or back flashovers. These flashovers can damage terminal equipment. The 3 effects produced by those events which do not cause flashovers can also be severe. There have been many interesting shielding-failure events captured in the AEP (American Electric Power) 765 kV system, which created moderate overvoltages [2]. Figure 1.1 depicts the phase-to-ground waveforms at a 765 kV transformer terminal, which resulted from a mul-tiple shielding-failure event, captured during heavy thunderstorms. This event is very unique since three consecutive different strokes, within a 110-ms period, terminated directly on the phase A conductor approximately 149 miles away; flashover occurred on the last stroke. The expanded views at the stroke locations are exhibited in Figures 1.2-1.4. The time between the first and the second, and the second and third strokes are approximately 65 and 45 ms, respectively. Figure 1.5 also shows another extraordinary shielding-failure event, which occurred near the negative peak of the system 60 Hz voltage. The expanded view is shown in Figure 1.6. Phase C is the struck phase in this case; it registered a maxi-mum magnitude of 1.80 p.u.! These overvoltages, with non-standard waveshapes, are fully seen by the system equipment at the substations. The resulting effects of these non-standard waveforms impinging across the substation transformers, for example, are not yet well known to the transformer manufacturers. Many unexplained E H V transformer failures [2] may partially have been caused by these effects. Thus, it is important to realistically simu-late and fully understand all the possible events created by lightning in power systems. The ensuing outcome can then be helpful for the future development of improved protection schemes. 4 1.4 Thesis outline Power system transients, especially those resulting from lightning, are by and large computed with general-purpose programs such as the EMTP [44]. The goal of this thesis is to develop new and improved transmission line and transmission tower models for lightning surge simulations for use in such computer programs. Overhead lines and transmission towers are both entirely exposed to lightning activities, and therefore their behaviour in-fluences greatly the travelling waveforms reaching the substations. For modelling other sys-tem components, guidelines are available which produce results with reasonable accuracy [6], [7], and [44]. The thesis gives a brief literature review in Chapter 2, where the lightning return-stroke current measurements and models for overhead transmission lines and towers are discussed. A nonuniform transmission line model, based on an exponential variation of the line para-meters, is presented in Chapter 3. This model is intended for the representation of trans-mission towers. For developing a time-domain frequency-dependent transmission line model, an eigenvector/eigenvalue routine based on the Newton-Raphson technique is in-vestigated and validated in Chapter 4. A phase-domain frequency-dependent transmission line model is described in Chapter 5, along with various frequency and time domain test re-sults. In Chapter 6, simulation results of common lightning events, along with field test comparisons, are presented. Possible for future research are outlined in Chapter 7, and con-clusions about the research findings of the thesis are made in Chapter 8. 5 1.4 Thesis contributions The original contributions from the work in this thesis are: • Implementing and validating the Newton-Raphson method for evaluating the modal parameters of overhead transmission lines as smooth functions of frequency [37]. • A time-domain algorithm for modelling single-phase nonuniform transmission lines for the representation of transmission towers [22, 46]. • A time-domain algorithm for direct phase-domain modelling of frequency-dependent overhead transmission lines [47]. • Modelling techniques for simulation of typical lightning events. 6 7 at o > o Q. O O a. o at o > m a> in in t o Q. o o a: O) 2 o > t o a O o CC - 6 "3~n3 3 0 . 5 Figure 1 . 2 - Expanded view of stroke # 1 of Figure 1 . 1 ir^r 3 1 . 5 Time (ms) x 1 0 5 Time (ms) - 2 . 52. - 3 . > - 3 . 5 . Time (ms) Figure 1.4 -Expanded view of stroke #3 of Figure 1.1 10 12 CHAPTER 2 Literature Review 2.1) Lightning waveform measurements in power systems The overvoltages caused by lightning in power system are influenced by various char-acteristics of the lightning flashes. These characteristics include the stroke leader mechan-ism, the length-of-path, current polarity, time duration, charge distribution in the stroke channel and the total number of strokes. Many researchers have attempted to measure light-ning currents from strokes to high towers. For example, the C N tower in Toronto is equipped with instrumentation to measure lightning stroke currents. Others have used lightning-flash counters to obtain statistics of the frequency of lightning strokes to earth, while others have used the power system itself as a means to record lightning transients, which impinge upon terminal equipment. Recordings from a monitoring system, which has been installed in the 765 kV system of American Electric Power Service Corporation, will be used in this thesis to compare simulation results against field measurements. Each type of transient phenomenon occurring in a power system requires a specific measuring technique, which depends primarily on the frequencies involved in the phenom-enon. Lightning phenomena are of the fast transient type (10 kHz - 3 MHz [6]). Its wave-forms in power systems have been captured with various monitoring systems [2,8,9,10]. Such monitoring systems must have a wide frequency bandwidth and fast sampling rates. 13 The system described in [2] captures lightning events in the form of voltages, through volt-age dividers connected to transformer bushing capacitance taps. This monitoring system is of particular interest to the author since he managed and maintained it for a three-year period during his stay at the American Electric Power Service Corporation. During this time, it captured many interesting lightning events, such as those shown in Figures 1.1-1.6. 2.2) Lightning stroke representation in engineering calculations One of the most critical factors in lightning surge simulations is the representation of the • lightning flash itself. It has been accepted that when lightning strikes a power line, the light-ning stroke injects a current source into the power system [3,7]. The injected current repre-senting the lightning return-stroke is generally modelled as a current source, with the stroke-channel impedance often being ignored because it is relatively large. A single light-ning flash may contain several strokes, with the average being three. If each stroke has a certain waveshape and magnitude, it may be necessary to represent the flash as a series of current sources in simulations. To understand the significant effects created by a multiple shielding-failure event, this kind of multi-stroke representation with lightning current sources is needed. In most cases, however, a single current source for the representation of only the leader stroke is justified [7]. Some other factors may be important as well. For example, Greenwood [3] mentions that a more accurate representation of lightning must take the leader stroke and the prestrike into account. Some researchers have theorizing that the effect of charge distribution in the stroke-channel should be included in lightning calcula-tions too. Researchers have held different viewpoints about this particular parameter. 14 Wagner-Hileman [11] stated that the effect of the charge distribution in the stroke-channel on the line insulator voltage is as important as the effect of the injected current. Choy and Darveniza [12] re-examined this theory and confirmed the significant influence of this para-meter on the predicted flashover rates. In their analysis, they assumed that the length of an upward streamer is one half of the striking distance. They found that the resulting insulator voltage due to charge in the stroke-channel decreases significantly as the length of the streamer increases. Subsequently, Darveniza et al [12] reversed this opinion and recom-mended to ignore this 'presumably' important parameter in backflash calculations. Most re-searchers today seem to model only the effects of the injected current stroke. This approach is also adopted in this thesis. Thus, the effect of charge in the stroke-channel will not be considered here. As far as the representation of the injected lightning current is concerned, factors such as the waveshape and the rate-of-rise significantly influence the overvoltages which lightning can create. Most earlier studies assumed the return-stroke current as a triangular-ramp func-tion, because of the need for simplicity in the analyses involved. IEEE recommendations for surge arrester studies still use such simple waveforms. With the advent of digital simula-tions, which make more complicated source representations easy, more recent studies often used a double exponential current source. The most recent lightning waveform, proposed by O G R E Working Group 33 [7], gives a somewhat more realistic picture, since it has the con-cave wavefront, similar to that observed in measurements. It is therefore recommended that this new waveform be used in digital simulations. The latest UBC EMTP version Micro-. Tran has this function implemented. 15 2.3) Tower surge impedance and its response to lightning strokes The overvoltage created by a lightning stroke depends on the surge response characteris-tics of the object it strikes. Transmission towers are difficult to represent in lightning surge simulations. From a wave propagation standpoint, the tower is a non-uniform structure, al-though many models in the literature assume that it is uniform. Some non-uniform models have been proposed [19-24], but the mathematics involved in most of these models are com-plicated and not easy to implement in general purpose programs such as the Electromagnetic Transients Program (EMTP). With the EMTP, transmission towers have been modelled as pure resistances, as pure inductances, or as lossless constant-parameter transmission lines, but the accuracy of such representations has seldom been verified. This has motivated re-cent research on more accurate tower models [18]. It is relevant to mention that poor lightning performance has been reported on some transmission lines, even if they had very good shielding and low footing resistances [27]. Some lines with low footing resistances have actually displayed outage rates as high as the ones with high values [28]. These unexpected results may well be caused by wave propaga-tion effects up and down the tower. For example, Chisholm [28] pointed out that the influ-ence of the transmission tower on lightning performance of transmission lines may well have been overlooked in many lightning surge simulations. In the past, it has often been said that it is not practical to model transmission towers in detail, because there are already so many uncertainties involved in lightning surge simula-tions (see, for example, the discussions to [17]). Recent line performance data indicate that the actual lightning outage rates may be underestimated with simplistic models. It can also 16 be shown that the wave propagations up and down the tower after a back flashover have an important influence on the waveform which reaches the equipment in the substation. Better tower models are therefore needed. The transmission tower is a non-uniform transmission line structure, whose surge impe-dance (Z) varies as the" surge travels through it. Wagner and Hileman [13] initially derived and proposed a cylindrical tower model with the equation, Z=601n 2 l f (2.1) where t is the travel time counted from the instant the current has entered the tower top, c is the velocity of light and r is the radius of the cylinder. Equation (2.1) illustrates that the tower surge impedance is not constant; it is lowest at the top and increases as the wave progresses down to the base. Kawai [14] performed measurements on isolated towers (without ground wires connected) and obtained a trend similar to the Wagner-Hileman equation, although the magnitudes were appreciably lower. However, the change in surge impedance in [13] and [14] is opposite to that found in most subsequent experiments [15]. Recognizing the inconvenience of the time-variation in equation (1), Sargent and Darve-niza [16] developed a conical tower formula with a constant surge impedance, Z=601n (2.2) where S is the sine of the half-angle of the cone. Chisholm et al [17] later performed experiments as well as analytical calculations, and discovered that the position-dependent characteristics of the tower surge impedance are 17 dependent upon the direction of the injected current source. For a vertical stroke (current terminating at the tower top), the trend follows that of Kawai's and Wagner-Hileman's re-sults. The trend is the opposite for the horizontal situation - where the stroke current termin-ates somewhere on the ground wire and progresses horizontally as a travelling wave towards the tower. To accommodate the case of the horizontal wave, they devised another equation which characterizes that of an inverted cone, where 0 is the half-angle of the cone. Ishii et al [18] recently performed excellent low current measurements and ascertained that the presence of the connected ground wires alters the change of the measured impe-dance. Their results confirmed that the surge impedance is higher at the top and lower at the bottom when the tested tower has ground wires connected to it. These results will be used to verify the proposed nonuniform tower model, which is described in detail in Chapter 3. 2.4) Transmission line characteristics for lightning studies The modelling of transmission lines in digital simulations is well described in the litera-ture; e;g., [30] to [41], to mention just a few. Most studies in the past have relied heavily on single-phase representations.: Due to. major improvements of the line models in the last two decades, detailed.multi-phase line models, are now used more and more. However, despite these improvements, limitations still exist. There are some guidelines for choosing the exist-ing line models for typical Ughtning studies [5, 7, 31, 32, ...], but further improvements are (2.3) 18 still required in order to take into account the complete frequency-dependent nature of un-transposed overhead transmission lines. There are presently four main overhead transmission line models available in the EMTP, namely, the nominal Pi-circuit, the exact Pi-circuit^ the constant distributed parameter, and the frequency-dependent distributed parameter models [31], [36, pp. 4-66 to 4-107]. The as-sumptions used in these line models are clearly described in the EMTP Theory Book [36]. Both Pi-circuit models cannot be used for lightning or fast transient studies because they are valid only for steady-state or low-frequency transient analysis. If the constant-parameter model is used in lightning studies, the line should usually be modelled as untransposed, and the transformation matrix for transforming phase to mode quantities should be evaluated in the hundreds of kHz range. Unfortunately, the option of lumping R in a few places, which is used in many studies to take the line losses into account, is not usable for long lines in light-ning studies, since the total line series resistance at high frequencies is no longer much less than the line surge impedance. This is one of the important restrictions which the con-stant-parameter line model has. Ignoring the line losses altogether might be a legitimate ap-proximation for lightning studies, but this implies that all modes travel at the same speed of light. In reality however, the ground mode does travel somewhat slower (-5%) than the aer-ial modes in the lightning frequency range. If the studies are concerned with short line dis-tances (a few tower spans), this model is considered to be accurate enough, since the mode separation over such a short distance is not noticeable [36, p. 4-57]. A better high frequency lossless constant-parameter model has been suggested [42], which starts from a line para-meter calculation (possibly at 500 kHz) with losses included, but the series resistances are 19 set to zero just before the modal parameters are calculated. This gives a realistic wave vel-ocity for the ground mode. One limitation should be kept in mind, though, namely that the resistance for the ground return mode in lightning transients is highly frequency dependent [32]. For this reason, it is best to use the frequency-dependent model. However, even this model has the limitation of using constant transformation matrices. Constant transformation matrices are usually accurate enough for single and flat circuit lines. They are less accurate for double-circuit lines and probably worse for triple circuit lines. In general, the constant transformation matrix assumption is questionable if the asymmetry of the line conductor configuration is strong [33]. Thus, a full wideband frequency-dependent transmission line model is desired [43]. This line model should provide correct answers in both steady-state and transient studies. This main objective has been achieved in this thesis. The modelling details of this new transmission line model are presented in Chapter 5. It is worth mention-ing that a full frequency-dependent modal-domain model was developed earlier to achieve the same goal [35], but primarily for underground cables. In parallel with the work reported in this thesis, other phase-domain frequency-dependent transmission line models have also recently been proposed and published as well [38 - 40]. 20 CHAPTER 3 Modelling of Single-Phase Nonuniform Transmission Lines for Representation of Transmission Towers 3.1 Overview A line model with exponentially varying parameters is introduced to represent nonuni-form transmission lines in this Chapter. From the travelling wave standpoint, the trans-mission towers are viewed as nonuniform transmission line structures. When the line parameters are assumed to vary exponentially, a set of two-port equations can be formed in the frequency domain, which contain frequency-dependent functions. These functions are then synthesized with rational functions of the minimum-phase-shift type. Utilizing a fast recursive convolution technique, the time-domain equations of the proposed model reduce to a form similar to those in Bergeron's method. Thus, the model is compatible with general electromagnetic transients programs such as the EMTP [1]. In the early stage of this thesis work, a tapered line model [22] was developed for the representation of transmission towers. It had shortcomings, however, inasmuch as it was only a high frequency approximation. Subsequently, improvements were made with the ex-ponential variation assumption of [23]. With this assumption, the non-uniformity of the line parameters can effectively be included in the time domain simulation, as described in the following sections. 21 zc ( x ) = Z 0 e^ * -> <-length = 1 Figure 3.1 - Single-phase exponential line 3.2 Exponential transmission line Notations: The uppercase represents frequency domain quantities, while the lowercase indicates their time domain correspondents. 3.2a Equations in the frequency domain The basic equations of a nonuniform transmission line, expressed in the frequency do-main, are (3.1a) (3.1b) V and I are the voltage and current phasors, and Z(x) and Y(x) are the space-dependent per-unit-length series impedance and shunt admittance, respectively. Following the procedure of [23], Eq. (3.1a) is differentiated again, (3.2) d2V y ( . d l , dZ(x) dx2 v dx dx Substituting (3.1a) and (3.1b) into (3.2) gives 22 For the exponential line shown in Figure 3.1, with losses ignored, Z(x) and Y(x) are Z(x) = j(oL(x) = j®L0eqx, (3.4a) Y(x)=j(oC(x)=j(oC0e-"x, (3.4b) where L(x) and C(x) are the per-unit-length inductance and capacitance, respectively; L 0 and C 0 are the values at x = 0. These parameters are related to the high-frequency approximation of the line characteristic impedance, Z h i A x ) = l m =z°eqx> (3-5) where Z 0 = f e . Substituting Eq.'s (3.4a) and (3.4b) into Eq. (3.3) gives ^ - ^ + ^ F = 0 ' ( 3 - 6 ) where c = ^ is the wave speed. J LoC o Eq. (3.6) is a second-order differential equation with the roots x - f - J S ' - d ) 1 : x , . | + J f i ) , - ( f ) J . (3.7) Then, and from Eq. (la), V(x) = Cie^x + C2e^x, (3.8) I(x) = -^[^Cle^x + X2C2exn- (3-9). The constants C, and C 2 depend on the boundary conditions. They can be found by evaluating Eq.'s (3.8) and (3.9) at nodes k and m. At node m we have 23 Vm = C]e^l + C2e^', x2i 1 Im=-^=-[XiCle^l + X2C2ex% where Zm =jwL0eqx. Thus, from Eq. (3.10), Cx=\Vm-C2e^l\e-^1. Substituting Eq. (3.12) into Eq. (3.11) and collecting terms, A.) Vm +ZmIm (3.10) (3.11) (3.12) C2 X\ — X2 (3.13) Then C, is obtained from Eq. (3.12), X2 Vm 4~ ZmIm [_ X2 — X\ Evaluating Eq.'s (3.8) and (3.9) at node k leads to vk = c,+c2, 1 (3.14) where Zk = j®L0. h = -—[X\C\ + X2C2], (3.15) (3.16) Multiplying Eq. (3.16) with t— and adding it to Eq. (3.15), the voltages and currents at K2 both ends can be related as [Vk + Zck((o)Ik]A = Vm+Zcm((o)Im, (3.17) where A = eXxl is the propagation function. Z c k and Z c m are the characteristic impedances at both ends of the line, which can be expressed as j<nL0eq' Zck((a) =J^~ and Zcm((o)=- X2 (3.18) By letting co go to infinity in Eq. (3.7) and inserting the result into Eq. (3.18), the char-acteristic impedances become the high-frequency approximation of Eq. (3.5). In [22], this high-frequency approximation was used. 24 3.2b Equations in the time domain Reversing the current direction at node m to make it flow into the line, Imk = -Im, Eq. (3.17) in the time domain becomes M O + zck(t) * ikm(t)] * a(t) = vm(f)-zcm{t) * imk(i), (3.19) where the symbol * denotes convolutions. If the characteristic impedances Z c k and Z c m , together with the propagation function A, are synthesized with rational functions, their corresponding expressions in the time domain will become simple sums of exponential functions. The function a(t) will also have a time delay, which approximately equals the time it takes for the fastest frequency component to travel along the line. Accordingly, the convolutions in Eq. (3.19) can be evaluated with a fast recursive algorithm [34]. Thus, the voltage at node m can be expressed in a simpler form as vm(t) = zmeqimk(i) + vhRcm(i) + VhProm(f), (3.20) where z m e q is a constant. The last two terms on the right- hand side are evaluated from the known values of previous time steps. The first term comes from the RC-network which approximates the characteristic impedance, and the second term comes from the propagation of conditions at the remote end k. Combining these two terms into a single history voltage source eh.m(t), Eq. (3.20) becomes vm(t) = zmeqimk(f) + eh-m(f). (3.21) The same procedure is used to evaluate the voltage at node k, but with the wave direction from node m to node k instead. Factor q has opposite sign now, and the roots in Eq. (3.7) must therefore be re-evaluated. The characteristic impedances and propagation function of 25 km k ^ k e q -> I W v ~ e h - k ' + ^meq m V\AA——i <— h - m m Figure 3.2 - Equivalent circuit for the exponential line the line are also recalculated, using these new roots. Carrying out the necessary steps, the voltage at node k can then be expressed in the same form as Eq. (3.21), vk(t) = zkegikm(t) + eh-k(t). (3.22) Eq.'s (3.21) and (3.22) are very similar to those in Bergeron's method [1]. They are com-patible and can easily be interfaced with the EMTP, with the equivalent circuit of Figure 3.2. 3.3 Synthesis of exponential line functions Recursive convolutions leading to Eq. (3.21) and (3.22) are possible only if the line char-acteristic impedances and propagation functions can be synthesized successfully. The syn-thesis results are discussed in this section. The line parameters are synthesized using rational functions of real negative poles and zeros. A 50 m long exponential line with high-frequency characteristic impedances Z h i g h k and Z h i g h m of 220 and 150 Q at node k and node m, respectively, is used for demonstration purposes. 3.3a Travelling from node k to node m The factor q for the example is evaluated as 1 Zhighm 1 , 150 q I Z h i g h k 50 220 = -0.0076 26 Figures 3.3 and 3.4 show the magnitude and phase angle of Zck(co)'. Note that it is only necessary to calculate and synthesize Zck(co) because Zcm(co) can be found by multiplying Zck(co) with eql. Figures 3.5 and 3.6 show the magnitude and phase angle of the propagation function A. The fitting accuracy is quite reasonable in these figures, and the magnitude and phase of the synthesized functions trace those of the exact ones very closely. The fitting accuracy can be further checked by comparing the results directly in the time domain. Figure 3.7 shows the exact and synthesized propagation function a(t) obtained from taking the inverse Fourier transform of A(oo), using the Fourier transformation program of [55]. To filter out high frequencies, the Lanczos filter of the form s ^ n ( 7 i a ) / ^ ( ° m ax ) w a s a ppii e c j tQ the (710)/© max) synthesized and exact functions, where comax is the truncated frequency and is equal 27tl08 for this case. The two curves in Figure 3.7 are practically identical, which illustrates that the fitted a(t) is an excellent representation of the exact one. The propagation function a(f) has a physical meaning [36, p. 4-95]: If the line were to be energized with a voltage from node k which has an amplitude of 1.0 at all frequencies, through an impedance matching the line characteristic surge impedance Zck(co), and with node m open, the voltage V m at the receiving end would be the propagation function A(co). This, in the time domain, is interpreted as applying an impulse voltage (infinite magnitude with an area of 1.0) at node k through the characteristic impedance to the line. The voltage arriving at the open end of the line will then be the propagation function a(t). Note that the line under consideration is 50 m long, which means that it should take approximately 16.67 ps for the impulse to reach the other end of the line. This is indeed seen in Figure 3.7. 'it was found experimentally that the fitting of Zck(co) is better in the low frequency region when a small conductance G is included in Eq. (3.4b). In that case, the term (co/c)2 in Eq.'s (3.6) and (3.7) becomes - Z Y . 27 x 1 0 • i • i Exact function 400 1 — • 1 • 5 350 a -4 k \ 300 250 o o o -3 1 200 150 -2 1 00 i o 4 1 0 5 6 7 10 10 1 o8 1 > E n l a r g e m e n t n 1 0 - 2 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 3.3 - Magnitude of Z c k in Q 5 0 0 5 0 - 1 0 0 2 0 1 0 1 0 Exact function Synthesized function 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 3.4 - Phase of Z c k in degrees .28 0 . 8 5 0 . 8 0 . 7 5 0 . 7 0 . 6 5 - 2 1 0 Exact function Synthesized function 1 0 1 0 1 0 Frequency (Hz) Figure 3.5 - Magnitude of A 1 0 1 0 2 0 0 - 2 0 - 4 0 - 6 0 - 8 0 1 0 0 1 2 0 1 0 Exact function Synthesized function 1 0 E n I a r g e m 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 3.6 - Phase of A in radians 29 . 2 l . —. .—- .—. 1 1 5 1 6 1 7 1 8 1 9 2 0 Time (us) Figure 3.7 - Propagation function a(t) 3.3b Travelling from node m to node k As mentioned previously, when the wave travels in reverse direction, a different factor q must be used in the equations, which results in different characteristic impedances and propagation function. The prime symbol will be used to identify these functions, namely Z ' c k , Z ' c m , and A'. They are used in Eq. (3.22) to evaluate the voltage at node k. The factor q in this case is the negative value of what it was for the wave travelling from k to m, namely l l n | ^ ! L = 0.00766. Figures 3.8 and 3.9 show the magnitude and phase angle of Z ' c m . The function Z ' c k is re-lated to Z ' c m with a multiplication factor of eql. On the other hand, A' can be evaluated di-rectly from the propagation function found for the other direction. 2 0 0 1 5 0 . 1 0 0 5 0 0 1 0 Exact function Synthesized function 1 0 1 0 1 0 Frequency (Hz) Figure 3 . 8 - Magnitude of Z ' c m in Q 1 0 1 0 1 0 0 8 0 6 0 4 0 2 0 0 - 2 1 0 Exact function Synthesized function 1 0 1 0 1 0 Frequency (Hz) 1 0 1. 0 Figure 3 . 9 - Phase of Z ' c m in degrees 31 Zh ighm (3.23) It can be seen that A' differs from A only by a constant factor. Thus, once A is syn-thesized, the approximate function for A' is also known! The good agreement between the fitted and the exact A functions applies to A' as well. 3.4 Reflection factors To gain confidence in the proposed line model, it is useful to calculate its input reflec-tion factor and compare it against previously published results. Eq.'s (3.8) and (3.9) are the starting point for this comparison. The impedance looking into the line at any point x, Z in(x), can be computed by simply taking the voltage-to-current ratio from these two equa-tions. Therefore, the impedance looking into node k, at the sending end, is evaluated with x = 0 as Z,„(0) = -zk Ci+C2 _ X\C\ + X2C2 Using the constants from Eq.'s (3.13) and (3.14), (A.2 Vm + ZmIm)e-^1 -(XxVm+ ZmIm)e-^1 Zin(0) = Zk kifri Vm + ZmIm)e-*1 -Xi(X2Vm+ ZmIm)e Z c(x)= z0eqx length = 1 Figure 3.10 - Single-phase exponential line terminated at node m with Z L 32 Letting the exponential line be terminated at node m with a load impedance Z L as shown in Figure 3.10, then V m = Z L I m , and.Zin(0) becomes Dxe'x'1 -D2e~X2' Zin(0) = Zk where _X2D2e-k2' -XtDie-*-1'. (3.24) Di=k2ZL + Zm and D2=XXZL + Zm. (3.25) Note that Z k and Z m are the per-unit-length series impedances evaluated at node k and node m, respectively. Figures 3.11 shows the magnitude of the input impedances evaluated from Eq. (3.24) for two different values of Z L . It is interesting to note that this input impedance is exactly the same as the load impedance Z L , at frequencies less than about 105 Hz. This is because the propagation function is purely real below this cut-off frequency, which presents only attenu-ation in the electromagnetic field without any progression (no travelling wave). Under this condition, the load impedance Z L is seen as if it were connected right at the sending end of the line. When travelling waves enter node k, the reflection coefficient seen at this node is de-fined as . ™ = t f ^ i ? <326> with the assumption that a uniform line of characteristic impedance Z 0 is connected to the sending end at node k. In [29], the reflection coefficient is approximated as 33 where (3 = Eq. (3.27) was derived from a first-order approximation for the reflection coefficient be-tween two infinitesimally short sections. It is therefore approximate and not as accurate as Eq. (3.26). In addition, Eq. (3.27) was established only for the case when the load impe-dance Z L and the high-frequency characteristic impedance at the receiving end m are nearly equal to each other [29]. For the exponential line considered in section 3.3, with the load impedance Z L = Z h i g h m = 150 Q, Figure 3.12 shows the magnitudes of the reflection coeffi-cients evaluated from Eq.'s (3.26) and (3.27). Both methods produce close answers for the line under consideration. The differences become larger when the Z h j g h k / Z h i g h m ratio increases to greater than 2. If Z L is set to 100 Q instead, Eq. (3.27) is no longer valid, as shown in Fig-ure 3.13. 3 5 0 8 0 1 0 1 0 . 1 0 1 0 Frequency (Hz) Figure 3.11 - Magnitude of Z i n in Q from E q . (3.24) 34 0 . 2 1 0 1 0 1 0 1 0 1 0 Frequency (Hz) Figure 3.13 - Magnitude of Reflection factor for Z L = 100 Q 35 3.5 E M T P simulations To show the behaviour of the proposed nonuniform line model, EMTP simulations of the line responses under open- and short- circuit conditions were performed. The results are shown in Figures 3.14 and 3.15. For testing purposes, the equations of the exponential line were interfaced with the UBC EMTP version MicroTran through subroutine CONNEC. The interface is described in Appendix C. Additional time-domain simulations for comparisons with field measurements will be shown in Chapter 6. An exponential line can be approximated as a cascaded line model [29], where the line is divided into many sections of uniform transmission lines. Figures 3.14 and 3.15 show si-mulation results for open and short circuit tests for this approximation as well. The expo-nential line, with parameters described in section 3.3, is energized with a 1 p.u. step voltage at node k. Three models are used for comparison purposes: 1) the new model, 2) a cascaded model consisting of 50 short uniform lines, and 3) an ideal uniform line having a uniform characteristic impedance equal to the average value (185 Q). Note that the proposed new model agrees well with the cascaded model for the open circuit voltages. It clearly shows the tapered effect created by the continuous reflections along the line, in contrast to the square waves obtained from the uniform transmission line model. The voltage obtained from the new model has larger peaks, compared to those produced by the cascaded model, because the latter model does not have enough short sections to represent these peaks. The short-circuit currents agree more or less among the three models during the first 1.5 ps. The current from the new model then starts to deviate from the others because of the continuous reflections along the line. 36 3 . 5 3 2 . 5 2 1 . 5 1 0 . 5 0 - 0 . 5 - 1 - 1 . 5 , Cascaded model Uni form model 0 . 0 7 0 . 0 6. 0 . 0 5 . 0 . 0 4 . 0 . 0 3 . 0 . 0 2 . 0 . 0 1 . 0 New model 0~5 1 1 . 5 2 2 T Time (u.s) Figure 3.14 - Receiving end open-circuit voltage (p.u.) Cascaded and Uni form models N e w model o: 5 1—\— r^ T" Time (us) Figure 3.15 - Receiving end short-circuit current (p.u.) 37 CHAPTER 4 Frequency-Dependent Transformation Matrices for Untransposed Transmission Lines using Newton-Raphson Method 4.1 Overview The development of a complete frequency-dependent overhead transmission line model in the time domain is one of the major objectives of this thesis. In the frequency-dependent overhead transmission line model, presently in the EMTP, constant and real transformation matrices are used for conversion between modal and phase domain quantities. As men-tioned in section 2.4, constant transformation matrices are only accurate enough for single flat circuit lines. For the case of double- or more circuit, lines, the accuracy decreases noti-ceably, especially in studies which cover wide frequency ranges [33]. Lightning is a phe-nomenon which, in most cases, involves a wide range of frequencies. For example, in the simulation of the shielding-failure flashover and back flashover events, the 60 Hz steady-state as well as the transient voltages are needed. Thus, a wideband frequency-dependent line model is needed for such situations. One such model is the full frequency-dependent modal-domain model of [35]. It includes the complete frequency-dependent nature of the transmission line and synthesizes not only the modal characteristics but the elements of the transformation matrices as well, with rational functions. To achieve this, the elements of these matrices must be continuously and smoothly evaluated over a wide frequency range 38 (1 Hz - 1 MHz). A stable eigenvector/eigenvalue routine based on the Jacobi algorithm was developed and used in [35]. The smoothness and minimum-phase requirements for the transformation matrix elements were satisfied, but only for the case of underground cables. For multi-circuit overhead lines, however, these elements could not be fitted with rational functions of minimum-phase shift. In addition, they did not always vary as smoothly as one would wish [45]. Realizing this limitation, a diagonalization routine was developed in this thesis project, together with an effective scaling method, which produces well-behaved transformation matrices. The time-domain implementation procedure can then follow the same formulation as described in [35]. The diagonalization method is based on the New-ton-Raphson iteration technique. It was first proposed by L . M . Wedepohl [48], and initially programmed by G.D. Irwin at the Manitoba HVDC Research Centre. The method was then additionally improved and implemented successfully with the U B C Line Constants Program by the author [37]. This algorithm produces well-behaved transformation matrices for over-head transmission lines. Although some elements of the transformation matrices may pres-ent non-minimum phase behaviour, they could still be synthesized with rational functions which have real negative poles and zeros, with reasonable engineering accuracy. 4.2 Equations for multiphase transmission lines Notations: Boldface represents matrix quantities, whereas uppercase depicts frequency-domain quantities. Considering a multiphase transmission line system as shown in Figure 4.1, the basic transmission line equations in the frequency domain are: 39 Node k Node m Ik legnth = l I m > > n • > Phase #1 • n V Phase #2 o —> Phase # n Figure 4.1 - Multiphase transmission line system - f = ZI, (4.1.). - ^ = Y V , (4.1b) ax where Z and Y are the per-unit-length series impedance and shunt admittance matrices, re-spectively. V and I are voltage and current vectors (or column matrices). dV Differentiating Eq. (4.1b) again and replacing — with Eq. (4.1a), we get for the line ax currents, V = [YZ]I. (4.2) ax It is well known that Eq. (4.2) is solved by transforming the n coupled equations into n decoupled equations. To achieve this, a modal transformation matrix T„ which diagonalizes the Y Z product, is found from the solution of the eigenproblem, T I " 1 Y Z T i = A,, (4.3) where A, is the diagonal eigenvalue matrix. The matrix T, varies with frequency since Y and Z are frequency dependent. 40 In order to understand fully the frequency-dependent characteristics of the matrix T„ it is best to examine the frequency-dependent nature of Y and Z first. The matrix Y has the fol-lowing form [36, pp. 4-16 to 4-17], Y = G+yooP-', (4.4) where co = 2nf is the angular frequency, G is the per-unit-length conductance matrix, and P is the constant Maxwell potential-coefficient matrix, having each element defined as s0 = permittivity of free-space, Djj = distance from conductor i to image of conductor j , dy = direct distance form conductor i to conductor j , or conductor radius when j = i. The conductance matrix G is typically represented as a diagonal matrix. G is small com-pared to jcoP"1, and is therefore often ignored. In [34] and [36, p. 4-99], values of 0.03 uS/km were used. The matrix Z has the form of Z=y'cop0s0P + Ze + Z c , (4.6) where is the permeability of free-space. Subscripts c and e refer to the internal impedance of the conductor and earth impedance corrections, respectively. Note that Z c is a diagonal matrix. For solid conductors, its elements can be evaluated with a simple but sufficiently accurate formula in the frequency range of interest here [49], (4.7) -f PO Zcii = [ ^ - j I ( r 0.777A ^ n , c A ' i —— icothi — I +0.356 \2de\J V de\ J 41 where: p, = conductor resistivity, r = conductor radius, del = I , conductor complex depth of penetration, p = conductor permeability. At very low frequencies, d . approaches infinity, and Zcu = -—r, which is the correct d.c. Pi resistance. At very high frequencies, d e l approaches zero, and Zcu —> -——, which is the 2nrde] correct skin effect impedance at high frequencies. The formula for evaluating the earth correction impedance matrix Z e was developed by Carson [50]. This formula has been used by power engineers for many years. Deri etal evaluated the elements of this matrix with a simpler formula [51], which was originally pro-posed by Dubanton [52]. It has the following form [36, pp. 4-4 to 4-5] and [41], 2 * - ^ * $ . (4.8) lit JJii where D'j = ^(y,+y, + 2afe2)2 + (*,• -Xj)' The term dei = / ~— is the earth complex depth of penetration, x and y are the hori-zontal and vertical co-ordinates of the conductors, respectively, and p 2 is the earth resistivity. Eq. (4.8) is accurate enough for engineering purposes, and the results agree very closely with those obtained from Carson's integral formula. The maximum error is only a few percent in the kHz range. This is of very little consequence in transient calculations where the entire frequency spectrum contributes to the results anyhow. Thus, Eq. (4.8) is preferred in this thesis because of its simplicity. It should be mentioned that the transformation matrices 42 evaluated in this thesis, using Carson's formula or Eq. (4.8) in the Y Z product, were practi-cally identical. 4.3 Eigenproblem solution with Newton-Raphson (NR) method For convenience, the following additional variables are introduced: i) Sis the Y Z product, ii) T I ( k ) is the column k of matrix T , , iii) U is the unit matrix. Then, equation (4.3) can be expressed as S T , = m . (4.9) Eq. (4.9) contains n 2 equations for an n-phase system, because each of the n eigenvalues Xk is related to n equations. Considering one set of n equations, which belongs to a particu-lar eigenvalue Xk [41], i.e.: (S-XkU)Tm = 0. . (4.10) Eq. (4.10) is a set of homogenous simultaneous equations. It may be shown that the sol-ution is trivial unless the determinant of (S - A,*U) is zero. Thus, the matrix S can be diago-nalized by equating this determinant to zero. This leads to a polynomial in X of order n, yielding n values of X. Sometimes one or more values of X may be equal. Once the values of X are known, the elements of T I ( k ) may be found from Eq. (4.10). There are various software packages available for determining the eigenvalues and ei-. genvectors. However, for the modal-domain transmission line model of [35], the eigenva-lues and eigenvectors must be found at many frequencies. The elements of the eigenvectors 43 have to vary smoothly as a function of frequency throughout the entire frequency range. This is very important because rational approximations for these elements are required. If a normal algorithm is used, the eigenvalues will switch places at certain frequencies, which creates jumps in the eigenvectors. A complicated "tracking procedure" was used in [35] to avoid jumps. The NR algorithm proposed here overcomes this limitation. It solves Eq. (4.10) itera-tively for each value of A,k together with each set of T , ( k ) . This system of non-linear equa-tions has "n +1" unknowns. The number of variables could be reduced by one, by specifying an arbitrary value for one of the elements of T I ( k ) . For example, the element with the largest value at the initial frequency could be chosen and set equal to 1.0 for all fre-quencies. However, this could result in overflow or excessive errors, because the value of this element may sometimes become zero or very small at some frequencies. A safer way is to introduce a constraint equation. One such constraint is to specify that the sum of the squares of the elements of the eigenvector be unity. This constraint method has been chosen in this thesis. With this constraint, the Jacobian matrix is well behaved, and the rational ap-proximations of the resulting smooth transformation matrices are found with ease. For example, for a three-phase line, one set of 3 equations belonging to A,, and column 1 of the matrix T, from Eq. (4.10) are (5 , i -Xi )T m +SuTm +SnTDl =0, (4.11) 5217-/1,+(522-^)7/21+523^1=0, (4.12) 53ir/„+532r/2,+(533 -Al)T/31 =0. (4.13) 44 There are four unknowns but only three equations. The extra equation is the constraint equation, which specifies that the sum of the squares of the elements of T I ( 1 ) be 1.0, (4.14) T2m+T2m+T2m-l=0. The Jacobian matrix for the above four equations then has the form, J = (S]\-X\) Sn S\3 - T / i i 521 (^22-A-l) 523 -721 531 532 (533 - ^ l ) -T31 27/n 27/2i 27/3i 0 (4.15) Let x represent the column vector of the four unknowns Tn\,Ta\,Tn\ and V Then x in iteration step n is found from x ( n ) _ x ( " - l ) . j(X(«-i)) ' F ( x ( " - 1 ) ) (4.16) where J(x)" 1 is the inverse of the Jacobian matrix and F(x) represents the left-hand sides of equations (4.11)- (4.14), which eventually become zero as the iterations converge. The iter-ation procedure starts from an initial guess x ( 0 ) . As starting values of x, the solution values from the previous frequency are used. The NR method generally converges very fast, provided that reasonable starting values of x are known and that J(x)" 1 exists. After convergence is reached, A, and column 1 of the matrix T, are found. The same procedure is then followed for X.2 and column 2 of T , , and then for X3 and column 3 of T , : It should be noted that this NR method is not self-starting. Another method must be used to supply the solution for the initial frequency first. In the program written for this the-sis, the diagonalization routine based on the power method provided by L . M . Wedepohl [41] is used to start at the initial frequency. In addition, in order to avoid floating-point overflow 45 and to ensure continuity of the modal parameters at the high frequency region, the Y Z prod-uct is scaled before the NR diagonalization algorithm is called. The scaling procedure con-sists of dividing each element of this product by a constant -(co2p0s0), and subtracting a unit matrix (1.0 for diagonal and 0.0 for off-diagonal elements) from the resulting matrix. As seen in Appendix A, this scaling procedure does not alter the eigenvectors of the original Y Z product. This scaling procedure was suggested in [54], and is further explained in [36, pp. 4-52 to 4-53]. 4.4 Modal parameters for multiphase transmission lines from Newton-Raphson method Modal parameters were found with the NR method for different line configurations, ranging from single- and double- circuit lines for horizontal and vertical constructions to un-common circuit lines with the phase conductors arranged in various asymmetric settings, with or without ground wires. An example of an uncommon double-circuit line is the case where the left three phases are vertical while the right ones are horizontal. Another case with strong asymmetry, which is common in practice is the case where two or more circuit lines of different voltage ratings share the same right-of-way and may be located one under-neath the other. The NR method had no difficulty in diagonalizing such uncommon line configurations, and the transformation matrices were all smooth and well behaved. The fitting routine used in these tests, referred to as FIT-S, is a modified version of the one described in [35]. This fitting routine was developed by J.R Marti [34]. The main modification involved the calculation of the time delay constant for the rational approxima-tion of the modal propagation function. Since the angle of this function can be found 46 directly, and since it is smooth and continuous, only one time delay at the fastest frequency was used, whereas the average time delay over a frequency range was used in [35, 45]. The synthesis and fitting criteria for all modal functions are explained clearly in references [34] and [35]. They will not be repeated here. It is significant to mention that the NR method produces modal parameters which are smooth and vary only slowly as a function of frequency, because of the inherent eigenva-lue-eigenvector tracking nature of the NR algorithm. This is the reason why syntheses using minimum-phase-shift rational functions have been successful in this thesis (FIT-S has not failed in fitting the parameters produced by the NR method for all the line configurations tested in this thesis). 4.4a M o d a l transformation T , matrix There are n modes in an n-conductor transmission line. In the simulation of lightning tower/ground wire surges and back flashover events, ground wires are included as conduc-tors, whereas in shielding failure and shielding-failure flashover events they may be elimin-ated from the Z and Y matrices by Gaussian elimination. For the double-circuit line of Figure 4.2, assuming that the ground wires have uniformly zero potential and can thus be eliminated, there are six modes describing this line. Figures 4.3 shows the magnitude and phase angle of the elements of column 1 of the 6 x 6 transformation matrix. Note that there are actually twelve curves superimposed in each of these, plots. It is seen that the fitting functions match the exact functions closely. The elements in other columns also match quite reasonably, but some differences in the angles showed up in the low frequency region. 47 11.1 m 12.7 m 12.7 m 7fT _ \ /_ A 26.0 m , 20.8 m . . 19.0 m . r n 1 O 6 4 c 19.0 m 2 O O 5 20.8 m 3 O Ground wire: 1.75 cm diameter length = 100 km Phase conductor: 4 X 3.84 cm diameter at 50 cm spacing o- o ! I o—o • i i 50 cm 100 Ohm-m ground resistivity o 6 Figure 4.2 - Vertical double-circuit line These differences, however, have very little bearing on the overall transient solution, as proven in the open/short -circuit tests in section 4.5. Because of symmetry, the elements of each column of the matrix T, occur either in equal or in equal and opposite pairs. The phase conductors are arranged such that elements 1 and 4, 2 and 5, and 3 and 6 of each column are equal in magnitude. It is shown in Appendix A that the elements of the matrix T, become real and approach asymptotic values at the high and low frequency ends. This characteristic is seen in Figure 4.3. 4.4b Modal characteristic admittance matrix Y, C M The modal shunt admittance and series impedance matrices are defined as Y M = T7'YT7 r and Z M = T [ Z T i , where T,"T is the inverse of the transposed matrix T, T . It follows from Eq. (4.9) that (4.17) Y M Z M = A,. (4.18) 48 It has been proved that both Y M and Z M are diagonal [53]. The modal characteristic ad-mittance Y C M is also a diagonal matrix and can be found from the scalar relationship, YcM-k (4.19) where Y C M . k is the diagonal element k of the Y C M matrix. Figure 4.4 shows the six elements of the modal characteristic admittance matrix and their corresponding synthesized functions for the line in Figure 4.2. Note again that the two fitted and exact curves are superimposed in the figure. 4.4c Modal propagation matrix A M The modal propagation matrix A M is also diagonal. Its elements have the form of AM-k = e-J1*1 = e-ale-W, (4.20) where Xk is the k-th eigenvalue and 1 is the,length of the line. The term e"al is the attenuation factor and e"jpl the phase shift. Figure 4.5 shows the superimposed (exact and synthesized) curves of these matrix elements for the line in Figure 4.2. There are only five distinct el-ements shown in each of the plots since two of the eigenvalues are very close or identical for this line. The phase angles, without time-delay phase shifts, are shown as well. Again, the fitting is very good. In order to ensure continuity in the travel time delays, it is important to obtain the phase angle of each element A M . k directly as arg{AM-k} = Im{JTk}l. (4.21) 49 With this equation, the angle is smooth and continuous, as shown in Figure 4.6. If it were evaluated as the phase angle of the function 1 , using the intrinsic function atan in FORTRAN, it would show highly oscillatory behaviour in the high frequency region, as shown in Figure 4.7. It would be difficult to evaluate the time delays needed for the fitted functions with such a behaviour in the phase angles. This aspect will be elaborated further in Chapter 5. FIT-S usually allocated fewer poles for the elements of Y C M and T, than for those of A M . On the average, twelve poles were used for the former and 28 for the latter. 50 0 . 8 <D 0 . 6 3 CD ro 0 Exact function Synthesized function Maximum error < 1.0 % 3 & 6 1 & 4 1 0 1 0 1 0 Frequency (Hz) 1 .0 2 0 1 0 in T3 D) c ro (U (/> ro 0 E 1 0 2 0 1 0 1 0 1 0 1 0 Frequency (Hz) Figure 4.3 - Elements of column 1 of the matrix T , 51 , X 1 0 b , , —. 1— # 2 0 1 0 1 0 1 0 Frequency (Hz) Figure 4.4 - Elements of the matrix Y , 52 1 0 1 0 1 0 1 0 Frequency (Hz) Figure 4.5 - Elements of the matrix A , 53 x 1 0 m c ro ~o ro _0) TO C ro w ro 1 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 4.6 - Phase angle of element A M _, calculated from E q . (4.21) ro a> a) c ro <u w ro o k 1 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 4.7 - Phase angle of element A M . , calculated with atan function 54 4.5 Open/Short -circuit tests in the frequency domain An efficient and simple way of comparing the solutions produced by the exact and ap-proximate modal parameters is to calculate the line responses under unfavourable condi-tions, such as short- and open- circuit tests, in the frequency domain. In order to perform these tests, a column voltage matrix E s is applied at the sending ends (node k) with the re-ceiving ends (node m) either open or shorted. The single-phase open-circuit voltage and short-circuit current for each mode at the receiving ends are evaluated as follows [45], ^ ( « » = T ^ # 1 . (4-22) ji , x 2YcM-iAM-iE's_i /m_/(co) = — — , (4.23) 1 ~Am-I where Y C M . j and A C M . j are the i-th elements of the modal characteristic admittance and propa-gation function matrices, respectively. The prime symbol indicates that the voltages and currents are in the modal domain. The conversion between modal and phase domain for the voltage and current column matrices is made with the following relationships, Im(co) = T,(co)I'w(co) and Vm(co) = T7r(coVVUco). ( 4 . 2 4) The exact line responses are calculated using the exact T , , A M and Y C M matrices, where-as the approximate responses are obtained using their respective approximations by rational functions shown in section 4.4. In the following, two forms of excitation will be used to test the double-circuit line of Figure 4.2. 55 4.5a Zero sequence excitation For this form of excitation the source voltages all have 1 p.u. magnitude and zero angle. Figures 4.8 and 4.9 show the magnitude of open- and short- circuit responses for all phases at the receiving end. Note that the approximate and exact responses agree very closely. The angles are not shown but they also agree closely. 4.5b Unsymmetric excitation The accuracy of the line responses with approximate modal parameters depends on the degree of imbalance in the voltages and currents. An unfavourable case of voltage source imbalances is used to test the line considered here. This case was suggested in [33], with these values: Es\ =ES4 = IZ0°, Es2 =Ess =0, and E& =ES6 = 1Z1200. Figures 4.10 and 4.11 show the magnitude of the responses for all phases at the receiving end. Note again that the approximate and the exact solutions match one another almost perfectly. These tests were repeated for many other line configurations. The degree of accuracy is generally very good for the cases of single- and double-, flat and vertical, circuit lines. For lines with strong asymmetry, the approximate and exact solutions still agree quite well when the excitation sources are unsymmetric or symmetric. For the zero sequence excitation, small differences appear at some frequencies in the medium frequency range. For all practi-cal purposes, these small differences are still considered to be acceptable in the simulation of power system transients. Exact function Synthesized function -A. 1 0 8 6 4 2 -r-~ n r 1 0 1 0 Frequency (Hz) 1 0 1 0 Frequency (Hz) 1 0 1 0 Frequency (Hz) 1 0 - 2 0 2 4 6 1 0 1 0 1 0 1 0 1 0 1 0 1 ( Figure 4.8 - Receiving end open-circuit voltages (p.u. magnitude) 57 Figure 4.9 - Receiving end short-circuit currents (p.u. magnitude) °3 co a) CO ro CL 7 0 6 0 5 0 4 0 3 0 2 0 1 0 L O CN w w ro 1 0 5 0 4 0 -3 0 -1 0 5 0 4 0 -CO oa C O CO (U CO ro Exact function Synthesized function 1 o L l 1 o 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Frequency (Hz) 1 0 Figure 4.10 - Receiving end open-circuit voltages (p.u. magnitude) 1 0 1 0 . 1 0 1 0 1 0 1 0 Frequency (Hz) 0 . 2 : : \ Figure 4.11 - Receiving end short-circuit currents (p.u. magnitude) 60 4.6 Validation of the Newton-Raphson algorithm with the power method An alternative method was used to validate the results obtained from the NR algorithm. This alternative method uses a diagonalization routine based on the power method. Since the eigenvalues will most likely switch places at certain frequencies, a correlation technique is used to keep track of the eigenvectors throughout the frequency range. This correlation technique is based on the fact that the eigenvectors belonging to the same set of eigenvalues are nearly orthogonal from one frequency to the next. This means that if an eigenvector at the present frequency is to be correlated with each and every eigenvector obtained at the previous frequency (the most recent frequency), it will be found that the strongest correlation factor will be produced by the two eigenvectors belonging to the same set of eigenvalues. The correlation procedure for tracking the eigenvectors can be summarized as followed: i) Obtain the eigenvector matrix T I 2 for the present frequency. ii) Transpose T I 2 and obtain its complex conjugate T I 2 . iii) Find the matrix product T I 2 T n , where the eigenvector matrix T„ is the one from the previous frequency. iv) In each row of the product in iii, identify the largest element. The row number of this element indicates the column number of the eigenvector from the previous frequency, where as its column number represents the column number of the eigenvector at the present frequency. For example, in the case of a 3 by 3 matrix, if eigenvalue switching has not oc-curred, the correlation product of iii would have the largest elements in the diagonal loca-tions as indicated by the X's in Figure 4.12a, whereas if the eigenvalues of elements 1 and 3 have been switched, there will also be a switch between the eigenvectors of columns 1 and 3, 61 and the correlation matrix in iii would have the largest elements in the locations shown in Figure 4.12b. At every frequency, after the eigenvectors are obtained from the diagonalization routine, they are correlated with those obtained at the previous frequency. When switchings among eigenvectors are detected, they have to be rearranged in the order of the eigenvectors from the previous frequency. This guarantees their continuity throughout the frequency range. Tests confirmed that the modal transformation (eigenvector) matrices obtained from the NR method and from the above alternative method agree perfectly. Thus, the validity of the proposed NR method was proven. X X X X X X (a) (b) Figure 4.12- Location of largest elements of the correlation matrix 62 CHAPTER 5 Frequency-Dependent Overhead Transmission Line Model in the Phase Domain 5.1 Overview As an alternative to the modal-domain transmission line model described in reference [35] and in Chapter 4, the model presented here works directly in the phase domain. This approach has the advantage of not having to synthesize the transformation matrices, and also uses fewer numerical convolutions at each time step. Thus, it is computationally more effi-cient than the modal-domain models of untransposed lines. Direct phase-domain modelling has been proposed in [38] - [40]. The approach devel-oped in this thesis differs from these references in that it synthesizes the frequency-depend-ent line functions directly in phase co-ordinates, using stable rational function approximations with real negative poles and zeros. This new wideband transmission line model has been tested extensively for various overhead lines with different and asymmetri-cal geometries. It is stable and provides answers which are accurate enough in both steady-state and transient conditions. 5.2 Phase-domain transmission line modelling Notations: Bold and normal faces represent matrix and scalar variables, whereas upper-case and lowercase are used for frequency-domain and time-domain quantities, respectively: 63 Considering again the multiphase transmission line system of Figure 4.1, the classical solution for Eq. 4.2 in the frequency domain is Ix = e-rxIf+erxlb, (5.1) where T = J\Z is the propagation function, and I f and Ib stand for the vectors of forward and backward travelling currents, respectively. The voltage vectors of the line can be found as V* = - Y " 1 ^ = Y c V r * I / - <?rxh), (5.2) where the characteristic admittance Y c is Y c = J(YZTl Y (5.3) Pre-multiplying Eq. (5.2) with Y c and adding the result to Eq. (5.1) gives Y c V , + I x = 2e-r jcI/. (5.4) Applying boundary conditions at both ends, At node k, x = 0: Y c V k + I k = 2I/, (5.5) At node m, x = l: Y c V m - I m = 2e-rlIf. (5.6) In Eq. (5.6), I m is the current flowing into the line, or !(*=[) = - I m . Substituting Eq. (5.5) into Eq. (5.6), Y c V m - I m = A ( Y c V k + I k), (5.7) where A = e~rl. (5.8) Eq. (5.7) can be rewritten as I m = Y c V m - A ( Y c V k + I u). (5.9) Analogously for the currents at node k, Ik = Y c V k - A ( Y c V m + I m ) . (5.10) 64 Y c and A in Eq.'s (5.3) and (5.8) are the line phase-domain characteristic admittance and wave propagation matrices, respectively. They can be evaluated using modal parameters ob-tained from the eigenproblem of Eq. (4.3) [41, 35], Yc = T , Y c M T , r , (5.11) A = T , A M T 7 ' . (5.12) Note that the modal characteristic admittance matrix Y C M and the modal propagation ma-trix A M are calculated from Eq.'s (4.19) and (4.20). Transforming into the time domain, Eq.'s (5.9) and (5.10) become i « ( 0 = y c ( 0 * v « ( 0 - a ( 0 * f * ( 0 . (5-13) i*(0 = yC(0 * v*(f) - a(0 * fm(t), (5.14) where f*(0 = y c(0*v*(0 + u(0, (5-15) f«(0 = y c(0 * v „ ( 0 + im(t), (5.16) and the symbol '*' denotes matrix-vector convolutions. If the elements of Y c and A are synthesized using rational function approximations, then the corresponding time-domain functions yc(t) and a(t) become matrices whose elements contain simple sums of exponential functions [34, 35]. The convolutions in Eq.'s (5.13) and (5.14) can then be evaluated using fast recursive methods [30, 34, 35]. In terms of the syn-thesis functions, the nodal equations at both ends of the transmission line can be expressed in the more convenient form (Appendix B), im(f)=yiqVm(t)-ihta-m(t), . (5.17) i*(f) = YeqVk(t) ~ his,-k(t), (5.18) 65 k m Figure 5.1 - Multiphase transmission line equivalent circuit in the time domain where ye q is a real, constant symmetric matrix. The second term on the right-hand side of Eq.'s (5.17) and (5.18) is the history current source evaluated from the known values of the variables at previous time steps. These equations are represented with the equivalent circuit of Figure 5.1, which can easily be incorporated into the EMTP. 5.3 Synthesis of the propagation and characteristic admittance matrices The NR diagonalization algorithm described in Chapter 4 was used to evaluate the phase-domain characteristic admittance and propagation matrices as functions of frequency. The syntheses were then performed using the routine FIT-S. 5.3a Synthesis of the characteristic admittance matrix Y c The characteristic admittance matrix Y c in the phase domain is a square symmetric ma-trix. Therefore, only n(n+l)/2 elements need to be synthesized. Experience in synthesizing this matrix shows that all of its elements can be fitted using minimum-phase-shift functions, that is, using only real negative poles and zeroes. The only detail is that the off-diagonal el-ements need to be multiplied by -1.0 before the synthesis is performed. The reason is that, 66 with this multiplication, the value of the modified function asymptotically becomes a real constant with zero phase angle at the high frequency end (Figure 5.2), which allows FIT-S to allocate the constant R 0 (the function's high frequency limit), the poles and the zeros with better accuracy. Figure 5.2 shows the magnitude and. phase angle of the elements in column 1 of Y c for the vertical double-circuit line shown in Figure 4.2. Because negative functions of the off-diagonal elements were synthesized, the phase angle of the off-diagonal elements shown there differs by 180 degrees from the original functions. Note that the syn-thesized and exact functions match very closely, with errors typically under 2% within the required frequency range (1 Hz - 1 MHz). If the off-diagonal elements were fitted directly without the multiplication with the factor -1.0, noticeable differences would show up be-tween the exact and synthesized functions at the low frequency region. 5.3b Synthesis of the propagation matrix A The intention is to synthesize all the elements of the propagation matrix A with minimum-phase rational functions, and associate a different time delay with each of them. To achieve this, some function manipulations are required. First, the modal diagonal matrix A M is evaluated from Eq. (4.20). Secondly, a phase shift is then extracted to give P M = ^ A M ; (5.19) where x is the travel time of the line based on the speed of light. Finally, a triple matrix product is calculated, P = T i P M T 7 ' . (5.20) 67 The resulting matrix P is then synthesized using rational functions. Each of the syn-thesized elements contains an additional (small) phase shift which corresponds to an addi-tional time delay in the time domain. As in the case of the matrix Y c , the off-diagonal elements of the matrix P must be multiplied by -1.0 before the synthesis is performed. Figure 5.3 shows the fitting of the magnitude and phase angle for the elements of column 1 of the matrix P, for the double-circuit line of Figure 4.2. It is seen that the exact and fitted magnitudes agree closely, and the phase angles are in good agreement as well. Because.of the phase shift extraction in Eq. (5.19), the angles of the elements of the ma-trix P are smooth throughout the required frequency range (Figure 5.3), and the additional time delay can be accurately determined. The approximation for the propagation matrix A is found as &fitted = e-JW?fitted- (5.21) Note that each element of the synthesized matrix A has a single time delay. This time delay is the sum of the delay x of the speed of light and the delay from the fitted matrix P. The multiplication by -1.0 of the off-diagonal terms of the matrix P affects very much the synthesized functions. Without performing this, the phase angles of these elements would be 180 degrees from the minimum-phase behaviour, as shown in Figure 5.4a. With this, it is extremely difficult to find the correct time delays for the fitted functions. There may be a physical meaning which justifies the multiplication of the off-diagonal elements with the factor -1.0: by observing the propagation functions in the time domain, as shown in Figure 5.5, one can see that when one phase is energized with an impulse, the other phases have negative surges through coupling. As explained in section 5.3c, the impulse 68 waveforms arriving at the other end of the line correspond to the elements of a column of the propagation matrix a(t). If the matrix A were to be calculated directly from the triple product in Eq. (5.12), with-out first extracting the phase shift e"j(0T, its elements would show highly oscillatory angles (between +180 and -180 degrees) in the high frequency region, as illustrated in Figure 5.4b. This would make it nearly impossible to relate the fitted function to a single time delay. For example, if A,,(co) were to be synthesized with FIT-S, a minimum-phase function, say FH(co), would be obtained. The two functions are related to each other as ^ii(ro) = e-- /' f f lwFn(co), 9—0 where imax = co ^' a n c * w n e r e a n d @F a r e the phase angles of A,,(co) and F,,(co), respectively. xm a x is the.time delay calculated at the fastest frequency comax, within the considered frequency range. In order to evaluate this time delay with acceptable engineering accuracy, the phase angle of the A,, element has to be smooth and continuous throughout the entire frequency range. This condition would not be met if the phase angles were obtained directly from Eq. (5.12), as shown in Figure 5.4b. 69 Figure 5.2 - Elements of column 1 of matrix Y c for a double-circuit line 70 1 0 1 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 1 0 1 0 Frequency (Hz) Figure 5.3 - Elements of column 1 of matrix P 71 2 0 0 1 0 0 0 1 0 0 2 0 0 L — z 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 5.4a - Phase angle (degrees) of elements P 2 , and P 3 , before mult iplying with -1.0 2 0 0 1 0 0k 0 - 1 0 01 2 0 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 5.4b - Phase angle (degrees) of elements A , , and A 2 , from E q . (5.12) 72 5.3c Inverse Fourier transform of the propagation matrix A The accuracy of the synthesized functions can be further checked by comparing their in-verse Fourier counterparts in the time domain with those of the exact functions. Figure 5.5 shows the exact and synthesized functions of the first three elements of column 1 of a(t) (time-domain form of A(co)) for the line of Figure 4.2; the fitting accuracy is very high. The program of [55], with the Lanczos filter described in section 3.3a applied to the exact and synthesized functions, were used to obtain these inverses. Analogously to the suggestion in [36, p. 4-95], the physical meaning of a(t) can be ex-plained (Figure 5.6) by connecting a set of currents I s o u r c e to node m, together with a shunt characteristic admittance matrix Yc(co), and shorting all phases at node k. In Figure 5.6, summing the currents at node m gives Isource = Yc(co)Vm + I m . x 1 0 5 - [ Exact function -\ Synthesized function -L ^ a i i ( t ) \ \ ^ - - " a 2 1 ( t ) J s ^ T ^ ' a 3 i ( t ) 0.331 Time (msec) 0.34 Figure 5.5 - Time-domain propagation functions: a,,(t), -a 2 1(t), and-a 3 1( t) 73 I k I m > k m S h o r t e d s o u r c e Figure 5.6 - Current source with matching admittance Yc(co) connected to node m Then, from Eq. (5.10), with the direction of the currents in node k the same as that from node m, If I s o u r c e is 1.0 for phase #1 and zero for the other phases, for all frequencies, then Ik be-comes the first column of the matrix A. This situation corresponds, in the time domain, to the injection into phase #1 of a unit current impulse at node m, with the other phases open. Therefore, the inverse Fourier transform of A,, must be the current impulse that arrives at node k. The inverses of the other elements in the first column of the matrix A are the cur-rents of the other phases arriving at node k, induced through coupling between phases. It is seen in Figure 5.5 that the term a,,(t) has positive polarity and a relatively high magnitude, while a21(t) and a3!(t) have negative and smaller magnitudes (-a21(t) and -a31(t) are plotted). These functions exhibit the typical attenuation and distortion of impulse waveshapes. 5.4 Synthesis of Y c and A for asymmetrical transmission line configurations Figure 5.7 shows a railway running parallel to a 1 mile long three-phase 138 kV line. The rails are represented as conductors and the ground wire is assumed to have uniformly zero potential. Thus, the line configuration is a five-phase system. This line configuration Ik = A l source (5.22) 74 6 ft G.W.: 1-159 MCM ACSR _4-^>, Dia. = 0.576 in. Phase conductor: 1-1033 MCM A C S R b Q<^ A 15 ft i - i u j j mum ftivSK D o< -Dia. = 1.212 in. 7 ftf 1 ? f t 5 f t Rails #1 & #2: I I Dia. = 9.19 in. #1 #2 28.5 ft o o < 10 ft Q C 14 ft a V Length = 1 mile 40 ft 100 Ohm-m ground resistivity Figure 5.7 - Railroad near a three-phase 138 kV line is extremely asymmetrical. Figures 5.8 and 5.9 show the elements of column 1 of its char-acteristic admittance and propagation matrices, respectively. The agreement between the fitted and exact functions is very good, similar to the case of the double-circuit line of Figure 4.2. Table 5.1 lists the average number of poles that the fitting routine FIT-S typically allo-cated for the synthesized functions. Table 5.1- Average number ofpoles allocatedfor synthesized functions Function Diagonal elements Off-diagonal elements A 24 30 Y c 14 17 . 75 x 1 0 Elements order from top down: 1, 2, 3, 5 & 4 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 (D <U L a 0> w "5b c < a in <a O H 1 0 0 8 0 6 0 1 4 0 2 0 0 2 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 5.8 - Elements of column 1 of matrix Y c for an asymmetrical line 76 0 . 8 L 0 . 6 ^ -a 3 J2 0 . 4 0 . 2 L 0 1 0 - T Exact function Synthesized function D i a g o n a O f f - d i a g o n a l s Maximum error < 3.0 % 1 0 1 0 1 0 1 0 Frequency (Hz) 1 0 0 5 0 u 60 (D w C < 0 a - s o - 1 0 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 5.9 - Elements of column 1 of matrix A for an asymmetrical line 77 5.5 Open/Short -circuit tests in the frequency domain As described in section 4.5, open- and short- circuit calculations are performed to check the accuracy of the synthesized functions. The open-circuit voltages and short-circuit cur-rents are directly evaluated in phase co-ordinates here. Applying boundary conditions to Eq.'s (5.9) and (5.10), these voltages and currents at the receiving end m become V m - < ¥ « „ = Yc 1 [A + A - , ] " , V k , (5.23) I m - ^ / = 2 [ A - A - 1 ] " 1 Y c V k . (5.24). If a voltage column matrix E s is applied at the sending end of the line (node k), then V k = E s , o r V m - ^ ^ Y c ' t A + A - ' f ' E s , (5.25) I m - ^ = 2 [ A - A - 1 ] " 1 Y C E S . (5.26) 5.5a Zero sequence excitation The zero-sequence source voltage phasors in this test are real and all have 1.0 p.u magni-tude. Figures 5.10-5.13 show the responses obtained for both lines (Figures 4.2 and 5.7) in this test. The approximate and exact responses are in very good agreement. The responses of phases b and c (not plotted) for the line of figure 5.7 also agree very well. 78 5.5b. Unsymmetric excitation For the double-circuit line of figure 4.2, the unsymmetric source voltage phasors (E s) are the same as in section 4.5b, namely Esl = Es4 = \Z0°, Es2=Es5=0, andEs3=Es6 = 1Z120 0, whereas for the line of Figure 5.7, the following values were used: Esa =.\Z0°, Esb = lZ- 120°, Esc = 1Z1200, andEs-raim = Es-raim = 1Z0 . Figures 5.14 - 5.17 show the responses obtained for this test. The agreement between the exact and synthesized functions is very good. 79 2 5 2 0 . 1 5 . 1 0 . 5 . 1 0 1 0 1 2 Exact function Synthesized function ^ L"2 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 1 Frequency (Hz) 1 0 1 0 Frequency (Hz) Figure 5.10 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 4.2 80 o3 w 0) in ro °a CM C/) <D in ro C D 08 CO in <D w ro 0 . 6 . 0 . 4 0 . 2 . 1 0 0 . 6 0 . 4 . 0 . 2 1 0 0 . 6 0 . 4 0 . 2 Frequency (Hz) 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 Exact function Synthesized function Figure 5.11 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 4.2 2 5 ro cu to ro 2 0 1 5 . Lx 1 0 - 2 1 0 1 0 ' r o CM ' r o a: 1 0 8 -6 -1 0 Exact function Synthesized function 1 o Aim, 1 o 1 0 Frequency (Hz) Frequency (Hz) 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Figure 5.12 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 82 1 2 a) ro 1 0 1 0 Frequency (Hz) 'ro 0C 2 0 1 5 . 1 0 1 0 2 0 1 5 -1 0 1 0 1 0 Frequency (Hz) C M 2 1 o ro 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 Exact function Synthesized function Figure 5.13 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 5.7 83 CO 0) CO 03 LO o 3 CM CO cu CO co CO o3 oo CO CO CO CD 2 5 2 0 1 5 1 0 5 -- 2 1 0 1 0 4 . 1 0 1 2 1 0 . Exact function Synthesized function 1 0 1 0 Frequency (Hz) 1 0 1 0 Frequency (Hz) Frequency (Hz) Figure 5.14 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 4.2 84 0 . 8 0 . 6 08 co 0 4 CO CD ° 3 CM CO CO CO CD x: Q. 0 . 2 0 I 1 0 0 . 2 0 . 1L 1 0 1 0 Exact function Synthesized function 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 1 0 Frequency (Hz) co o3 co CO CO CO CD sz CL 0 . 8 0 . 6t 0 . 4L 0 . 2L 1 0 1 0 1 0 Frequency (Hz) Figure 5.15 - Receiving end short-circuit currents (p.u. magnitude) for line of Figure 4.2 85 ro 1 4 1 2 .1 0 8 cu CO c ro o. Exact function Synthesized function 1 o 1 0 CN ro or 1 0 J 1 0 1 0 1 0 Frequency (Hz) 1 0 . 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 Figure 5.16 - Receiving end open-circuit voltages (p.u. magnitude) for line of Figure 5.7 86 1 2 1 0 |-8 CD 8! 6 l-CD ^ 4 2 I-1 0 2 0 1 0 2 0 1 0 Exact function Synthesized function 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 1 0 Frequency (Hz) 1 0 1 0 . 1 0 Frequency (Hz) 1 0 1 0 1 0 Figure 5.17 - Receiving end short-circuit currents (p.u. magnitude) for line o f Figure 5.7 87 5.6 Time-domain simulation tests To demonstrate the accuracy of the proposed frequency-dependent line model, time-domain simulations for transient as well as steady-state conditions were performed. The re-sults are presented in this section. For testing purposes, the equations of the developed fre-quency-dependent transmission line model were interfaced with the UBC EMTP version. MicroTran through subroutine CONNEC. This interface is described in Appendix C. Addi-tional tests and field measurement comparisons for typical lightning events will be shown in Chapter 6. 5.4a Open-circuit test with zero-sequence energization In this test, all phases of the double-circuit line of Figure 4.2 are energized simultaneous-ly with a 1.0 p.u. step voltage at time zero, with the receiving end of the line open. Figure 5.18 shows voltages at the receiving end of phase 1. Note that the response produced 0 5 1 0 1 5 Time (msec) Figure 5.18 - Receiving end phase 1 open-circuit voltage (p.u.) 88 by the new line model agrees very well with those simulated with the frequency-domain model [41], and with the existing EMTP frequency-dependent line model [34]. It is worth emphasizing that the EMTP frequency-dependent line model using a constant transformation matrix at 5 kHz also gives a very reasonable response in this test, because the zero sequence modal transformation vector has very little, if any, frequency dependence. 5.4b Magnetically induced voltages from a three-phase 138 kV line Figure 5.19 shows a railway running parallel to a three-phase 138 kV line. As men-tioned earlier, the rails were modelled as conductors. In order to investigate the possible magnetically induced voltages on the rails and to show the performance of the new line model,- a test condition is set up as illustrated in Figure 5.19. In this test, all three phases are energized simultaneously at time zero with a symmetric 60 Hz voltage source through a 190 Q source impedance. All phases, including the two rails, are shorted at the receiving end. Magnetically induced voltages will appear at the open end of the rails. Figure 5.20 shows the induced transient and steady-state voltages of rail #1. Table 5.2 lists the RMS va-lues obtained with the compared models, as well as the exact solution obtained from a pha-sor solution. The results in Figure 5.20 and Table 5.2 show that the existing EMTP frequency- dependent model, with constant transformation matrices evaluated at two differ-ent frequencies, produces reasonable answers during the transient period, however, it leads to inaccurate results in the steady-state condition. The steady-state induced voltages pro-duced by the new line model agree much closer with the exact solution. Note that even though the transformation matrix was calculated at 60 Hz, the existing frequency-dependent 89 t = 0 G.W. L e n g t h = 1 m i l e 190 Ohms ^ J P w ) f \ j ) i 3 8 k v M Rail #2 i Rail #1 Figure 5.19 - Magnetically induced voltages on rails CU s o > c CD '</> c CD 1 - New model 2 - Existing model w. constant Ti at 5 kHz 3 - Existing model w. constant Ti at 60 Hz 0 . 0 2 0 . 0 4 0 . 0 6 Time (msec) 0 8 1 0 1 5 2 0 2 5 3 0 Time (msec) 3 5 4 0 Figure 5.20 - Magnetically induced voltage on rail #1 90 line model, which uses a constant transformation matrix, still provided wrong steady-state solution at this frequency. The reason is because this line model used only the real part of the transformation matrix in its solution algorithm. This assumption is not required by the new line model. Table 5.2 - Steady-state induced voltages (RMS) Line Model Rail #1 Rail #2 Exact Solution 16.41 V 17.49 V New Model 16.27 V 18.19 V Frequency-dependent Model w. constant T, at 5 kHz 20.75 V 23.03 V Frequency-dependent Model w. constant T, at 60 Hz 11.32 V 15.98 V 91 CHAPTER 6 Simulation of Lightning Surges and Field Test Comparisons 6.1 Overview Unexplained system equipment failures often occur during heavy thunderstorms, for which power system engineers try to find out the causes with lightning surge studies. Most engineers use computer programs today such as the EMTP for such studies. Two new mo-dels have been developed for such programs, to represent transmission lines and trans-mission towers. This Chapter covers the application of these two models for the simulation of four typical lightning surge types, as outlined in Chapter 1. Sections 6.2 and 6.3 empha-size the improvements accomplished with these models, by comparing the simulation results with field tests and with field recordings. 6.2 Tower/ground wire surges When lightning strikes the tower or the ground wires, overvoltages are developed across the insulator strings. If these overvoltages are above the voltage which the insulators can withstand, flashovers to the phase conductors will occur. Simulations of these scenarios have indicated that the line insulator voltages are influenced by the models used in repre-senting the transmission towers. In the past, the transmission tower has often been repre-sented as a transmission line with uniform surge impedance. As pointed out in section 2.3, 92 the transmission tower surge impedance does vary as the wave travels from the top to the base. Thus, a nonuniform transmission line model is required. The simulations of the measured insulator voltages from the published experiments of Ishii et al [18] are analyzed in this section. They performed measurements on 500 kV double-circuit towers. Their experimental set-up is similar to the case where lightning strikes the transmission towers, except that the current source was low in magnitude. The circuit parameters for EMTP simulations are displayed in Figure 6.1. The circuit consists of three sections of a double-circuit line with the geometry of Figure 6.2, of three transmission towers, and of a low current source with a 400 Q source impedance. A 300 m cable with a surge impedance of 50 Q is connected between the source and the tower top. All of the con-ductors are terminated with 350 Q resistors at both ends. The rest of the parameters used are similar to those described in [18], with the exception of the line model: the authors of [18] used the Semlyen frequency-dependent line model, while a constant-parameter line model is used here. This difference has no noticeable effect on the insulator voltages in this case. The middle tower No. 7 is modelled as 4 transmission line sections (Figure 6.3), three uni-form and one exponential -line, with the crossarms ignored. In concurrence with the measured impedance values in [18], the uniform sections have a uniform 220 Q surge impe-dance, while the exponential one has a nonuniform surge impedance of 220 Q at the upper end and 150 Q at the bottom end. Capacitances of the insulator strings and stray capaci-tances from the phase conductors to the tower are also modelled. The simulated and measured waveforms of the three crossarm insulator voltages (potential difference between crossarm and phase conductor) at tower No. 7 are shown in Figure 6.4. A time delay of 93 approximately 1 ps is seen in these voltages because of the connected cable between the cur-rent source and the tower top. The voltages are very close compared to the measurements. Thus, the simulation results supports the validity of representing transmission towers with the exponential transmission line model! It should be mentioned that the capacitances (from the insulator strings and from the phase conductor to the tower) influence and retard the electromagnetic field propagation in the transverse direction. This aspect is not taken into account by the conventional time-domain transmission line models because they only represent the T E M mode in the horizon-tal direction [56]. Without including these capacitances, the lower-phase insulator voltage would have the highest magnitude. This would be the opposite of what had been measured. 400m Constant parameter line < . . ' Uniform line model 17 Q 449 m Constant parameter line Uniform line model 17 n No. 8 Exponential line model Figure 6.1 - Circuit for E M T P simulations 13-2 m (No. 56) If 11.1 m(No. 6-7) 133m(No.7«) y 127m A 127m 40X)m(No.56) 26X)m(No.6-7) 36.0 m (No. 7-8) ^ 20.8 m ^ 19.0 m 1 o 2 O 3 6 19.0 m O 4 -5* 20.8 m Ground wire, 1.75 cm diameter Phase conductor 4 X 3.84 cm diameter at 50 cm spacing o—o ! o o ^ i 50 cm 500 Ohm-m ground resistivity O 5 O 6 94 Figure 6.2 - Conductor geometry at average height 1 CR#1 CR#2 CR#3 L i n e #1 Line #2 Line #3 Line # 4 Uniform Uniform Uniform Nonuniform (Exponential line) / / / / / / / / / / / ~~" "~ " Figure 6.3 - Tower representation 95 Measurement Simulation 1 2 0 40 V / d i v . 50 « (89%) Time (us) 1 2 0 > 0 1 2 3 4 5 Time (us) Figure 6.4 - Simulation of Ishii et al experiment and measurement from [18] 96 6.3 Simulation of a recorded shielding failure event in a 765 kV system A simplified diagram of the American Electric Power (AEP) 765 kV system surrounding the Rockport area is shown in Figure 6.5. A monitoring system has been installed at Rock-port since 1986 [2]. This system has captured many interesting shielding failure events dur-ing heavy thunderstorms. Shown in Figure 6.7 are the recorded waveforms for one such event. The expanded waveforms around the strike location are shown in Figure 6.8. A care-ful examination of these waveforms revealed that a shielding failure had occurred on phase A of the Rockport-Sullivan line, approximately 31 miles away from Rockport. The travel-ling surge reflections denoted by a-d in Figure 6.8 are identified on the Rockport phase A voltage as: • (a) - Initial lightning surge arrives at Rockport, which also produces coupled surges on phases B and C. • (b) - First reflection from Sullivan comes. The travel time between this and the first surge of (a) suggests that the stroke hit the Rockport-Sullivan line, approximately 31 miles from Rockport. • (c) - Second reflection from Sullivan arrives. • (d) - First reflection from Jefferson arrives. In the simulation of this shielding failure event, all of the transformers are ignored, ex-cept for the one at Jefferson station, The transmission lines, system linear reactors, and sys-tem 60 Hz sources behind their appropriate short circuit impedances are represented. The struck line (Rockport-Sullivan) was first represented with the new transmission line model and then with the existing frequency-dependent model. Because of insignificant 97 contributions to the simulation results, all other secondary lines were modelled with untrans-posed constant-parameter lines, calculated at 5 kHz. With the exception of the Rockport-Jef-ferson line, these secondary lines do not have to be included. The reason for including them in this simulation is because the author wanted to investigate if there are any significant re-flections coming back to Rockport from far away stations. It was concluded from observing the results that no significant reflections arrived back to Rockport from these far away stations. According to calculations performed in Appendix D, the maximum current for lightning strokes which can escape the shielding system (the two ground wires) and terminate directly onto the phase conductors is 8.13 kA. The event was simulated by injecting a current source with this magnitude into phase A at the strike location, on the Rockport-Sullivan line. The injected current waveshape had the concave wavefront recommended in [7]. Preliminary re-sults showed that the struck phase voltage was noticeable lower in magnitude than that of the recording. Subsequently, the current source was increased to 11.2 kA. The other character-istics of the current source are listed in Figure 6.6. The equations for the current source are: i(t) = -{At90(t/t90) + Bt90(t/t90)n} for t<t90 (front), / ( ? ) = - { / ] e - ( ' - '9o ) / ' ._ / 2 g - (W9o) / '2 ) for t>t90 (tail), where the front part has At^ = 5324.225, BV = 4755.775, t^ = 31.034483 [is and n = 8.0163697, and the tail part has I, = 11571.533,12 = 1491.5327, t, = 112.98168 \xs and t2 = 0.99275381 |xs. With this current source, the simulated waveforms are shown in Figure 6.9. As can be seen from these plots, the waveforms are in very good agreement with the recordings of 98 Figure 6.8, with respect to the magnitudes, travelling times and polarities. It can also be ob-served that the initial surges coupled into phases B and C produced by the new transmission line model match the recordings closer than those given by the constant T, frequency-de-pendent line model. It should be noted that a current magnitude of 11.2 kA is as probable as the value of 8.13 kA calculated with the electrogeometric theory in Appendix D. This electrogeometric model is only an approximation, and has a margin of error. In addition, it has been reported that subsequent strokes in a single flash sometimes have higher magnitudes than the first stroke (with a factor of up to 200%) [44]. The event described in Figure 1.5 is evidence of this, where the third stroke actually caused a flashover, which most likely had a higher magnitude than the earlier strokes. 99 5 ^ " 36.5' A 1 8 " ^ K— 6 6 o o o o L 44.5' 36.5' .W. A A 34.75" X 4 4 5 ' > h c =57.5' Line geometry at average height hq = 99.5' D U M O N T S T R I K E L O C A T I O N AAAA R O C K P O R T H A N G I N G R O C K Figure 6.5 - Simplified diagram of 765 kV system surrounding Rockport - 6 0 0 0 - 8 0 0 0 1 0 0 0 0 1 2 0 0 0 ACrest "11-2 kA Tfront = 20 LIS T t a i l=100Lis \ Sm a x = -1.4 kA/us V Note that T f l o n , is defined through-\ the 30% and 90% points [7]. 2 0 4 0 6 0 Time (us) 8 0 1 0 Figure 6.6 - Lightning current source (A) 100 x 1 0 o O) s? o > •e o Q. J< O O o > m © ra o a o o a: 0) o> o > (J a> in re t o a O o a. 0 . ' 5 U 0 . 5 1 k 0 . 5 U 0 . 5k 1 V 0 . 5 L 0 : 5 Figure 6 . 7 - Lightning shielding failure waveforms recorded on Oct. 7 , 1990 102 103 6.4 Shielding-failure flashover If currents in the order of tens of kA or less escape the shielding and terminate directly onto the phase conductors, it is most unlikely that these shielding failures will flash over. This is not always the case, however. Shielding failures can result in flashovers, provided that there is a combination of worst scenarios. For example, considering the 765 kV system of Figure 6.5 again, if a subsequent stroke of approximately -16 kA in magnitude is injected at a negative peak of the system steady-state voltage, a voltage larger than the line standard BIL would appear across the line insulators. For a normal 4.4 m 765 kV line insulator, a value of 2600 kV critical flashover voltage is assumed for lightning strokes having a 20 ps wave front [5]. The lightning current source in this case has a similar waveshape with the one in Figure 6.6, with A c r e s t = -16 kA, T f r o n t = 20 ps, T t a i l = 100 ps, Sm a x= -2.0 kA/ps. A vol-tage-dependent switch is used between phase A and ground to model the flashover. This switch closes as soon as the voltage across it is greater or equal to the 2600 kV BIL. Figure 6.10 shows the simulated phase-to-ground voltages at Rockport when the current source was applied 2 miles away, on the Rockport-Sullivan line. Note that Phase A shows a voltage close to zero during the fault period (about half a cycle). The purpose of the foregoing analysis is to identify situations where shielding-failure fla-shovers could occur. It is not within the scope of the present thesis to investigate the details of the breakdown mechanism of the line insulation. There is a large volume of experimental work on this topic [5, 7, 44]. To assume that flashover occurs just when the voltage exceeds the breakdown volt-time characteristic is an approximation, which is useful for preliminary studies. More detailed line insulation co-ordination studies can then be done as outlined in 104 the above references. Experimental volt-time characteristics are only valid for standard lab-oratory impulses. For non-standard lightning waveshapes, as typically seen in field record-ings, these characteristics are of limited use in predicting the performance of lines. 6.6 Back flashover The circuit of Figure 6.11 is used for simulating back flashover events. It represents a 765 kV line of the geometry shown in Figure 6.5, and consists of 13 towers with each span being 0.25 miles. To avoid reflections coming back, the line is terminated with its own surge impedance at both ends. The untransposed constant-parameter line model, calculated at 500 kHz, is used for the line spans between the towers. The struck tower (tower #7) is modelled as an exponential transmission line, while the others are modelled as constant-pa-rameter lines. With 60 Hz voltages connected to both ends of the line, the lightning current source shown in Figure 6.12 is injected at the tower top of tower #7. Since phase A was as-sumed to be at its maximum 60 Hz positive voltage at time zero, its line insulator will see a higher voltage compared to other two phases. Thus, a back flashover will most likely occur to this phase first. Figure 6.13 shows the tower top, phase A, and phase A insulator voltages obtained at tower #7, for this incident. Note that the phase A insulator voltage reaches a maximum value of -1.7 M V . With the volt-time characteristics of [5], for a 4.4 m 765 kV insulator, a critical flashover voltage of 3.1 M V is obtained for standard lightning wave-forms having a 3 us wavefront. Therefore, for a flashover to develop, the injected current must be increased by a factor of (3.1/1.7). Thus, a back flashover could occur if a current of 45(3.1/1.7) or 82.6 kA magnitude, with the same front time, hits tower #7. The same back 105 flashover scenario can be established by increasing the tower footing resistance to a larger value, without having to change the magnitude of the lightning current. Figure 6.14 shows the respective voltages when the tower #7 footing resistance is increased to 180 Q. Thus, one can conclude that back flashovers could occur on E H V transmission systems provided that the lightning current magnitude is large enough, or that the tower being hit has a high footing impedance. It must be noted, however, that EHV towers are usually designed with low footing impedances, with the exception of those in rocky mountainous locations. For transmission lines over this type of terrain, simulations are important to predict the line out-age rates, and to analyze the influence of transients on protection systems. Summary Lightning studies for tower/ground wire surges and back flashover events are only con-cerned with a few tower spans close to the strike location, because the surges attenuate quite fast in the ground (zero sequence) mode as the distance increases. For shielding-failure and shielding-failure flashover events, however, the entire struck line has to be modelled because the waveforms may still be significant at a far distance down the line. For example, refer-ring to Figure 6.5, it has been observed that a reflection coming back from Sullivan rein-forced the oscillations which were initially created by the shielding failure surge; in the high-voltage transformer tertiary windings at Rockport [57]. Thus, even after travelling the length of the line, a reflected surge could still excite the transformer natural frequency os-cillations, and result in overvoltages in the windings. For this reason, it is best to represent all incoming lines to the substation with their entire length, in investigating the effects of 106 shielding failure and shielding-failure flashover waveforms impinging on substation equipment. For most practical purposes, the struck tower can be represented with a one-section ex-ponential line model, as shown in Figure 6.11, which gives good enough results and is easy to use in simulations. For the situation where the insulator voltages of different crossarms need to be known more accurately, the transmission tower must be broken up into and mo-delled as a number of exponential line sections. The number of sections depends on the number of crossarms which the tower has. For the situations of tower/ground wire surges and back flashover, the line representing each individual span between towers should include the ground wires as additional conduc-tors, and represent the towers along the path, whereas in shielding-failure and shielding-fai-lure flashover cases, the line can be modelled with the ground wires eliminated. The accuracy of a particular line model affects greatly the results in those events which involve the phase conductors directly, such as the shielding failures and the shielding-failure fla-shovers, whereas it has less influence in the tower/ground wire surges and back flashover cases. The reverse holds true for the transmission tower representation. 107 108 Exponential line Z c = 300 and 100 n at top and bottom, respectively. # Uni form line Z c = 3 0 0 f i Figure 6.11 - E M T P circuit for simulation of a back flashover event X 1 0 A c r e s t = -45.0kA Tfront = 3 U.S T,ail = 75 JXS Sm a x = -37.5 kA/ns Figure 6.12 - Lightning current source (A) 109 Figure 6.13 - Simulated voltages at tower #7 with a 15 Q footing resistance 110 x 1 0 6 0 . 5 I Time (us) Figure 6.14 - Simulated voltages at tower #7 with a 180 Q footing resistance I l l CHAPTER 7 Future Research in Transmission Line and Transmission Tower Modelling Although the main goals set for this thesis have been accomplished, the author believes that improvements can be made if the following additional research is carried out. 1) Extend the direct phase-domain modelling technique of overhead transmission lines to underground cable modelling. The method described in this thesis is general and could po-tentially be applicable for the case of underground cables. 2) Include corona effects in the simulation of lightning surges, if reliable corona models are available. This aspect has been ignored in this thesis due to the lack of field test verifica-tions. It is appropriate to note here that the upper part of the wavefront is most likely being reduced by corona, as the surges travel along the line. This wavefront reduction has been noticed in the recorded waveforms in reference [2]. This system is in a remote location, and has been experiencing frequent lightning incidents 50 miles or farther away. The recorded wavefronts have been somewhat slow (> 20 ps), compared to typical lightning data in the literature. The corona effect is thought to have a major influence on this. 3) Develop equations for the calculation of the surge impedances of the exponential line model, which can be used to represent typical high voltage transmission towers. The tower surge impedance is a function of its dimensions. When the conductor geometries are speci-fied in line constants programs, the tower dimensions are indirectly given as weil, and it may 112 therefore be possible to modify line constants programs to produce tower surge impedances as well. Built-in models in such programs should cover most of the typical towers presently used in industry. The common classes are described in [5, 15], which includes the lattice, lattice and steel pole, H-frame, guyed-Y and typical 735-800 kV towers. 4) Modify the EMTP to reduce storage requirement for cases involving long trans-mission lines and A C steady-state initial conditions. The derivation of this modification is shown in Appendix E. 113 C H A P T E R 8 Conclusions EMTP simulations have become widely accepted for lightning transient studies, which are done to prevent equipment failures and transmission line outages. In such studies, it is important to realistically simulate and understand all possible events created by lightning. From preliminary EMTP simulations for these events, it was seen that the representation of the transmission tower and transmission line influences greatly the simulation results. It was also observed that the existing models used in representing towers and lines have limitations. For this reason, new and improved models were developed in this thesis. An exponential tapered line model was derived for the tower. Its time-domain simula-tion results showed good agreement with those produced by a cascaded model where the line is divided into many uniform sections. The model is numerically stable because it syn-thesizes the line functions in the frequency domain with stable rational functions of minimum-phase-shift type. It includes these functions in the time domain by means of re-cursive convolutions. The model is suitable for the representation of transmission towers because the tower, from a travelling wave standpoint, is a nonuniform structure. It can also be used to model general nonuniform transmission lines with complicated space-dependent characteristic impedances, as described in [21]. A simple, fast and accurate method based on the Newton-Raphson iteration technique was implemented for finding transmission line modal parameters. It overcomes the . • • . • • 1 1 4 eigenvector and eigenvalue switching problems experienced by most other methods. It pro-duces modal parameters which can reasonably be fitted with rational functions of minimum-phase type. This should improve the modal transmission line formulation of [35]. A new wideband time-domain transmission line model was developed. It includes the complete frequency-dependent nature of untransposed overhead transmission lines directly in phase co-ordinates. It. avoids the assumption of a constant transformation matrix, which the existing frequency-dependent transmission line model uses. Extensive simulation results were presented to show that the new model provides reasonably accurate answers for both steady-state and transient studies. The model is numerically stable and can be used for most types of transient studies. In terms of simulation time, the proposed model was, on average, about 25% slower than the existing frequency-dependent model using a constant T,. How-ever, this could be improved if the programming code is optimized. Typical lightning surges were simulated to compare the solutions produced by the devel-oped models with actual measurements. With these new modelling approaches, simulations matched field recorded transient waveforms closer than existing models. 115 APPENDIX A Asymptotic Behaviour of Modal Transformation Matrices at Very Low and High Frequencies It is shown in this appendix that the elements of the matrix T, of Eq. (4.3) become real and approach the asymptotic values at the high and low frequency limits. A. l Behaviour at very low frequencies When the frequency approaches zero, the conductor internal impedance matrix Z c be-comes the d.c. resistance matrix Rj,.. Then YZ=;coP-'[/-coLi 0s 0P + Z e + R dc], (Al) ignoring the shunt conductances of the line. Since the first two terms in the right-hand side of the above equation approach zero faster at very low frequencies, it follows that Y Z —> ycoP_ 1RdC as co - » 0. With P"1 is the capacitance matrix C, the Y Z product of Eq. (Al), at low frequencies, is Y Z « / c o C R d c . (A2) Because the eigenvectors stay unchanged if the matrix to be diagonalized is scaled with a constant, Eq. (A2) can be reduced to just the product of CRdc- If all of the conductors have the same dc resistance value, the matrix to be diagonalized is then reduced further to just the matrix C. Thus, the eigenvectors of the Y Z product at very low frequencies can be found by diagonalizing only the matrix C, and Since C is real and symmetric, its eigenvectors are real and constant. This agrees with the solution method described in [33]. 116 A.2 Behaviour at very high frequencies At very high frequencies, the elements of the two matrices Z c and Z E in Eq. (4.6) become pi j(o\ifdei 2nrd, el 2n V r )> ycop0 2n 2de2(yi +yj) Dl (A3) (A4) Eq. (A4) is obtained because \de2\ << Dy, for co -> oo, and In WA « l n ^D]j + Ade2iyi+yj) = ln 1 + Ade2{yi +yf) Dl 2de2(yi +yj) 4 because for small 0, ln(l + 0) —» 0. Then, Y Z = - w 2 u 0 s 0 U + 1 p - . B 2TIZ0 (A5) where U is the unit matrix, and the elements of the matrix B are: Bii = ^ + ^ a n d B i J \irde 2de2(yi +yj) D (A6) where u.r = 77- is the conductor relative permeability. Because de\ and de2 -> 0, as (a ->• 00, the elements of the matrix B of Eq. (A6) event-ually become negligible at very high frequencies. Thus, from Eq. (A5), Y Z - » - c o 2 p 0 s 0 U , asm -> 00 , (A7) Therefore, the eigenvectors of the Y Z product approach constant and real values as the frequency increases to infinity. 117 Eq. (A7) is also obtained in [36, p. 4-58] by making the two assumptions: 1. Conductor resistances and ground return resistances are ignored. 2. At very high frequencies, all currents flow on the surface of the conductors, and on the surface of the earth. This is because the complex depth of penetration constants d e l and d e 2 approach zero at high frequencies. In actual calculations, losses are small, and that the complex depth of penetration con-stants are not exactly zero at the 1 MHz limit. Therefore, the eigenvectors have very small angles (close to zero) at this limit. This has been seen in all the eigenvectors obtained for the lines tested in this thesis. 118 APPENDIX B Recursive Convolutions Bl . In matrix form If the elements of Y c are approximated by rational functions with real negative poles and zeros, then the numerical convolution of the yc(t) matrix and the voltages vector v(t) has the following form, y c ( 0 * v ( 0 = g(0 = y e q v ( 0 + M O , (Bi) where is the history current vector evaluated from g(t- At) and \(t- At), and y e q is a real, constant symmetric matrix. In a similar fashion, the elements of the A(co) matrix can also be synthesized with ra-tional functions, but with a time delay associated with each element. Numerical convolution of the corresponding time-domain function a(t) with a forcing vector f(t) is evaluated as a(0 * f(0 = r(7) = W ( 0 - (B2) The currents vector \pr0p(f) depends exclusively on past values of f(t) and r(t). Introduc-ing (Bl) and (B2) into (5.13) and (5.14) leads to (5.17) and (5.18), with = iprop(t) - hi (0 • (B3) B2. Integration constants and scalar convolutions When rational function syntheses are performed, each element of the matrix yc(t) be-comes a simple sum of exponential functions and contains a constant, namely 119 m yc-ij(t) = k0d(t) +1 kie-P*u(t). (B4) i=\ In a slightly different form, each element of the a(t) matrix contains a unique and differ-ent time delay Xy associated with it, and the constant k0 becomes zero in this case, namely m = k,e-P>^ (B5) i=\ Consider the matrix convolution of Eq. (Bl). This convolution operation is carried out as a normal matrix multiplication, but each multiplication involves a scalar convolution [30, 34, 35]. For example, the Convolution between element 11 of the yc(t) matrix and element 1 • of the v(t) matrix can be performed recursively as yc-i i (0 * v, (0 = g(t) = dv , (0 + hi(t), (B6) where m d = k0 + zZdi, (B7) /=i m hi(t) = T,cigi(t-At) + eivi(t-At), (B8) , i-\ ci5 di and e; are integration constants which depend on the integration method. For the trap-ezoidal rule of integration, they are defined as c, = ! ^ £ ; .*=., = - ^ £ - . - (B9) 2+piAt 2+pjA It should be noted that the stability of the trapezoidal rule of integration was proven and has been used in most of the present EMTP models [36, 34, 35]. The scalar convolution of an element of the matrix a(t) and a forcing function from equa-tion (B2) is also evaluated recursively, but based solely on past values, m aij(t)*fx(i) = r(t) = zZri(i), (BIO), 120 where n(f) = cir,(t- At) + dfi (t- Tij) + e(f\{t- Ty - At). (B11) The integration constants for Eq. (Bl 1) are the same as those from Eq. (B9). 121 APPENDIX C Interfacing User-defined Models with the EMTP version MicroTran The equations of the developed models in this thesis have been interfaced with the EMTP version MicroTran, in the time domain, through the subroutine CONNEC. This is a convenient way to test new models. The final implementation will be done directly in the main time step loop of the program, in a much more efficient manner. Consider a transmission line connecting node k to node m in Figure CI. MicroTran can describe the behaviour of the network seen by both terminals of the line, in the form of a Thevenin equivalent circuit, with the equation Vk Vko R THEV-k 0 hm vm Vmo 0 RTHEV-M Imk (CI) The Thevenin equivalent resistance matrix has zero off-diagonal elements if the rest of the network has no lumped-element connections between nodes k and m, which was true in all cases studied in this thesis. k ex km THEV-k m 'mk ^ T H E V - m 9 m m Figure C I - Equivalent circuit seen by the transmission line terminals 122 At each time step, the subroutine CONNEC supplies values for the open-circuit voltages and the Thevenin equivalent resistance network. These parameters, together with the equa-tions of the line are used to solve for the line terminals' voltages and currents. The equations are expressed as follows: C.l Single-phase Exponential transmission line equations From Eq.'s (3.21) and (3.11) derived in section 3.2b of Chapter 3, (C2) ^mep vm <^h-m Substituting Eq. (C2) into Eq. (CI), ikm Zkeq 0 -1 Vk - eh-k imk 0 Z mep Vm @ h-m . vk Vko R THEV-k 0 Zkeq 0 -1 Vk - eh-k vm Vmo 0 RmEV-m 0 Zmep Vm eh—m Hence, 1 0 0 1 RmEV-k 0 0 RmEV-m Zkeq 0 0 z mep Vk Vm Vko + Vmo R THEV-k 0 0 RmEV-m Zkeq 0 0 z, mep eh-k e/i-m (C3) The only unknown variables in Eq. (C3) above are the voltages vk and vm. Thus, these voltages can be solved and obtained at each time step. The currents can then be found from Eq. (C2). C.2 Multiphase frequency-dependent transmission line equations For a multiphase transmission line system as in Figure 4.1, Eq. (CI) contains matrix quantities instead, v* 0 RTHEV-IC 0 v m 0 R'77/£K-»i From Eq.'s (5.17) and (5.18) derived in section 5.2 of Chapter 5, yeq 0 - l 1km + lhist-k . V w -0 y e ? Imk + ^ hist-m Substituting Eq. (C5) into Eq. (C4), 123 (C4) (C5) yeq 0 0 y e q Hence, Urn imk + y.eq 0 0 Yeq lhist-k lfiist-m Vko Vmo R-THEV-k 0 0 ^THEV-m Hm imk J yeq 0 -1 _L RmEV-k 0 0 yeq _ 0 Rj-ffiF-m Vyt0 -1 imk 0 yeq _ lhist-m (C6) Thus, the currents i k m and i m k are solved first with Eq. (C6), the voltages can then be found from Eq. (C5). 124 APPENDIX D Estimation of Maximum Lightning Shielding Failure Stroke on a Typical 765 kV Line Based on Electrogeometric Theory G.W. o\o Conductor Bundle ~5~oo A hg = 129.5' - (40')2/3 = 99.5' = 30.33 m hc = 99.5'- (50')2/3 = 57.5'= 17.53 m a = 2.44 m d = 1 3 . 0 3 m a=l0.£> Figure Dl - Shielding angle of the right conductor bundle Figure D l shows the location of a ground wire and the closest phase conductor bundle at their average heights for a typical three-phase 765 kV flat transmission line. The average height is equal to the height at the tower, minus two thirds of the sag. In this conductor ar-rangement, shielding failures can only occur with the two outer phases. According to refer-ence [5], the following formulas are used, s=\or°-65, [3 = 0.8, where S is the striking distance and [3 is the correction factor. The meanings of these para-meters are illustrated in Figure D2. It is observed that, based on the striking distance theory, 125 the conductor bundle is vulnerable when lightning enters the uncovered area X s . The strik-ing distance S is determined by the lightning stroke magnitude I. Thus, if h g, h c and a are held fixed, the uncovered distance X s would be diminished when the current I increases until reaching its maximum value. With this condition, Figure D2 is modified to D3. From Figure D3, PSmax = 0.8S , m a x = hc +S , m a xsin9 , -> sinG = 0.8 17.53 Also sin 0 = sin 1 0 0 . 8 ° - cos"1 ( ^ - ] Expanding by using the trigonometric identity, sinO = sin 100.8° • cos _ / 6.5 ^ c o s \s~J v " J m a x y - cos 100.8° •sin v " - > m a x ' -> sine 6.38 -cos 1 0 0 . 8 ° - s i n cos _ , f 6.5 "l — ; v > - > m a x ' Combining Eq.'s (Dl) and (D2), 6 7 ^ - c o s l 0 0 . 8 ° s i n < J m a x Rearranging, cos 100.8°s inicos Dividing Eq. (D3) by coslOO.8' cos _i( 6.5 ^ i - > m a x ' = 0.8- 17.53 _ / 6.5 23.91 0.8 = 0. —> sin cos if 6.5 } + ^ ^ - 4 . 2 7 = 0. Solving for S„ (Dl) (D2) (D3) Smax = 38.83w. 126 Since S m a x = 1 0 » / ° 6 5 -> /= 0.029 • SlM, Therefore, 7 m a x = 0.029(38.85)154 = 8.13 It is observed from the above procedure that the maximum lightning current which can penetrate this particular line shielding is 8.13 kA. This value is comparable to that obtained from the modified equations suggested in [44], which were also derived from the electrogeo-metric model of [5]. The concept of attractive radius was included in these modified equa-tions. For helpful purposes, they are carried out here. The equation used for the current is still the same as that listed above, which is 7 m a x = 0.029^ 5 4 , but with S calculated as f S=Ya -BS - JB2 +ASCS V where: n = ^ =23.93, Xc-xg 2.44 n , o n A m = h^r 30.33- 17.53 = 0 - 1 9 0 6 ' p = 0.36 +0.168 ln(43-F 0 ) = 0.85 , As = m2-m2p2-$2= -0.7213, CS = (m2 + \)= 1.036, BS = pC, = 0.886. Hence 5=35.84 -> /max = 0.029(35.84)154 = 7.18 kA. Striking distance Figure D2 - Incomplete Shielding, Width X s is uncovered Maximum striking distance GW Conductor Bundle a = 79.2° 6 = 1 8 0 ° - a - c o = 100.8° - cos-1(z/2Smax) ps, max SlL. Figure D3 - Effective shielding, Unprotected width X s is reduced to zero 128 APPENDIX E Representation of a Long Line with Linear AC Steady-State Initial Conditions When the maximum simulation time t ^ is less than the travel time x of a long line, the line will' only be represented by its surge impedance matrix at the two ends if the case starts from zero initial conditions in the MicroTran version of the EMTP [36]. Computer time and storage space are then significantly saved because there is no need for recording the history terms for such a long line. This simplified model is most often used in lightning studies. However, one often wants to include the system steady-state voltages in lightning simula-tions as well. In that case, the simulation has to start from linear A C steady-state conditions. This can be accommodated by a simple modification, by adding A C steady-state history cur-rent sources at both ends of the line. This simplified modification can be illustrated as follows: For simplicity, consider a single-phase transmission line in Figure E l . Its current and voltage relationships at node k can be expressed as ikm(J) = v*(0 + histk, (El) Z km mk ) m Figure El - Single transmission line 129 where histk = - V m ( f ~ i ) -imk(t-x). (E2) Similarly, the same expressions can be applied for node m, but with k and m subscripts interchanged. With zero initial conditions, the line can then be replaced by shunt resistances with a value of Z at both ends, as shown in Figure E2. Figure E2 - Equivalent circuit o f a long line without initial conditions For linear A C steady-state initial conditions, additional current sources are needed, thus this equivalent circuit is modified to the circuit of Figure E3, where histk is the A C current source of Eq. (E2). The same expression with subscripts k and m interchanged is applied for hist™. m Figure E3 - Equivalent circuit o f a long line with linear A C initial conditions 130 R E F E R E N C E S [I] H.W. Dommel, "Digital Computer Solution of Electromagnetic Transients in Single-and Multiphase Networks," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-88, No. 4, pp. 388-399, 1969. [2] J.M. Schneider, E.N. Fromholtz, D.K Nichols, B.J. Ware, "The Rockport Transient Voltage Monitoring System," CIGRE Report 23-04, 1988 Session. [3] A. Greenwood, Electrical Transients in Power Systems, 2nd edn., Chapter 10, Wiley (1991). [4] K. Ragaller, Current Interruption in High-Voltage Network, Plenum Press, New York (1978). [5] J.D. Anderson, "Lightning Performance of Transmission Lines," Chapter 12 of Trans-mission Line Reference Book-345 kV and Above, 2nd edn., EPRI (1982). [6] CIGRE WG. 33.02, "Guidelines for Representation of Network Elements When Calcu-lating Transients," CIGRE Technical Brochure, 1990. [7] CIGRE WG 33.01, " Guide to Procedures for Estimating The Lightning Performance of Transmission Lines," CIGRE Technical Brochure 63, 1991. [8] P.P. Barker, et al, "Characteristics of Lighting Surges Measured at Metal Oxide Dis-tribution Arresters," IEEE Trans, on Power Delivery, Vol. 8, No. 1, pp. 301-309, 1993. [9] M.D. Perkins, et al, "Summary of Transient Data Obtained from Long-Term Monitoring of a 138-kV and a 500 kV Transmission Line," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-103, No. 8, pp. 2290-2298, 1984. [10] K. Nakamura, et al, "Artificially Triggered Lightning Experiments to an E H V Trans-mission Line," IEEE Trans, on Power Delivery, Vol. 6, No. 3, pp. 1311-1318, 1991. [II] C F . Wagner and A.R. Hileman, "A new Approach to The Calculation of The Lightning Performance of Transmission Lines III-A Simplified Method: Stroke to The Tower," AIEE Trans, on Power Apparatus and Systems, Vol. 79, pp. 589-603, 1960. [12] M . Darveniza, et al, "Modelling for Lightning Performance calculations," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-98, No. 6, pp. 1900-1908, 1979. [13] C F . Wagner and A.R. Hileman, "A new Approach to The Calculation of The Lightning Performance of Transmission Lines III-A Simplified Method: Stroke to The Tower," AIEE Trans, on Power Apparatus and Systems, Vol. 79, pp. 589-603, 1960. 131 [14] M . Kawai, "Studies of The Surge Response of a Transmission Line Tower," IEEE Trans, on Power Apparatus and Systems, Vol. 83, pp.30-34, 1964. [15] P.G Buchan and W.A. Chisholm, "Surge Impedance of HV and E H V Transmission Towers," Final Report for Canadian Electrical Association Contract #78-71, 1980. [16] M.A. Sargent and M . Darveniza, "Tower Surge Impedance," IEEE Trans, on Power Ap-paratus and Systems, Vol. PAS-88, No. 5, pp. 680-687, 1969. [17] W.A. Chisholm et al, "Lightning Surge Response of Transmission Towers," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-102, No. 9, pp. 3232-3242, 1983. [18] M . Ishii et al, "Multistory Transmission Tower for Lightning Surge Analysis," IEEE Trans, on Power Delivery, Vol. 6, No. 3, pp. 1327-1335, 1991. [19] L.O. Barthold, "An Approximate Transient Solution of the Tapered Transmission Line," AIEE Trans., PAS 77, pp. 1556-1561, 1958. [20] C. Menemenlis and Z.T. Chun, "Wave Propagation on Nonuniform Lines," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-101, No. 4, pp. 833-839, 1982. [21] M . M . Saied, A.S. Alfuhaid, M.E. ElShandwily, "S-Domain Analysis of Electromag-netic Transients on Nonuniform Lines," IEEE Trans, on Power Delivery, Vol. 5, No. 4, pp. 2072-2083, November 1990. [22] H.V. Nguyen, H.W. Dommel, J.R. Marti, "Tower Models for Lightning Surge Simula-tions," IEEE Paper No. 045-5 PWRD, Winter Power Meeting, New York, New York, Feb. 1994. [23] E.A. Oufi, A.S. Alfuhaid, M . M . Saied, " Transient Analysis of Lossless Single-Phase Nonuniform Transmission Lines," IEEE Trans, on Power Delivery, Vol. 9, No. 3, pp. 1694-1701, July 1994. [24] M.T. Correia de Barros, M.E. Almeida, "Computation of Electromagnetic Transients on Nonuniform Transmission Lines," IEEE Paper No. 397-0 PWRD, Summer Power Meeting, Portland, Oregon, July 1995. [25] W.A. Chisholm and Y.L . Chow, "Travel Time of Transmission Towers," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-104, No. 10, pp. 2922-2928, 1985. [26] Working Group on Lightning Performance of Transmission Lines, "A Simplified Method for Estimating Lightning Performance of Transmission Lines," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-104, No. 4, pp. 919-927,1985. [27] J.T. Whitehead, "Lightning Performance of TVA's 161-kV and 500-kV Transmission Lines," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-102, No. 3, pp. 752-768, 1983 132 [28] W.A ChisKolm, "Lightning Surge Response of Transmission Towers and Lines," Ph.D. Thesis, University of Waterloo, Department of Electrical Engineering, 1983. [29] R.S. Elliott, An Introduction to Guided Waves and Microwave Circuits. Prentice-Hall, Inc., [13] M.A. Sargent and M . Darveniza, "Tower Surge Impedance," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-88, No. 5, pp. 680-687, 1969. [30] A. Semlyen and A. Dabuleanu, "Fast and Accurate Switching Transient Calculations on Transmission Lines with Ground Return using Recursive Convolutions," IEEE Transac-tions on Power Apparatus and Systems, Vol.PAS-94, pp. 561-571, March/April 1975. [31] H.W. Dommel, "Line Model in EMTP," EMTP Review, October 1988 & January 1989. [32], J.R. Marti and H.W. Dommel, "Line Models for Lightning Studies," CANADIAN ELECTRICAL ASSOCIATION, March 1989. [33] J.R. Marti, et al, "Approximate Transformation Matrices for Unbalanced Transmission Lines," 9th Power System Conference, Lisbon 1987. [34] J.R. Marti, "Accurate Modelling of Frequency-Dependent Transmission Lines in Elec-tromagnetic Transient Simulations," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-101, pp 147-157, 1982. [35] L. Marti, "Simulation of Transients in Underground Cables with Frequency-Dependent Modal Transformation Matrices," IEEE Trans, on Power Delivery, Vol. 3, No. 3, pp. 109-1110, 1988. [36] H. W. Dommel, EMTP THEORY BOOK, 2nd edn., Microtran Power System Analysis Corporation, Vancouver, B.C., Section 4, 1992. . [37] L. M . Wedepohl, H. V. Nguyen and G. W. Irwin, "Frequency-Dependent Transform-ation Matrices for Untransposed Transmission Lines using Newton-Raphson Method," Paper 95SM602-3 PWRS, IEEE PES Summer Meeting, Portland 1995. [38] G. Angelidis, A. Semlyen, "Direct Phase-Domain Calculation of Transmission Line Transients Using Two-Sided Recursions," IEEE Transactions on Power Delivery, Vol. 10, No. 2, pp. 941-949, April 1995. [39] T. Noda, N. Nagaoka and A. Ametani, "Phase Domain Modelling of Frequency-De-pendent Transmission Lines by Means of an A R M A Model," Paper 95WM245-1 PWRD, IEEE PES Winter Meeting, New York 1995. [40] F. Castellanos, J. R. Marti, "Phase-domain Multiphase Transmission Line Modelling," IPST'95 - International Conference on Power System Transients, Lisbon, September 3-7, 1995, pp. 17-22. 133 [41] L . M . Wedepohl, Frequency Domain Analysis of Wave Propagation in Multiconductor Transmission Systems, E E 552 Class Notes and Computer Program, Dept. of Elec. Engr., The University of British Columbia, 1993. [42] V. Larsen of EFI in Trondheim, Norway; personal communication with H.W. Dommel. [43] Personal communication with H.W. Dommel. [44] The Fast Front Transient Task Force of the IEEE Modelling and Analysis of System Transients Working Group, "Modelling Guidelines for Fast Front Transients," Drafted version undergoing review process, No. FFT-01/02/95-40. [45] L. Marti, "Simulation of Electromagnetic Transients in Underground Cables with Fre-quency-Dependent Modal Transformation Matrices," Ph.D. dissertation, The University of British Columbia, November 1986. [46] H.V. Nguyen, H.W. Dommel, J.R. Marti, "Modelling of Single-Phase Nonuniform Transmission Lines in Electromagnetic Transient Simulations," Accepted for presenta-tion, 1996 IEEE Summer Power Meeting, Denver, Colorado. [47] H.V. Nguyen, H.W. Dommel, J.R. Marti, "Direct Phase-Domain Modelling of Frequen-cy-Dependent Overhead Transmission Lines," Accepted for presentation, 1996 IEEE Summer Power Meeting, Denver, Colorado. [48] Personal communication with L . M . Wedepohl. [49] L . M . Wedepohl, "Transient Analysis of underground power-transmission systems System-model and wave-propagation characteristics," Proc. IEE, vol. 120, No. 2, pp. 253-260,1973. [50] J.R. Carson, "Wave propagation in overhead wires with ground return," Bell Syst. Tech. Journal, vol. 5, pp. 539-554, 1926. [51] A. Deri et al, "A Simplified Model for Homogeneous and Multi-Layer Earth Return," IEEE Trans, on Power Apparatus and Systems, Vol. PAS-100, No. 8, pp. 3686-3693, 1981. [52] C. Dubanton, "Calcul approche des parametres primaires et secondaraires d'une ligne de transport. Valeurs homopolaires" ("Simplified calculation of primary and secondary parameters of transmission line. Zero sequence values", in French), E.D.F. Bulletin de la Direction des Etudes et Recherches, Serie B, pp. 52-62, 1969. [53] L . M . Wedepohl, "Application of matrix methods to the solution of travelling-wave phe-nomena in polyphase systems," Proc. IEE, vol. 110, pp. 2200-2212, Dec. 1963. 134 [54] R.H. Galloway, W.B. Shorrocks and L . M . Wedepohl, "Calculation of electrical para-meters for short and long polyphase transmission lines," Proc. IEE, vol. I l l , pp. 2051-2059, Dec. 1964. [55] H. W. Dommel, Fourier Transformation Program, Microtran Power System Analysis Corporation, Vancouver, B.C. [56] H.V. Nguyen, Discussion of "Computation of Electromagnetic Transients on Nonuni-form Transmission Lines," Paper 95 SM 397-0-PWRD, IEEE PES Summer Meeting, Portland, 1995. [57] A.J.F. Keri, Y.I. Musa, and J.A. Halladay, "Insulation Coordination for Delta Con-nected Transformers," IEEE Trans, on Power Delivery, Vol. 9, No. 2, pp. 772-780, April 1994.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Simulation of lightning surges on transmission lines
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Simulation of lightning surges on transmission lines Nguyen, Huyen V. 1995
pdf
Page Metadata
Item Metadata
Title | Simulation of lightning surges on transmission lines |
Creator |
Nguyen, Huyen V. |
Date Issued | 1995 |
Description | Computer simulations with programs such as the EMTP are widely accepted for lightning transient studies, to design the system in such a way that equipment failures and transmission line outages are prevented. In these studies, an accurate representation of the transmission lines and towers is required. Two new and improved time-domain models for the representation of these components are presented in this thesis. They are the frequency- dependent and the exponential transmission line models. The former is based on synthesizing the line functions directly in phase co-ordinates, and is used to represent multiphase transmission lines, while the latter is a single-phase model, which also synthesizes the line functions, to be used for tower representations. These models are valid for the required frequency range of 1 Hz to 1 MHz. They can be used for lightning studies as well as for many other types of transient studies. Extensive time-domain simulations have been performed to compare the results with those produced by existing line models. The results indicate that the new models are more accurate than some of the existing models, and numerically stable. Simulations of typical lightning surges have also been performed to compare the solutions produced by the developed models with actual measurements. With these new models, simulations match field recorded transient waveforms very well. In addition, a simple, fast and accurate method based on the Newton-Raphson iteration technique is introduced for finding transmission line modal parameters. It overcomes the eigenvector and eigenvalue switching problems experienced by most other methods. With a proper constraint equation for scaling the eigenvectors, this proposed method produces modal parameters which can reasonably be fitted with rational functions of minimum-phase type. With this, the modal transmission line modelling can be improved. The routine is very stable and can be used to diagonalize the Y Z product for uncommon lines with strong asymmetry. It has also been used to obtain the phase-domain parameters for the new frequency- dependent line model. |
Extent | 4968449 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-02-19 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0065121 |
URI | http://hdl.handle.net/2429/4816 |
Degree |
Doctor of Philosophy - PhD |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 1996-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_1996-091414.pdf [ 4.74MB ]
- Metadata
- JSON: 831-1.0065121.json
- JSON-LD: 831-1.0065121-ld.json
- RDF/XML (Pretty): 831-1.0065121-rdf.xml
- RDF/JSON: 831-1.0065121-rdf.json
- Turtle: 831-1.0065121-turtle.txt
- N-Triples: 831-1.0065121-rdf-ntriples.txt
- Original Record: 831-1.0065121-source.json
- Full Text
- 831-1.0065121-fulltext.txt
- Citation
- 831-1.0065121.ris
Full Text
Cite
Citation Scheme:
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>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0065121/manifest