Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Calculation of frequency-dependent parameters of power cables with digital images and partial subconductors Rivas, Richard A. 2001

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

Item Metadata

Download

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

Full Text

Calcu la t ion of Frequency-Dependent Parameters of Power Cables w i th D ig i ta l Images and Pa r t i a l Subconductors by R I C H A R D A. RIVAS Elec. Eng., Universidad Simon Bolivar, Caracas, Venezuela, 1988 M.Sc.,Univci'sidad Simon Bolivar, Caracas, Venezuela, 1993 A THESIS S U B M I T T E D IN P A R T I A L F U L F I L L M E N T O F T H E R E Q U I R E M E N T S FOR T H E D E G R E E O F D O C T O R OF P H I L O S O P H Y in T H E F A C U L T Y O F G R A D U A T E STUDIES (Department of Electrical & Computer Engineering) We accept this thesis as conforming to the required standard T H E U N I V E R S I T Y OF BRITISH C O L U M B I A July 2001 © Richard A. Rivas, 2001 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 £ l £ c - T R i CA>— g Co Me UT&fZ £ M 0 ' r V E B R J W £ The University of British Columbia Vancouver, Canada -f-Date AUGUST 1 ^ - Z.OO \ DE-6 (2/88) Abst rac t This thesis work presents the development of a general-purpose cable-parameter algorithm which uses digital images to discretize the cable geometry and the partial subconductor equivalent circuit method to estimate the cable parameters. The idea of the proposed geometry discretization technique is to draw the cable geom-etry, or scan its photograph, and to use this digital image (pixel map) to automatically identify the spatial coordinates of the square-shaped (pixel-shaped) partial subconduc-tors into which the different conductors can be subdivided. Image resolution, penetration depth, edge detection, as well as area and geometric mean distance error reduction tech-niques are then used to reduce the dimensions of the problem and improve the accuracy of the results. The idea of the partial subconductors method is to transform the system of conductors into an equivalent network of subconductors. A set of coupled circuit equations represents the equivalent network, and "bundling techniques" are used to obtain the parameters of the cable system. These "bundling techniques" transform the equations of the subconductors into the equations of the conductors by row and column operations. Since the number of subconductors increases noticeably with the frequency, an algo-rithm to partition the impedance matrices is proposed. The proposed methodology adapts itself to the physical memory of the computer, thus allowing the program to partition the partial subconductor impedance matrices when their sizes exceed the available physical memory. Coaxial cables, buried cables, and sector-shaped cables are studied with the proposed technique, and the results are compared with those obtained from the analytic method, the finite element method, and Ametani's approximate method. 11 A B S T R A C T 111 Spanish version: Este trabajo de tesis presenta el desarrollo de un algoritmo para el calculo de parametros de cables de potencia que utiliza imagenes digitales para discretizar la geometria del cable y el metodo de los circuitos equivalentes de subconductores parciales para determinar los parametros. La idea de la metodologfa propuesta para la discretizacion de la geometria es dibujar la section transversal del cable, o escanear su fotografia, y utilizar dicha imagen digital para identificar en forma automatica las coordenadas espaciales de los subconductores parciales cuadrados (pfxeles) que dividen los diferentes conductores. Para mejorar la exactitud de los resultados y reducir las dimensiones del problema, la metodologfa propuesta considera aspectos tales como la reduccion de la resolution de la imagen, la profundidad de pene-tration de las ondas electromagneticas, la detection de los bordes de los conductores y la reduccion de los errores de discretizacion asociados con las areas y las distancias medias geometricas de los conductores. La idea del metodo de los subconductores parciales es transformar el sistema de con-ductores en una red equivalente de subconductores. De este modo ecuaciones de cir-cuitos acoplados describen la red equivalente y tecnicas de "bundling" permiten obtener los parametros. Estas tecnicas de "bundling" transforman las ecuaciones acopladas de los subconductores en las ecuaciones acopladas de los conductores operando filas y columnas. Como el numero de conductores incrementa notablemente con la frecuencia, la tesis propone un algoritmo para particionar las matrices de impedancia. E l algoritmo propuesto se adapta a la memoria ffsica del computador, particionando las matrices de impedancia de los subconductores parciales cuando sus tamanos exceden la memoria fisica disponible. La tecnica propuesta se utiliza para calcular los parametros de cables coaxiales, cables con forma de sector circular y cables enterrados, y los resultados se comparan con aque-llos del metodo analitico, el metodo de los elementos finitos y el metodo aproximado de Ametani. Contents Abstract ii List of Tables v i i List of Figures v i i i Acknowledgements x 1 Introduction and Overview 1 1.1 Background 1 1.2 Thesis Motivation 2 1.3 Summary of the Thesis 5 2 Literature Review on Cable Parameter Calculations 10 2.1 Introduction 10 2.2 Classification of Power Cables 10 2.3 Literature on Cable Parameter Calculations 12 2.3.1 Analytic Method 13 2.3.2 Finite Elements 17 2.3.3 Partial Subconductor Equivalent Circuits 19 2.4 Summary 23 3 M e t h o d s for the Calculat ion of Frequency-Dependent Parameters of Power Cables 25 3.1 Introduction 25 3.2 Basic Assumptions 26 3.3 Analytic Method 27 3.4 Approximate Method 35 iv C O N T E N T S _v 3.5 Finite Element Method 36 3.6 Earth-Return Impedances 37 3.6.1 Uribe's Methodology 40 3.6.2 Wedepohl's Methodology 42 3.6.3 Comparison of Methods 43 3.7 Summary 45 4 Geometry Discretization: Proposed Methodology 47 4.1 Introduction 47 4.2 Analysis of the Modeling Problem 48 4.3 Spatial Discretization 50 4.4 Developed Algorithm 60 4.4.1 Processing and Segmentation of Input Image 61 4.4.2 Rediscretization for High Penetration Depths 63 4.4.3 Rediscretization for Low Penetration Depths 67 4.4.4 Storage of Data 72 4.5 Summary 72 5 Calculat ion of Frequency-Dependent Parameters of Power Cables: P r o -posed Methodo logy 73 5.1 Introduction 73 5.2 P S E C Method 74 5.2.1 Method Formulation 74 5.2.2 Equations Modification 77 5.2.3 Equations Reduction 78 5.3 Partition Methodology 79 5.3.1 Partition for Matrix Construction 80 5.3.2 Partition for Matrix Modification 81 5.3.3 Partition for Matrix Reduction 84 5.3.4 Subpartition for Storage and Retrieval of Elements 87 5.4 Summary 88 6 Case Studies 89 6.1 Introduction 89 C O N T E N T S yj. 6.2 Coaxial Cable 90 6.2.1 Simulation a 90 6.2.2 Simulation b 93 6.2.3 Simulation c 95 6.2.4 Summary of Results 98 6.3 Buried Cable System 100 6.4 Sector-Shaped Cable I l l 6.5 Summary 116 7 Conclusions and Recommendations for Future Research 118 7.1 Conclusions 118 7.2 Recommendations for Future Research 120 List of Symbols 124 Bibl iography 132 A Derivat ion of the Frequency-Dependent Parameters of a C y l i n d r i c a l C o n -ductor 139 B Wedepohl's Formulae for the Frequency-Dependent Parameters of T u b u -lar Conductors 146 C Capacit ive Susceptances Associated with Underground Tubular Conduc-tors 148 D Fini te Element M e t h o d D . l Calculation Procedure D.2 Derivation of the Basic Formulae 150 150 152 List of Tables 4.1 Resolution as a function of frequency for core of coaxial cable 71 4.2 Resolution as a function of frequency for sheath of coaxial cable 71 4.3 Resolution as a function of frequency for skin layer of core (coaxial cable). 72 5.1 Hard disk time versus record length 87 6.1 Simulation time versus number of subconductors (coaxial cable, Simulation c) 98 6.2 Loop resistance as a function of frequency for coaxial cable 99 6.3 Loop inductance as a function of frequency for coaxial cable 99 6.4 Differences in percentage terms between corrected impedances (Wedepohl's corrections with respect to Pollaczek's corrections) 104 6.5 Self resistance as a function of frequency for sector-shaped cable 112 6.6 Self inductance as a function of frequency for sector-shaped cable 115 6.7 Mutual resistance as a function of frequency for sector-shaped cable 115 6.8 Mutual inductance as a function of frequency for sector-shaped cable. . . . 115 vii List of Figures 3.1 Cable geometry for analytic (classical) method 28 3.2 Equivalent network for single-circuit a (similar for single-circuits b and c). . 29 3.3 Mutual earth-return impedances 29 3.4 Cable geometry for earth-return impedances 39 3.5 Self earth-return impedance 44 3.6 Mutual earth-return impedance 45 4.1 Cross section of single-core coaxial cable 51 4.2 Effect of changing the spatial resolution 52 4.3 Coordinate axes for discrete image 53 4.4 Integer matrix for discrete image 53 4.5 Discretization errors for core 57 4.6 Discretization errors for sheath 58 4.7 G M D errors after corrections 59 4.8 Developed algorithm for geometry discretization 60 4.9 Block of processing and segmentation of input image 61 4.10 Segmentation of coaxial cable 62 4.11 Block of rediscretization for high penetration depths 63 4.12 Subblock of rediscretization according to thickness of conductor 63 4.13 Subblock of iteration for G M D error minimization 64 4.14 Block of rediscretization for low penetration depths 68 5.1 Subconductors enclosed by ring 74 6.1 Resistance and inductance of coaxial cable (Simulation a) 92 6.2 Resistance and inductance of coaxial cable (Simulation b) 94 6.3 Resistance and inductance of coaxial cable (Simulation c) 96 viii LIST O F F I G U R E S _ix 6.4 Resistance and inductance of coaxial cable (Simulation c, log scale) 97 6.5 Geometry of buried cable system 101 6.6 Self core impedance Z c c of buried cable system 102 6.7 Mutual core-sheath impedance Z c s of buried cable system 103 6.8 Self sheath impedance Z s s of buried cable system 105 6.9 Mutual impedance Z a b = Z b c of buried cable system 106 6.10 Mutual impedance Z a c of buried cable system 107 6.11 Self core impedance with return through ring (coaxial cable) 108 6.12 Self sheath impedance with return through ring (coaxial cable) 109 . 6.13 Mutual core-sheath impedance with return through ring (coaxial cable). . . 110 6.14 Geometry of sector-shaped cable I l l 6.15 Self resistance and self inductance of sector-shaped cable 113 6.16 Mutual resistance and mutual inductance of sector-shaped cable 114 A . l Solid-wire skin effect quantities compared with dc values 145 Acknowledgements I have had help from many people. My apologies to anyone I miss. I gratefully acknowledge the financial assistance provided by "La Universidad Simon Bolivar (USB)" and "El Consejo Nacional para el Desarrollo de las Investigaciones Cientifi-cas y Tecnologicas (CONICIT)", both in Caracas, Venezuela. I also acknowledge the support I received from my colleagues of the "Departamento de Conversion y Transporte de Energia (CTE)" at USB, the staff of the "Direction de Desarrollo Profesoral (DDP)" at USB, and the staff of Human Resources at CONICIT. Big thanks to my supervisor Professor Jose Marti for his excellent guidance and support along the last four years. His continuous advice and his numerous ideas greatly improved this research work. I also appreciate that he trusted me and gave me the opportunity to pursue my Ph.D. studies here at U B C . Special thanks to Professor Hermann Dommel for supplying his valuable expertise, technical literature, and support routines (especially the one for matrix reduction and the code lines for bundling of conductors; both were very helpful). I am honored he is on my committee. Thanks to Professor Matthew Yedlin for bringing back my knowledge of Electromag-netism with his lectures on Fields and Waves. I am also grateful for the ideas he gave me for my research and for his participation in my committee. Thanks to Dr. Takahide Niimura for introducing me to the world of Fuzzy Logic and Artificial Intelligence as well as for his continuous interest in the progress of my work. Thanks to Dr. Jose Luis Naredo and Felipe Alejandro Uribe, from C I N V E S T A V in Mexico, for providing me with their papers on earth-return impedances. I also appreciate they kindly gave me their routines for earth-return impedance calculations and clarified all my doubts and questions. Thank you to my friend Mazana Lukic for proofreading the first draft of this thesis work. I am truly thankful for her kindness and useful suggestions, and also glad of having her friendship and trustfulness. Thanks to my friend Luis Linares, who gave me the idea of developing the partition algorithm in 1998. I also thank him for helping me during my admission process in U B C , for lending me his powerful dictionaries, for introducing me to the world of WPgX., and for the ski lessons. Thank you also to his wife Maria Luces for inviting me to sing carols. Thanks to my friend Meliha Selak for all her support here in Vancouver. Also to my friend and coworker Paloma de Arizon (and her husband Nelson Bacalao), and to my friends Fernando Moreira and Benedito Bonatto (and their respective families). Thanks to Jesus Calviho for giving me the idea of the bitmapped files in 1997, to Ting-Chung Yu for those useful discussions on cables, and to Daniel Lindenmayer for providing the template used to write this thesis. Thanks also to all the other faculty members, colleagues and friends in the U B C Power Group. I have learned a lot from all of them. x Acknowledgements xi Especially, I want to thank my parents (Oly Marlene and Marcos Arturo), grandparents (Mita, Celia and Maximo) and siblings (Roger and Carolina) for the love and support I have received from them during all my life. My achievements have always been theirs. Also I thank my aunts, my cousins, my wife's family and my friends in Venezuela for their usual willingness. In addition, I want to thank my elementary, high school and university teachers. Many of them made important contributions to my life. Last but not least, I would like to say thanks to my beloved wife Martha, who has always been very supportive. Without her love and organizational skills this adventure would not have come true. In fact, this thesis is a joint effort of the team Martha-Richard. I love her very much and am indebted to her for all she has done for the last nine years. Richard A. Rivas May 11, 2001 Chapter 1 Introduction and Overview 1.1 Background Power cables are commonly used to transmit and distribute electric energy in those situa-tions where overhead transmission lines are not applicable for technical and environmental restrictions. Typical examples are underground distribution networks in major cities and submarine interconnections. To develop and maintain power cable systems certain steps must be followed. The cable systems must be planned, designed, operated and maintained following rigorous standards. Those standards include electric analyses of the cable system operation both in steady state and during transient periods. In addition, adequate models of the conductors must be used to study the different operating conditions. The models can be either scaled-down physical prototypes or virtual models. If they are virtual models and the steady-state analysis is the objective, ampacity-, load flow- and short circuit-type programs are required. If the analysis of transient periods, either during switching operations or surge lightning is the target, EMTP-type programs [17] are employed. Significant effort has been made to incorporate power system component models into 1 1.2. Thesis Motivation 2 EMTP-type programs [19]. As far as the modeling of transmission lines and underground cables is concerned, J. Marti [38], L. Marti [39], Castellanos [11], Gustavsen [28] and Noda [47] have made important contributions. For instance, taking the frequency-dependent nature of the electric parameters into consideration, they have developed transmission line and power cable models suitable for time-domain EMTP-type environments. Before setting up any cable model, the cable parameters must be calculated over an extended range of frequencies. However, given the wide variety of geometric shapes, ma-terials used, types of arrangement, as well as the presence of strong skin and proximity effects among conductors, power cable parameter calculations can be a difficult task. For example, the cables can be single-phase or three-phase, stranded or solid, and circular (with or without sectors) or rectangular [34]. For power rails of traction systems the shapes are more irregular and the analysis becomes more complex. 1.2 Thesis Mot i va t ion The reasons that motivated the realization of the present thesis work were stated by Dommel in his E M T P Theory Book [19]. They are summarized as follows: • Given the broad range of cable shapes and designs, it has been considered very complicated, if not impossible, the creation of an algorithm for the calculation of the frequency-dependent parameters of arbitrarily-shaped cables. As a result, a general-purpose cable-parameter program for EMTP-type studies is not still available. • Available E M T P support routines, such as C A B L E C O N S T A N T S and C A B L E PA-R A M E T E R S [3], make it possible to study power coaxial cables and high-voltage pipe-type cables. Since the mathematical expressions used by those programs for pipe-type cables are approximate, proximity effects are neglected in eccentric cases with cradled or triangular configurations. • There is no E M T P support routine for the parameter calculation of distribution 1.2. Thesis Motivation 3 cables with concentric neutral conductors, pipe-type cables with sector-shaped con-ductors, busbar arrays, or power rails of traction systems. • Very frequently, cables are installed in underground tunnels, which makes it difficult to apply the traditional formulae for the earth-return impedances of cables directly buried. Two main approaches have been suggested in the past for the calculation of the frequency-dependent parameters of cables: analytical solutions and numerical ones. The numerical methods have also been divided into approximations based on: finite elements and partial subconductor equivalent circuits. Analytical solutions are typically used for cylindrical conductors, whereas numerical solutions are usually applied to arbitrarily-shaped conductors. The finite element method is a powerful tool widely used in the fields of civil, mechan-ical, and electrical engineering. Given the current development of standalone software in the area of mesh generation [49], [50], the complexity of the solution region is not a diffi-culty for the finite element method. Nevertheless, this method needs an explicit modeling of the entire solution region, which for the case of cable parameter calculations includes in addition to the conductors and the insulation, the earth-return path with infinite dimen-sions. The earth-return path usually has a poor conductivity and, as a result, a large penetra-tion depth for electromagnetic waves, even at high frequencies. Therefore, a disadvantage of the finite element method is that it requires a large number of elements for the repre-sentation of the earth-return path. In addition, finite element programs can be expensive given their wide range of possible applications. They also require skilled users with the knowledge of the problem setup, which involves a formulation based on electromagnetic field analysis, a discipline that practicing power engineers are less familiar with. In 1965, Graneau [25] proposed using square-shaped partial subconductor equivalent circuit techniques for the calculation of parameters and current distributions of conductors with any cross section. The finer the subdivision, the closer will the array be to the original 1.2. Thesis Motivation 4 cross sectional area, and the more accurate the current distributions and conductor pa-rameters will become. This, however, implies more computing time to obtain the solution, especially when the penetration depth is smaller than the dimension of the conductors. The formulation of the partial subconductor equivalent circuit method is simpler than that of the finite element method as the cable parameters are found using circuit analysis theory. The partial subconductors method only requires an explicit modeling (subdivision) of the conductors, and the cable parameters can be obtained by "bundling techniques" [29]. The "bundling techniques" transform the equations of the coupled subconductors into the equations of the coupled conductors by row and column operations. Circuit analysis theory, moreover, is an area power engineers feel more comfortable with. Nonetheless, given an arbitrarily-shaped configuration, the partial subconductor equiv-alent circuit method still requires considerable user intervention. This intervention is es-sential in the setup of the subconductor coordinates and in the definition of appropriate shapes and adequate resolutions for the subconductors. For instance, round, square, or concentric subdivisions could represent the best option in a particular case, with more or fewer subconductors needed depending on the frequency under study. A main advantage of the partial subconductors method over the finite element method is that the earth return does not have to be explicitly subdivided when the conductor parameters are calculated. After assuming a fictitious return path of infinite conductivity, the ground return can be incorporated into the analysis employing the traditional formulae for earth-return representation [10], [15], [54], [71]. Using the concept suggested by Graneau [25], [26], i.e., the subdivision of arbitrarily-shaped cross-sections into square-shaped subconductors, this thesis is aimed at developing a general-purpose cable-parameter algorithm based on the automation and enhancement of the geometry discretization process. The basic idea of the proposed methodology is to draw the cable geometry, or scan its photograph, and to use this digital image (pixel map) to automatically identify the spatial coordinates of the square-shaped (pixel-shaped) partial subconductors into which the different conductors can be subdivided. After that, 1.3. Summary of the Thesis 5 the conductor parameters can be obtained from the representation of the subconductors as sets of coupled circuits whose currents depend on the skin and proximity effects of the system. The input file stores the image sampling (discretization of the spatial coordinates) and the color quantization (color selection). Subsequently, basic digital image process-ing/computer graphics concepts such as connectivity among pixels and image frequency distribution histograms can be used to trace and quantify the different connected compo-nents [24] (e.g., conductors and insulation) and the number of subconductors (pixels) of each one. In summary, the main advantages of the proposed methodology for the study of arbitrarily-shaped cables are simplicity and economy. Simplicity of the mesh generation and of the parameter calculation method, and economy of software resources. The sim-plicity comes from using a square-shaped subconductor grid and an equivalent electrical network for all the possible geometries, whereas the economy comes from using off-the-shelf computer graphics formats such as bitmap [45], which can be created by many drawing and photography programs [2]. The simplicity of the proposed methodology, however, leads to implementation prob-lems such as high geometry discretization errors and large densities of subconductors. This thesis proposes solution algorithms for such difficulties. 1.3 Summary of the Thesis This thesis work deals with power cable parameter calculations. The cross sections of the conductors are treated as arbitrarily-shaped cable arrangements, and the frequency-dependent series parameters (resistance and inductance) are calculated from the subcon-ductor coordinates and the electric properties of the material (conductivity and perme-ability). These parameters can be used to synthesize the time-domain equivalent networks needed for the representation of cable systems in EMTP-type programs. 1.3. Summary of the Thesis 6 Chapter 2 reviews the technical literature on the topic of power cable parameter calcu-lations. The information is organized by date and type of method (analytic, finite elements, and partial subconductors), starting with the work of Carson [10] and Pollaczek [54] in 1926, and finishing with the work of Papagiannis [51] and Wang [69], [70] in 2000. Even though more information on the topic can be found in the technical literature, the sum-mary presented here provides a considerable amount of detail on the historical evolution and the state of the art in this field of knowledge. This chapter also presents a section on the classification of power cables. Chapter 3 develops the mathematical formulations of the alternative methods for the calculation of the frequency-dependent parameters of power cables. The methodologies described are the analytic (classical) method [3], [19], [63], [71], the approximate method [4], and the finite element method [75], [76]. The classical method formulates the loop equations of an underground cable system as functions of the surface impedances of the conductors and the impedances of the insulation layers. Then, it transforms the loop equations into conductor equations through row and column operations. The approximate method assigns the area and perimeter of arbitrarily-shaped conduc-tors to equivalent cylindrical conductors. As a result, the classical method can be applied to the equivalent system of conductors. The finite element method gives trial functions to the basic equations relating the current densities of the conductors with the magnetic vector potential at the element nodes, and then it forces a residual to be zero. Thus, the impedances can be obtained from the solution of a set of linear equations and from the ratios between electric field intensities and currents injected into the conductors. Chapter 3 also deals with the topic of earth-return impedances for underground cables. The methodologies described are Pollaczek's method [54] and Wedepohl's method [71], which assume a semi-infinite earth-return path extending itself from side to side and from the earth surface downward. Pollaczek's method uses Bessel functions and infinite inte-grals for the calculation of the impedances, whereas Wedepohl's method uses closed-form 1.3. Summary of the Thesis 7 formulae obtained from the expansion of Pollaczek's formulae into series. Uribe's method [68] for the numerical integration of Pollaczek's integral is also explained, and the results of the numerical integration are compared with those of Wedepohl's method. Chapter 4 introduces the proposed methodology for the discretization of power ca-ble geometries. The cable geometry is obtained from a bitmapped file, which provides a rediscretization program with the data on the dimensions and coordinates of the subcon-ductors (pixels). Later on, this data is used by a partial subconductors program to obtain the geometric mean distances of the subconductors. Information on the mathematical representation of the edges, however, is not contained in the bitmapped file storing the image discretization. As a result, staircase effects on the conductor peripheries, jagged shapes, as well as cross sectional areas and geometric mean distances different from the actual ones (phenomenon called aliasing in Computer Graphics) arise after carrying out the rasterization.1 The area errors are compensated for by dividing the actual area among the subconduc-tors, as proposed in [69], [70], whereas the geometric mean distance errors are corrected by an iterative procedure. The proposed iterative procedure keeps the number of sub-conductors constant and varies the resolution in discrete steps until detecting a change of sign in the geometric mean distance error. Then, the method uses linear interpolation to calculate a resolution closer to the crossing through zero of the error function. The use of this interpolated resolution reduces the geometric mean distance discretization errors noticeably. For simple geometries, the actual cross sectional areas and the actual self geometric mean distances of the conductors can be supplied by the user. Furthermore, more com-plicated shapes can be subdivided into simpler shapes, e.g., triangles, so that the actual areas can be approximated to the sum of areas of the triangles. Similarly, the actual self geometric mean distances can be obtained by bundling the triangles. To take skin effects into consideration, the subdivisions must have the same order 1 T h e process of transforming an image represented with mathematical equations into an image repre-sented with pixels. 1.3. Summary of the Thesis 8 of magnitude as the penetration depth associated with the frequency, the conductivity, and the permeability of the conductor. Hence, the number of subconductors increases noticeably as the penetration depth diminishes. To reduce memory requirements and simulation timings, Chapter 4 proposes segment-ing the input image so that the resulting images of the conductors can be rediscretized separately. Thus, low resolutions adapted to the thicknesses of the conductors, or higher resolutions adapted to the low penetration depths of high frequencies, can be assigned to the different images of the conductors. This chapter also develops a method to model only the boundaries (edges) of the conductors when the penetration depths are very small. As a result, the dimension of the problem is reduced considerably. In addition, criteria to choose the appropriate resolutions for the representation of rectangular and circular conductors are suggested. Chapter 5 discusses the details of the proposed partial subconductor equivalent circuit method and its implementation in a general-purpose cable-parameter program. The basic idea of the method is to transform the system of conductors into an equivalent network of subconductors. This equivalent network is then represented by a set of coupled circuit equations whose solution yields the current distributions in the subconductors and the parameters of the cable system. To reduce the number of arithmetic operations the program bundles the subconductors using a partial Gaussian elimination algorithm which exploits the symmetry of the system [7], [19]. To overcome the memory problems, an algorithm to partition the impedance ma-trices is proposed. The proposed methodology adapts itself to the physical memory of the computer, thus allowing the program to partition the partial subconductor impedance ma-trices when their sizes exceed the available physical memory. The proposed algorithm also reduces the interaction between hard disk and C P U and allows the program to calculate the frequency-dependent parameters in low R A M computers. Chapter 6 studies typical cable arrangements with the developed techniques and com-pares the results with those of the classical method, the finite element method, and the 1.3. Summary of the Thesis 9 approximate method. The cases under study are a coaxial cable, a buried cable system and a sector-shaped cable. It is shown that the results of the proposed method are in good agreement with those of the classical method and the finite element method. Finally, Chapter 7 presents the conclusions of the thesis. This last chapter also dis-cusses the contributions of the work and gives recommendations for future research. These contributions can be summarized as follows: • It is proved that it is feasible to perform an accurate calculation of the frequency-dependent parameters of arbitrarily-shaped power cables with digital images for hu-man vision. • A technique based on a correction of the geometric mean distances is introduced to decrease the inductance discretization errors. • A n adaptive resolution scheme is proposed to decrease the number of subconductors needed according to the penetration depth. At low frequencies the proposed scheme assigns low resolutions to the conductors, whereas at higher frequencies it assigns higher resolutions. For very low penetration depths at high frequencies, the adaptive scheme switches to a boundary representation of the conductors. The scheme reduces the requirements of simulation time and memory noticeably. • A matrix partition methodology is proposed that allows the partial subconductors program to trade solution time for memory when the sizes of the matrices exceed the available physical memory, even in low R A M computers. Chapter 2 Literature Review on Cable Parameter Calculations 2.1 In t roduct ion Much research has been done to address the problem of calculating the parameters of power cables. It is the purpose of this chapter to present a general classification of power cables and a literature review on the research done on cable parameter calculations. First, Section 2.2 classifies power cables taking into account the type of cable construc-tion and the method of installation. Then, Section 2.3 describes the developments in the field of cable parameter calculation since 1926. This review is divided into three areas: the analytic method, the finite element method, and the partial subconductor equivalent circuit method. Finally, Section 2.4 gives a summary of the chapter. 2.2 Classi f icat ion of Power Cables It is difficult to classify power cables given the wide variety of shapes, arrangements, material used, and voltage levels that exist in the power cable industry. In spite of that, classifications based on type of cable construction and method of installation have been proposed [5], [23], [59]. The aforementioned references give much information on this subject, so this section only presents a brief summary. 10 2.2. Classification of Power Cables 11 From the point of view of type of cable construction, cables can be grouped considering aspects such as conductor material (copper or aluminum), conductor construction (con-ventional stranding or stranding for compact design), type of electrical insulation, presence of sheath or concentric neutral conductors, and presence of armor. For example, considering the type of electrical insulation, cables can be divided into: • Paper-insulated cables. • Extruded-dielectric cables. • Gas-insulated cables. Paper-insulated cables can be split into two groups: High-Pressure Fluid-Filled ca-bles (HPFF), also called Pipe-Type cables (PT), and Self-Contained Liquid-Filled cables (SCLF). H P F F cables are usually enclosed by pipes, while S C L F cables are directly buried. In the latter, each single circuit, e.g., each circuit made up of core and sheath, is indepen-dent and sheath- or/and armor-protected. Hence, S C L F cables, although more expensive than P T cables, are preferred in submarine crosses. Extruded-dielectric cables have lower dielectric losses than paper-insulated cables and, as a result, higher ampacities. Extruded-dielectric cables can be divided into three types: X L P E (crosslinked polyethylene), P E (low- or high-density polyethylene), and E P R (ethyl-ene-propylene rubber insulation). Extruded-dielectric cables require less amount of acces-sories, use simple taped- or fused-joints, and need no pressurizing plants or reservoirs nor cathodic protection. Therefore, simpler design and installation, as well as lower cost and competitiveness, have turned the extruded-dielectric cable into an attractive choice. Gas-insulated cables use sulfur hexafluoride (SF 6) as the insulation between the con-ductors and the shields. They have high diameters and therefore higher ampacities. How-ever, they are restricted to short distances and high voltages for cost reasons. They are traditionally used in busbar arrays and bays of gas-insulated substations. From the point of view of the method of cable installation, cables can be grouped considering aspects such as laying conditions, bonding arrangements, and cooling method. 2.3. Literature on Cable Parameter Calculations 12 Considering the laying conditions, and considering the possible combinations, the fol-lowing installations result: • Cables surrounded by air, outdoors (e.g.: transitions between overhead lines and underground cables). • Cables surrounded by air, indoors (e.g.: cables in trays) • Cables surrounded by air, installed in an underground tunnel. • Cables in duct, and the duct is surrounded by air. • Cables directly buried. • Cables installed in an underground thermal backfill. • Cables installed in ducts in an underground duct bank. • Cables installed in a pipe and the pipe is directly buried. • Cables installed in a pipe and the pipe is in an underground thermal backfill. 2.3 L i terature on Cable Parameter Calculat ions Two main approaches have been suggested in the past for frequency-dependent parame-ter calculations: analytical solutions and numerical ones. The numerical methods have also been divided into approximations based on: finite elements and partial subconductor equivalent circuits. Techniques based on analytical solutions, on the one hand, are traditionally used to analyze simple geometries, e.g., concentric cylindrical conductors. Techniques based on nu-merical procedures, on the other hand, are usually employed to study cables with arbitrary shapes and configurations. 2.3. Literature on Cable Parameter Calculations 13 2.3.1 Analytic Method This method calculates the parameters solving analytically the field equations associated with the propagation along the conductors and the earth-return path. The solution typ-ically involves Bessel functions, which have to be evaluated numerically through series expansions. As far as the calculation of earth-return impedances, very useful work has been done by Carson [10] and Pollaczek [54] (1926), who developed the formulae for the earth-return corrections of overhead and underground conductors, respectively. This formulae stem from representing the conductors as infinitely long, infinitely thin filaments, and require the evaluation both of Bessel functions and infinite integrals. Later on (1934), Schelkunoff [63] derived the surface impedances of concentric coax-ial conductors, thus connecting the theory of propagation along cylindrical wires having external coaxial returns to circuit analysis theory. On the assumption of T E M propagation, Schelkunoff transformed the field equations in cylindrical coordinates into transmission line equations and expressed the solution of the electric and magnetic field intensities as Bessel functions with unknown coefficients. The coefficients are calculated from the magnetic field intensities on the conductor surfaces by visualizing the current injected into a tubular conductor as the sum of two components, a component returning through an inside return and another component returning through an outside return. Thus, the inner and outer surface impedances of a tubular conductor can be obtained from the terms in the equations relating the electric field intensities in the conductor surfaces with the return currents. For the calculation of the surface impedances associated with arbitrarily-shaped con-ductors, Schelkunoff recommended a field solution based on the integration of the complex Poynting Vector over the surfaces of the conductors. Wiseman [73] (1948), Meyerhoff and Eager [42] (1949), and the A I E E Working Group on ac resistance of pipe-cable systems [1] (1952) developed semi-empirical formulae for the 2.3. Literature on Cable Parameter Calculations 14 calculation of 60-Hz resistances and losses of high-voltage pipe-type cables. The main idea was to modify formulae for conductors in air so that the calculated values were similar to the test values, thus taking factors such as skin effect, proximity effects, shielding and pipe presence into consideration. Meyerhoff's formulae, for example, are still used to calculate the losses of steel pipes in ampacity studies [5]. Lewis and Tuttle [36] (1959) presented a thorough study on the calculation of the parameters of stranded conductors. Aside from explaining the geometric mean distance method and giving tables for the determination of the resistance and the reactance of Aluminum conductors, steel reinforced (ACSR), the main contribution of their paper was to use Bessel functions for the calculation of the frequency-dependent resistances of tubular conductors in the frequency range of Power Transmission Engineering. They also compiled information of the state of the art at that moment in the area of parameter calculations of stranded conductors. Their observations can be summarized as follows: • Parameters of stranded conductors are usually calculated assuming straight strands and a 2% allowance to take into account an increased length due to spiraling. The effect of the central steel core is considered negligible in A C S R since the outer layers of aluminum strands are spiraled in opposite direction. • Experimental results show that the skin effect for stranded conductors can be an-alyzed using equivalent solid tubular conductors with the same thickness and dc resistance as the group of layers of stranded conductors. For ACSR, for example, it is sufficient to use equivalent tubes enclosing only the aluminum strands. • Spiraling decreases ac resistance and increases ac reactance very slightly, as compared to the case of a stranded but unspiraled conductor. The correction of such differences is therefore negligible. • For A S C R having only one layer of aluminum strands the effect of the central steel core is not negligible. The magnetization of the core is not neutralized since no opposite layers are present and, therefore, hysteresis and eddy currents appear. In such a case, resistance and reactance increase noticeably and the parameters should 2.3. Literature on Cable Parameter Calculations 15 be corrected using the active and reactive losses of the steel wires obtained from experimental data. With respect to the automatic calculation of the parameters, Hesse [29] (1963) de-veloped a "bundling technique" to combine and eliminate subconductors in parallel and ground wires through matrix operations. Programs for conductor parameter calculations such as M T L I N E [44], for example, use a similar technique to transform multi-circuit overhead transmission line geometries into equivalent reduced systems. On the assumption that the current densities are uniform, Smith and Barger [65] (1972) studied the 60-Hz impedances of underground distribution cables using the geometric mean distance technique. They presented two methods for the calculation of the self, mutual and sequence impedances of underground distribution cables with multi-wire concentric neutral conductors. One method keeps the identity of each neutral wire when formulating the coupled equations, and the other replaces the neutral wires with an equivalent concentric sheath. The earth-return corrections employed were based on Carson's formulae [10], and the results obtained from both methods were similar at power frequency. For the computation of the sequence impedances of this type of cable, EPRI recommends Smith's formulae [21]. Wedepohl and Wilcox [71] (1973) applied Schelkunoff's formulae for surface impedances and Pollaczek's formulae for earth returns [54] to the calculation of the frequency-dependent earth-return impedances of buried cable systems. A n important contribution of this paper was to have approximated Schelkunoff's formulae (based on Bessel functions) and Pol-laczek's formulae (based on Bessel functions and infinite integrals) to simpler closed-form expressions. Bianchi and Luoni [8] (1976) analyzed impedances, losses and current distributions in submarine cables. The main contribution of their work was to visualize the sea re-turn as a medium with infinite dimensions. As a result, the self and mutual sea-return impedances can be calculated with Schelkunoff's formula for the inner surface impedance on the assumption that the sea is a cylindrical conductor of infinite outer radius. 2.3. Literature on Cable Parameter Calculations 16 Ametani [3] (1980) developed E M T P support routines [19] for the calculation of the series and shunt impedances of single-core coaxial cables and high-voltage pipe-type cables. For the series impedances of coaxial cables Ametani uses SchelkunofT's formulae, whereas for the series impedances of pipe-type cables he uses Tegopoulos's formulae and Brown's formulae [19]. Dommel (1986) gives a very good description of the two methodologies in his E M T P Theory Book [19], pointing out that for the case of pipe-type cables the formulae neglect proximity effects among conductors. Both methods are widely used by the EMTP-user community. More recently (1996), Ferkal, et al. [22], developed an analytic method to calculate current distributions and losses taking skin and proximity effects into consideration. Given a system with two neighboring cylindrical conductors, the method uses Maxwell's equations and finds the field distribution considering that one of the conductors is replaced by an infinitely thin filament which carries the return current. A similar technique has also been used by Kane, et al. [33] (1995), to calculate the self and mutual impedances of pipe-type cables. As a result of including the proximity effects in the formulation, the magnetic field intensity, the magnetic vector potential, and the volume current density become functions not only of the radial distance r, but also of the angular position 9. Then, the integration coefficients of the formulae can be found examining the continuity of the tangential com-ponents of H and the continuity of the radial components of B at the boundaries given by the outer and inner radii of the conductors. The final formulae for the current densities and the impedances are functions of r, 6, Bessel functions, and infinite series. The method incorporates more conductors into the analysis applying the principle of superposition, and its results are in good agreement with those of the finite element method. Parsi-Feraidoonian and Dommel [53] (1992) applied the "bundling technique" to the calculation of the impedances of cables with semi-conductive layers. The proposed method is based on treating the semi-conductive layers of core and sheath as additional conductors, and on writing the loop equation between conductors and semi-conductive layers as if an insulation layer were in the middle. 2.3. Literature on Cable Parameter Calculations 17 Since the semi-conductive layers are in parallel with the core conductors or conducting sheaths, and the impedance of the insulation between a semi-conductive layer and a con-ductor is equal to zero, the semi-conductive layer can be eliminated and combined with the conductor in a similar manner as it is done for ground wires and conductor bundles of overhead transmission lines. Given the low conductivities of the semi-conductive layers, the impedances obtained from their inclusion in the analysis are quite similar to those obtained from their treatment as insulation. Ametani and Fuse [4] (1992) developed an approximate method to transform arbitrarily-shaped conductors into equivalent circular conductors. Thus, an arbitrarily-shaped mul-ticonductor system can be studied with the existing impedance and admittance formulae given in [3]. 2.3.2 Finite Elements The finite element method transforms a continuous physical system into a discrete system. The basic idea of the method is to break the system into discrete elements interconnected at discrete node points, and to approximate the field distribution throughout each element by functions whose coefficients can be found through Galerkin techniques. The method is widely used in areas such as structural mechanics, electrical field theory, and fluid mechan-ics [52], and was applied by Chari, Konrad and Weiss ('70s and '80s) to the calculation of electric parameters [75]. Yin and Dommel [75], [76] (1989, 1990) solved the problem of calculating the frequency-dependent parameters of arbitrarily-shaped power cables with the finite element method. The proposed technique finds the magnetic vector potential and the current density distri-bution of the solution region (conductors, insulation and earth-return path) after assigning zero values to a faraway boundary. The parameters are then calculated with the losses or the current densities of the conductors. A similar technique was simultaneously developed in Italy by Cristina and Feliziani [14] (1989). Yin analyzed cable geometries such as coaxial cables, buried three-phase cable systems, 2.3. Literature on Cable Parameter Calculations 18 tunnel-installed three-phase cable systems, sector-shaped cables, and stranded conductors. To obtain accurate results, Yin came to the conclusion that in the case of deeply buried coaxial cables rb should be greater than 36, whereas in the case of shallowly buried coaxial cables rb should be greater than 125, where rb is the radius of the boundary and 6 is the penetration depth of the earth. Since the penetration depth of the earth is large as compared to the dimension of the cables, even at high frequencies, large solution regions are required for the earth return. As a result, the computing time increases noticeably and the mesh generator faces difficulties subdividing the region and the details around the cables [75]. To reduce the dimension of the earth-return path, Yin proposed solving the system with non-zero value boundary conditions. On the assumption that the field distribution in the cable does not affect the field distribution in the earth, the non-zero value boundary conditions can be calculated transforming the cables into buried filaments. Thus, the field solution of the filaments and the partial earth-return current can be obtained from the formulae given by Pollaczek [54]. Yin concluded that the results of the non-zero value boundary technique are accurate if rb/6 is small at small r e , or if re/6 is small at large re, where re is the outer radius of the outermost insulation of the cable (equal to the inner radius of the earth-return path). For example, when re = 24 mm the differences between the non-zero value boundary technique and the zero value boundary technique are lower than one percent \irb/6 < 0.2. However, Yin obtained discrepancies of up to 20 percent between the results of Pollaczek's formulae and those of the finite element method. Yin also proved that tunnel-installed cables can be analyzed with Pollaczek's formulae using an approximate re. More recently, Satsios, et al. [61] (1998), used the finite element method to calculate electromagnetic fields and eddy currents of geometries with faulted overhead transmission lines and buried pipelines. Satsios assumed a stratified earth-return path and studied the influence on the calculations of parameters such as the depth of the first earth layer, the resistivities of the different earth layers, and the number of mitigation wires. 2.3. Literature on Cable Parameter Calculations 19 Hill, et al. [30] (1999), used the methodology to compute the electric parameters associated with power rails of transit systems. Hill analyzed the effects on the impedances of aspects such as current, frequency, track material properties, ground conductivity, and skin effect. Satsios, et al. [62] (1999), applied the technique to study the inductive interference caused to telecommunication cables by nearby ac electric traction lines. Satsios evaluated parameters such as the separation distance between the electric traction line and the buried telecommunication cable, the earth resistivity, and the number (and material) of mitigation wires. The results obtained were compared with measurements. Triantafyllidis, et al. [66] (1999), employed the method to determine frequency depen-dent parameters and sequence impedances of overhead transmission lines. Triantafyllidis compared the results of cases with terrain irregularities with those of the classical method, e.g., a line next to a mountain of variable slope, a line inside a canyon, and a line near a water region. Papagiannis, et al. [51] (2000), also applied the methodology to calculate the parameters of overhead transmission lines. Papagiannis modeled conditions such as stratified earth-return paths and terrain irregularities. 2.3.3 Partial Subconductor Equivalent Circuits The partial subconductor equivalent circuit method calculates the parameters of the con-ductors using circuit analysis theory. The basic idea of the method is to subdivide the cross sections of the conductors into thin subconductors so that the system of conduc-tors is transformed into a multi-circuit equivalent network. The equivalent system is then represented with a set of coupled circuit equations whose solution yields the current dis-tributions in the subconductors and the parameters of the cable system. Because the subconductors are sufficiently thin, constant and uniform current den-sities can be assumed in the subconductors, and dc resistances and inductances can be assigned to each of them. The dc resistances are obtained from the conductivities and areas per subconductor, while the dc inductances are calculated from the self and mutual 2.3. Literature on Cable Parameter Calculations 20 geometric mean distances among subconductors. Since there is mutual coupling among subconductors, the currents vary from subconductor to subconductor according to the skin and proximity effects present in the system. As opposed to the finite element method, the partial subconductors method does not require the modeling of the non-conducting regions. Nor is an explicit modeling (subdivi-sion) of the earth return path required (on the assumption that the conductors are infinitely thin filaments as compared to the earth-return path, the earth-return impedances can be incorporated later on using the classical earth-return formulae). The partial subconductors method is also a practical alternative for the treatment of cables with unusual shapes, where a rigorous field analysis (analytic method) is complicated due to the fact that the boundary conditions must hold at the surfaces of the conductors [64]. The method, however, produces full matrices and demands more computer resources as the frequency increases and the subdivisions become thinner. The methodology was proposed by Graneau [25], [26] in 1965, who suggested applying it to the study of induction problems, eddy current losses, rotating electrical machines, busbars, overhead lines, and cable conductors. In [25] Graneau used it to determine the current distribution in a straight rectangular busbar of finite length, whereas in [26] he used it to calculate current distributions and parameters of a sector-shaped power cable. In both cases the cross sections were subdivided into square-shaped partial subconductors. Silvester [64] (1968) used the method along with a modal transformation to study the frequency-dependent parameters of solid and hollow square conductors and Comellini, et al. [13] (1973), employed it to estimate quantities such as skin effect ratios of stranded and sector-shaped conductors, proximity effect ratios and current distributions in neigh-boring conductors, eddy currents in cable sheaths, and current distributions in earth-return paths. Comellini subdivided the conductors and the earth-return path into circle-shaped subconductors and carried out all the studies at 60 Hz. Lukas and Talukdar (1978) [37] applied the subconductors technique to the calculation of the frequency-dependent parameters of coaxial cables. In order to improve accuracy, 2.3. Literature on Cable Parameter Calculations 21 reduce computing time, and minimize storage requirements they proposed subdividing the conductors into arc-shaped subconductors which they called elementals. The elementals are thicker toward the core and thinner toward the edges of the conductors. They also derived and verified the formulae for the calculation of the self and mutual geometric mean distances of the arc-shaped subdivisions. In addition, they introduced the concept of the cut-off resistance to explain the discretization error at high frequencies. For example, in a core conductor this resistance is given by the dc resistance of the outermost layer of subconductors. Later on (1979), Weeks, et al. [72], used the subconductors methodology to compute the parameters of intercircuit wiring of modern circuit packages. Weeks derived a closed-form formula for the calculation of the mutual geometric mean distances of rectangular subconductors and proposed subdividing the outer layers of the cross sections with subcon-ductors having thicknesses of the same magnitude as the penetration depth. The results of the proposed calculation method were compared with measurements made on a large-scale model of striplines over a ground plane. The partial subconductors method, also known as Weeks' method or P E E C method in the field of VLSI circuits, is used extensively in the study of interconnections of integrated circuits [6], [31], [32], [56], [67], [74], [77]. De Arizon and Dommel [7] (1987) used the subconductors method to calculate the frequency-dependent parameters of coaxial cables, pipe-type cables, and stranded conduc-tors. For the coaxial cable they compared the results of the analytic method with those of circle-shaped subdivisions, square-shaped subdivisions and arc-shaped subdivisions, ob-taining the best results with arc-shaped subdivisions. They also analyzed the parameters of two eccentric conductors. Very useful contributions of this paper were to have applied the "bundling technique" to the reduction of subconductor matrices of power cables and to have minimized the number of arithmetic operations using a partial Gaussian elimination technique. Another contribution of this work was to have explained how to incorporate the earth-return path impedances into the conductor matrices. Cao, et al. [9] (1988), employed the subconductors technique to study busbars mounted on the arms of arc furnaces. Cao calculated current distributions, phase imped-ances, 2.3. Literature on Cable Parameter Calculations 22 bus-voltage drops, and power losses for various bus positions and arm configurations, and compared the results of the partial subconductors method with those of laboratory experiments. Tsuk and Kong [67] (1991) used the subconductors methodology to obtain the low-frequency parameters of interconnections between integrated circuits. To analyze arbitrari-ly-shaped conductors, they proposed subdividing the cross sections into triangular and polygonal subconductors. To derive the formulae for the self and mutual geometric mean distances of triangle-shaped and polygon-shaped subconductors, they employed Green's identities. For high frequencies, they elaborated a surface integral equation approach in which the resistance and inductance of each conductor were expressed in terms of the current and its normal derivative only on the surface of the wires, and in which the solution was obtained through the method of moments. Since the partial subconductors method (circuit approach) and the surface integral equation method (field approach) are combined to find the parameters in the entire range of frequencies, Tsuk's method is called hybrid. Zhou and Marti [77] [78] (1993, 1994) utilized the subconductors method to study pipe-type cables. They proposed subdividing the conductors with arc-shaped subconductors whose thicknesses are chosen following a sinusoidal rule. The sinusoidal rule is a function of the penetration depth 5 and yields finer subdivisions at the peripheries of the conduc-tors. The proposed technique also assumes a linear variation of the current density with matching values in the boundaries between subconductors. This procedure improves the accuracy of the results. However, it requires the system to be solved twice: A first time assuming uniform and constant current densities in the subconductors and a second time with the correction factors obtained from the first solution. The results of the partial subconductors method were compared with those of the finite element method. Antonini, et al. (1999) [6], also used the subconductors technique to calculate the frequency-dependent parameters and current distributions of rectangular conductors. For low frequencies Antonini suggested using subdivisions smaller than the penetration depth 5 divided by two, while for high frequencies he suggested using subdivisions smaller than 5 divided by ten. The results obtained show that the current density along the periphery of a 2.4. Summary 23 rectangular conductor is not constant, which does not occur in isolated circular conductors. As a result, the assumptions of a propagation constant equal to (1 + j)/5 and a high-frequency internal resistance equal to the high-frequency internal reactance are no longer valid, and, therefore, the resistances and inductances of rectangular conductors must be calculated separately. Wang (2000) [69], [70] used the subconductors method to calculate the impedances associated with power rails of dc traction systems. These impedances are needed to analyze harmonics in those systems. To take into account the skin effect Wang proposes using finer subdivisions for the outer subconductors. In addition, Wang gives a simple method to correct geometry discretization errors. He assigns, to each subconductor, a dc resistance n times greater than the actual dc resistance of the cross section, where n is the total number of equal subconductors into which the cross section of the power rail is subdivided. The formulae developed for the calculation of the dc resistances of the subconductors allow the method to consider cases where the outer subconductors are finer than the inner subconductors. 2.4 Summary A general classification of power cables and a review on the research done on cable pa-rameter calculations have been presented. The review shows that the methods of analytic solutions, finite elements, and partial subconductors have been studied extensively in the technical literature and that each of them presents advantages and disadvantages. The analytic method is appropriate for circular cases where the equations and the boundary conditions can be conveniently formulated. The finite element method is appropriate for arbitrarily-shaped cases. However, it demands more computer resources as the number of elements increases with the frequency. In addition, it requires the use of large solution regions for the representation of the earth-return path. 2.4. Summary 24 The partial subconductors method displays a simpler formulation and is also appropri-ate for arbitrarily-shaped cables. Nevertheless, it also demands more computer resources as the number of subconductors increases with the frequency. Chapter 3 Methods for the Calculat ion of Frequency-Dependent Parameters of Power Cables 3.1 Int roduct ion Significant effort has been made to model power cables in EMTP-type studies [11], [28], [38], [39], [47]. However, to perform detailed transients simulations an accurate calculation of the frequency-dependent parameters of the conductors is required. This chapter summarizes the mathematical formulations of some of the available tech-niques for the theoretical calculation of the frequency-dependent parameters of power cable arrangements. The methods explained are the analytic (classical) method for cylindrically-shaped cables, the approximate method, and the finite element method for arbitrarily-shaped cables. The results of these methodologies are compared in Chapter 6 with the results of the proposed new methodology based on the partial subconductor equivalent circuit method, which is described in detail in Chapter 5. In addition, the topic of earth-return impedance calculations is addressed in the present chapter. Techniques such as the incremental inductance method and the perturbational method are also described in the technical literature [32]. Those two methods apply to cases where the thickness of the conductor is much greater than several skin depths, e.g., microstrip lines at very high frequency, and are not dealt with in the present work. 25 3.2. Basic Assumptions 26 The chapter is organized as follows: Section 3.2 lists the assumptions on which the frequency-dependent parameters of power cable parameters are traditionally calculated. Section 3.3 explains the analytic (classical) method for cylindrically-shaped cables, de-scribing the calculation of the loop impedances, the transformation of the loop quantities into conductor quantities, and the equations reduction process (elimination of the sheaths). Section 3.4 discusses Ametani's approximate method for arbitrarily-shaped cables, present-ing the equations for the impedances and the expressions relating the areas and perimeters of the conductors with those of the equivalent tubes. Section 3.5 summarizes the finite element method, mentioning the basic equations of the solution region, the base functions of the field distribution, the approximation of the basic equations to a residual, and the calculation of the impedances from the field solution. Section 3.6 presents Pollaczek's formula for the calculation of earth-return impedances, describing its solution through numerical integration (Uribe's method) and through closed-form approximations (Wedepohl's method). Finally, Section 3.7 summarizes the content of the chapter. 3.2 Bas ic Assumpt ions In general, the available techniques for the theoretical calculation of frequency-dependent parameters of multi-conductor systems are based on the following assumptions [19], [32], [76]: • The conductors are considered long enough to make the problem two-dimensional, i.e., only the cross section of the system is considered. • T E M propagation is assumed, i.e., magnetic and electric field are perpendicular to each other and occur only in the plane transverse to the direction of propagation. • The conductor cross sectional dimensions are smaller than the wavelength. • The conductors are rectilinear, parallel to each other, and longitudinally homoge-neous, so spiraling effects are not taken into account. 3.3. Analytic Method 27 • The space, all conductors, and insulation have homogeneous permeabilities. The permeabilities are also constant, so non-linear effects are neglected. 3.3 Ana l y t i c M e t h o d This method is described by Dommel in [19] and is based on calculating the impedances of cylindrically-shaped cables either with Bessel functions [63] or formulae [71]. Reference [57]. recommends the use of this methodology with the formulae proposed by Wedepohl and Wilcox in [71]. The problem can be stated as follows: 1) Consider a transmission system consisting of n buried single-circuit cables (Fig. 3.1), in which each single-circuit is made up of a high voltage conductor and a cable sheath, and in which, for simplicity, only three single-circuits a, b, and c have been depicted; 2) Determine the self and mutual series impedances of such a system taking into account the earth return. Calculation Procedure: 1. Calculate the Loop Impedances: The method and the formulae proposed in [19] and [71] are appropriate to describe the electric quantities associated with this cable system. Such a formulation takes into account the skin effect in the conductors, makes use of the traditional earth-return corrections proposed in the literature [10], [15], [54], [71], and allows system studies at frequencies higher than power frequency. By assuming far-end terminals short-circuited, six coupled equations describe the loop quantities associated with the equivalent networks depicted in Figs. 3.2 and 3.3. " dVu/dx " " Z U a Zl2a 0 0 0 0 Ila dV2Jdx Zl2a Z22a 0 Zab 0 Z a c l2a dVlb/dx 0 0 Zn b Z12b 0 0 dV2Jdx 0 Zab Z12b Z22b 0 Zbc h b dV1Jdx 0 0 0 0 Zii„ Zl2c dV2Jdx 0 Z a c 0 Zbc Zl2c Z22c . l2<= (3.1) where 3.3. Analytic Method 28 Figure 3.1: Cable geometry for analytic (classical) method. 3.3. Analytic Method 29 Z core-out„ I sheath. Z sheath-in, a_ V sheath. 2a Z sheath-out„ Z self earth-retur Z core/sheath-insulation I , ) Z sheath-mutual Z sheath/earth-insulation^ 2a^  core„ sheath„ earth x Figure 3.2: Equivalent network for single-circuit a (similar for single-circuits b and c). 3.3. Analytic Method 30 V2 = V s heath, II Icore, 1-2 Isheath "h Icore, ^11 Z c o r e _ o u t -f- Z c o r e y s n e a th—insulation ~r~ Z s neath—in, Z22 Zsheath—out "h Zsneath/earth—insulation "T" Z s e ] f earth—return, Zl2 — — Z s heath-mutual, ^ab Z m u t u a j earth—returnab, Z a c Z m u t u a ) earth—returnac, Zbc — Z m u t u a l earth—returnbc, and a,b,c — subscripts denoting quantities associated with the single-circuits of phases a, b, and c, respectively, V C ore = voltage drop across core with respect to ground, V s heath = voltage drop across sheath with respect to ground, Icore = current flowing through core conductor, Isheath — current flowing through sheath, Zcore-out = internal impedance per unit length of core conductor, calculated from the voltage drop on the outer surface of the core per unit current, when the current returns through the outer conductor. In this case, the outer conductor is the sheath, Zcore/sheath-insuiation = impedance per unit length of insulation between core and sheath, Zsheath-in = internal impedance per unit length of sheath, calculated from the voltage drop on the inner surface of the sheath per unit current, when the current returns through the inner conductor. In this case, the inner conductor is the core, Zsheath-out = internal impedance per unit length of sheath, calculated from the voltage drop on the outer surface of the sheath per unit current, when the current returns through the outer conductor. In this case, the outer conductor is the earth-return path, Zsheath/earth-insuiation = impedance per unit length of insulation between sheath and earth-return path, 3.3. Analytic Method 31 Zsheath-mutuai = mutual impedance per unit length of sheath. In this case, it is the mutual impedance between the inside loop "core/sheath" and the outside loop "sheath/earth" of one single-circuit, Z s e i f earth-retum — self impedance per unit length of the earth-return path, Zmutuai earth-return = mutual impedance per unit length of the earth-return path. In this case, it is the mutual impedance between the outermost loop "sheath/earth" of one single-circuit and the outermost loop "sheath/earth" of another single-circuit. The impedances in Cl/m associated with the insulation are frequency-independent and can be obtained from where Hi and fi2 are the magnetic permeabilities in H / m of the insulation between core and sheath and between sheath and earth, respectively, r c o r e is the radius of the core conductor, rSh-in is the inner radius of the conducting sheath, Tan-out is the outer radius of the conducting sheath, R is the outside radius of the outermost insulation of the cable. The other impedances are frequency-dependent and can be obtained from modified Bessel functions. If the subscript "tube" is used either for core or sheath, or for more layers of concentric conductors (e.g., armors), the impedances of the conductors in Cl/m can be calculated as follows [19], [63], [71] (3.2) (3.3) pm { l o (mr i n )K 1 (mr 0 U t ) + K o (mr i n ) I i (mr 0 U t )} (3.4) e—in 2nrinD pm {l 0 (mr 0 U t )Ki(mr^) + K o ( m r 0 „ t ) I i ( m r i n ) } (3.5) 27rr 0 U t D P (3.6) 'tube—mutual — 2irrinr0UtD 3.3. Analytic Method 32 where D = I 1 ( m r o u t ) K i ( m r m ) - I i ( m r i „ ) K i ( m r o u t ) , m = yj (j to/J,) / p is the reciprocal of the complex penetration depth of conductor in m - 1 , p is the resistivity of conductor in Q, • m, rin is the inner radius of conductor, r o u t is the outer radius of conductor, / i = pr • p0 is the magnetic permeability of conductor in H / m , with pr ^ 1 if the material of the conductor is magnetic, to is the angular frequency in rad/s, I 0 , I i , K Q , and K i are modified Bessel functions with complex arguments (first kind order 0, first kind order 1, second kind order 0, second kind order 1, respectively). The Bessel functions with complex argument can be represented through Bessel func-tions with real argument, and the Bessel functions with real argument can be calculated through polynomial approximations. The EMTP-support routine T U B E [19], [20], for example, uses those polynomial approximations to estimate the parameters with absolute errors lower than 10~ 7 when the argument equals or is lower than 8, and relative errors lower than 3 • 10~ 7 when the argument is greater than 8. Matrix-oriented programs such as M A T L A B ® [40] also carry routines for the evaluation of Bessel functions with complex arguments. For the core conductors depicted in Fig. 3.1, r i n = 0. Since l i (0 ) = 0, (3.5) becomes rj Pcore^core lo(mcore^core) / Q —\ ^core—out — ~Z Z ~, T V"-*''/ The derivation of (3.7) from Maxwell's equations is presented in Appendix A. The derivation of (3.4), (3.5) and (3.6) is given in [63]. Notes: • Alternative formulae for the calculation of the impedances (3.4), (3.5) and (3.6) are listed in Appendix B. They were proposed in [71] and use hyperbolic functions to represent the frequency-dependent nature of the parameters. 3.3. Analytic Method 33 • Formulae for the calculation of the earth-return impedances Z s e i f earth-return and Zmutuai earth-return are presented in Section 3.6. • Appendix C describes the methodology given in [19] to calculate the shunt parame-ters of underground tubular conductors. • If the single-circuits have additional conductors, e.g., armors, add three more coupled equations and include the corresponding impedances Z 2 3 and Z 3 3 in the single-circuits a, b, and c, as indicated in [19]. Move also the impedances Z a b , Z a c , and Z b c since the outermost loops will be the ones "armor/earth" and derive the formulae for the new impedances by analogy utilizing the right electric properties as well as the appropriate radii. 2. Transform the Loop Quantities into Conductor Quantities: To transform the loop quantities into conductor quantities, use the procedure recom-mended in [19] as follows. In (3.1) add row 2 to row 1, add row 4 to row 3, add row 6 to row 5, and replace I i with I c o r e and I 2 with I s h e a t h + Icore- By doing so, it is possible to prove that the system is described through the conductor quantities indicated below <9V c o r e a /dx ' Zcc a Zcs a Z a b Z a b Z a c Z a c 'core a <9V s heatha/da; Z C S a Z S S a Z a b Z a b Z a c Z a c ••-sheatha 9 V c o r e b / dx Z a b Z a b Zcc b Z c s b Zbc Zbc ^•coreb <9V s h eath b/d.T Z a b Z a b Z C s b S^Sb Zbc Zbc ^sheathb <9V c o r e c /dx Z a c Z a c Zbc Zbc Z C C C Z C s c •••corec . dVsheathc /dx Z a c Z a c Zbc Zbc Z C S C Zss c •••sheathe (3.8) where Z u + 2 Z 1 2 + •^ 22, Z c s = Z12 + Z22, Z s s = z 2 2 . The diagonal elements Z c c and Z s s are the self impedances of core and sheath with return through earth, respectively. The off-diagonal elements Z c s , Z a b , Z a c , and Z b c are the mutual impedances between core and sheath of one cable with return through earth, between sheath a and sheath b with return through earth, between sheath a and sheath 3.3. Analytic Method 34 c with return through earth, and between sheath b and sheath c with return through earth, respectively. As a result of the above-mentioned arithmetic operations, the system is represented in nodal form, with currents expressed as conductor currents, and voltages expressed as voltage drops across the conductors with respect to ground. Note: If armors are present, add rows 2 and 3 to row 1 and add row 3 to row 2. Add also rows 5 and 6 to row 4 and add row 6 to row 5. Similarly, add rows 8 and 9 to row 7 and add row 9 to row 8. 3. Equations reduction: By interchanging the corresponding rows and columns in the impedance matrix, move voltage drops across sheaths and current flows through sheaths to the bottom of the vectors of voltages and currents, respectively. Then, let the set of linear equations be split into subsets of core conductor equations and subsets of sheath equations, and the impedance matrix be divided into the core conductor submatrices [Z c c] and [Z c s] and the sheath submatrices [Z s c] and [Z s s ]. Similarly, let the far-end terminals be short-circuited and voltage drops per unit length and current flows be defined for core conductors and sheaths by the vectors [V c] and [Ic] and [V s] and [Is], respectively [Vc]" [Z c c] [Z c s] r tie] l [Vs] . [Z s c] [Z s s] . [Is] . (3.9) Since the vector of voltage drops across the sheaths [V s] = 0, assuming that both termi-nals of each sheath are grounded, the system can be written as — [V c] = [Z r ed] • [Ic], where the reduced matrix [Z r e d ] = [Z c c] - [Z C S ][Z S S ] - 1 [Z S C ] is the result of a Kron's reduction. The elimination is accurate enough provided that the cable section between grounding points is shorter than the wavelength of the maximum frequency under study [19]. At 60 Hz, for instance, such a condition is satisfied and a three-phase system of coaxial cables (six conductors) is reduced to a 3 x 3 system. Using the reduced 3 x 3 matrix the sequence impedance matrix can be obtained postmultiplying and premultiplying by the symmetrical component transformation matrix and by its inverse, respectively. 3.4. Approximate Method 35 3.4 Approx ima te M e t h o d This method is proposed in [4] and approximates the self and mutual impedances of arbitrarily-shaped cables to those of equivalent tubular conductors. Calculation Procedure: Given an arbitrarily-shaped conductor, the method proposes the calculation of its frequency-dependent impedance (Z c ) by Z c = v/Zde 2 + Z h f 2 (3.10) where Z d c is the impedance (resistance) at dc in Q/m, Z h f is the high frequency impedance in fl/m. Z d c can be calculated as Rdc and Z h f as p c / ( 4p c ) , where 4 is the perimeter of the cross section and p c is the complex penetration depth of conductor; p c is equal to T / p c / ( j wp c), where p,c is the permeability of conductor; and p c can be expressed as RdCS, where 5 is the area of conductor. Introducing the above relations into (3.10), the frequency-dependent impedance becomes A similar expression can be derived for a tubular conductor considering that Rdc = p c/[7r(r 2 u t — r2n)] and that Zhf can be obtained from (3.5) when the frequency approaches to infinity. Assessing such a limit, one obtains Z h f = ujpcpc/(27rrout), and the frequency-dependent impedance becomes J u»c(r2mt - r2in) n / i , J WH-cV out ' i ) (o 1 o\ Rdc\ll + D , 2 (3.12) Matching (3.11) and (3.12), one can prove that the radii of the equivalent tubular conductor are given by Tout = ^ (3.13) 3.5. Finite Element Method 36 ( 2TT ) 2 S 7T (3.14) For instance, if the arbitrarily-shaped cable is a pipe-type cable enclosing sector-shaped conductors, the cylindrical equivalents can be studied using the method explained in Sec-tion 3.3 and the formulae for Z m s u i a t i o n , Z p i p e _ i n and Z m u t u a i given in [19]. However, those formulae neglect proximity effects and are only accurate if there is symmetry or if the conductor radii are much smaller than the distances among conductors [19]. With the approximate method the spacing among equivalent tubular conductors is obtained from the geometric mean distances among actual conductors, and the equivalent conductors may overlap if the actual conductors are close to each other. This method was thoroughly studied by Yin in his Ph.D. thesis work [75] and obtains the cable series parameters from the magnetic vector potential distribution and the volume current densities of the conductors. The basic idea of the method is to subdivide the cable system and its surrounding regions into subregions or elements, and to find the current densities of the conductors and the values of the magnetic vector potential at the element vertices (nodes). If desired, additional nodes can be defined on the element sides or inside the elements. To solve the differential equations, this formulation employs numerical methods such as the Galerkin technique. The Galerkin technique approximates the field variables at a generic spatial point with "base" (also called "trial" or "shape") functions whose coeffi-cients become the unknowns of a set of algebraic equations. The method can be applied to arbitrarily-shaped cables and takes into account both skin and proximity effects. The principal equations of the solution region are a two-dimensional diffusion equation and the equations relating the magnetic vector potential with the current density of each conductor [75]. 3.5 F in i te Element M e t h o d 3.6. Earth-Return Impedances 37 The base functions are usually polynomials since they can be differentiated and inte-grated straightforwardly, and their selection depends on the type of shape assigned to the elements. As cables usually display circular geometries, the solution is obtained more efficiently using curved-sided elements (e.g., cylindrical shells or wedge-shaped subregions). Yin proved this by comparing for the same discretization error the results of meshes of straight-sided elements (e.g., rectangles or triangles) with those of meshes of curved-sided elements. The basic equations can be approximated to a residual, and such a residual can be forced to satisfy the integral equation of the Galerkin technique. The result is a set of algebraic complex equations where the unknowns are the field values at the nodes of the mesh and the source current densities of the conductors, and the forcing functions are the currents of the conductors and the known field values at the boundary nodes. The matrix is sparse and banded, and can be factorized using Choleski decomposition algorithms. The impedances can be obtained directly from the source current densities of the conductors with the Js method, or from the power losses and the stored magnetic energy with the loss-energy method. The Js method is usually preferred since it avoids the calculation of all the field variables. The detailed calculation procedure for this method is presented in Appendix D. 3.6 Ea r th -Re t u rn Impedances A cable system can consist of overhead conductors, underground conductors or a combi-nation of both. A case with underground and overhead conductors can be, for instance, a buried pipeline close to an overhead transmission line. On the assumption that the earth resistivity is homogeneous, series-based asymptotic approximations [10], [18], [19] allow the calculation of the infinite integrals (Carson's inte-grals) associated with the earth-return impedances of overhead transmission lines. Alter-natively, Deri, et al. [15], propose formulae to simplify the calculations. Deri's formulae are 3.6. Earth-Return Impedances 38 functions of the complex penetration depth of the earth-return path and of the distances among conductors and their images. The main advantage of Deri's method is that accurate results can be obtained by simple closed-form expressions. Dommel [19], for example, has compared the two methods and the results differ in no more than nine percent. For the calculation of the earth-return impedances of insulated underground cables, or for geometries involving both overhead conductors and insulated underground conductors, Pollaczek [54] derived a formula which models the earth-return path as a semi-infinite, homogeneous, half space extending itself from left to right and from the earth surface downward (Fig. 3.4). This formulation requires the evaluation of Bessel functions and infinite integrals, and obtains the field distribution in earth and air assuming that the conductors are infinitely thin filaments as compared to the earth-return path. Such an assumption is very reasonable inasmuch as the penetration depth of the earth-return path is generally much greater than the cross-sectional section of an underground conductor. Pollaczek's infinite integrals can be evaluated through numerical integration methods as indicated in [46] and [68], or through infinite series and closed-form approximations as proposed in [71]. A rigorous derivation of Pollaczek's formula from Maxwell's equations is given in [75]. Consider the cable system depicted in Fig. 3.4. According to Pollaczek's formulation the mutual earth-return impedance per unit length of such a system can be calculated as[19],[71] p m 2 { K 0 ( m d ) - K 0 ( m / J ) + P } (3.15) mutual earth—ret urn where + 0 0 exp {-(^  +W/3 2 + ™ 2 } exp{j /3x}d/3 (3,16) \/3\ + y/P + m* —oo and d is the distance between cables % and k, d = y/x2 + (hi — hk)2, (for the self earth-return impedance d is replaced with the radius of the outermost insulation layer R), 3.6. Earth-Return Impedances 39 image cable k Figure 3.4: Cable geometry for earth-return impedances. D is the distance between cable i and image of cable k, D = y/x2 + (hi + hk)2, (for the self earth-return impedance Z;;, D is replaced with 2/ij), K 0 is the modified Bessel function, second class, order zero, m = \J(J up)/ p is the reciprocal of the complex penetration depth of the earth return, x is the horizontal distance between cables i and k, hi and hk are the depths at which cables % and k are buried, respectively, \i is the permeability of soil and air, fj, = /i0, p is the soil resistivity, co is the angular frequency. Making a change of variable f3 = U y / a , where yfa. — l / | p | = |m| = y/cop/p, and after some algebraic manipulations, (3.16) becomes + ^ ° e x p { -2hy /ay /u 2 + j l P = 2 / ^ j = }- cos{xy/au}du (3.17) J U + Ju2 + j 0 where h = (hi + hk)/2. 3.6. Earth-Return Impedances 40 The above expression is given in [46] and resembles Carson's integral. The main dif-ference is that the exponent in Carson's integrand is real and non-oscillatory, whereas the exponent in Pollaczek's integrand is complex and oscillatory. As a result, Carson's inte-gral can be approximated with infinite series of rapid convergence, whereas Pollaczeck's integral can only be expressed as involved infinite series of difficult convergence, especially when |md| > 0.25 [71]. Hence, numerical integration methods are nowadays preferred for those cases. Nguyen [46], for example, uses the trapezoidal rule of integration to solve (3.17). Alternatively, Ametani has proposed evaluating Pollaczeck's integral by replacing the term exp{—2hy/a^u2 + j} with the term ex^{—2hy/au} [19]. Thus, Pollaczeck's integral becomes Carson's integral and the infinite series of rapid convergence can be used. This approximation, however, is valid only when |it| 3> 1. Since Pollaczek's integrand is highly oscillatory, difficulties can also arise when it is integrated numerically [47]. Uribe, et al. [68], for instance, have recently analyzed and identified the oscillatory terms of Pollaczek's integrand, and have proposed effective strate-gies for their numerical integration. 3.6.1 Uribe's Methodology The basic idea of Uribe's method is to detect the crosses through zero of the terms in the integrand associated with the oscillations, and then to subdivide the integration in-terval into subintervals whose limits are given by those zero crossings. The subintervals are then split into samples, and the integration is carried out by conventional numerical methods such as the trapezoidal rule or Simpson's method. Uribe also compares the re-sults of the numerical integration with those of Ametani's method and Saad's closed-form approximations [60]. To identify the oscillatory terms and proceed with the numerical integration, the inte-3.6. Earth-Return Impedances 41 grand in (3.17) is transformed into four factors +oo -2j J{F(u)-u+iG{u)}exp{-£F(u)}exp{-j£G{u)}cos{£riu)du (3.18) where F(u) = y V + Vu^+l/V^, G(u) = v7-"2 + \ftF+\/\/2, £ = 2h\fa. — 2h^ujfi/p, V = x/(2h). In the above expression, the first factor in the integrand is a fixed damping envelope independent from the cable system properties, the second factor is an additional damping envelope, and the last two factors are the oscillatory terms. These last three factors depend on the dimensions of the cable system and the electric properties of the earth-return path. The factor exp{— j £ G ( u ) } is identified as a complex function with irregular oscillations, whereas the factor cos(£ r]u) is identified as a real function with regular oscillations. The truncation criterion is suggested using the asymptotic approximation of the damp-ing factor e x p { - £ F ( u ) } . Since F(u) —» u for u > 1, then exp{—£F(u)} —> exp(£u) for u > 1. A semi-infinite integral can be approximated to a finite integral for numerical evalua-tion. For the asymptote associated with the factor exp{—fF(u)} such an approximation can be expressed as H~00 Umax J exp{-£ 'u}<iu ^ J exp{-£u}du (3.19) 0 0 where Umax = A / f (3,20) A = - ln(e r ) (3.21) i and er is the relative error of the numerical approximation. It is possible to prove; for example, that er = e x p ( £ « m a x ) . 3.6. Earth-Return Impedances 42 Uribe has determined empirically that a value of A = 6 (a relative error of 0.25 percent in (3.19)) is appropriate for the numerical calculation of Pollaczeck's integral provided that Umax > 2; otherwise umax = 2 and A = 2£ from (3.20). The zero crossings of the three oscillatory functions can be determined analytically from 0 to umax, and sorted from the lowest to the highest to subdivide the interval [0. umax] into the subintervals [0,ui] (first lobe), [tii,^] (second lobe), . . . , [un,umax] (last lobe); Ui V i = 1,... , n is a zero crossing of Pollazeck's integrand and n is the number of zero crossings between 0 and umax. The subintervals are then split into iV equally spaced values to represent each lobe with at least TV samples. The solution converges uniformly as N —v oo. Uribe presents results for a wide range of geometrical configurations, e.g., depths of burial between 0.2 m and 100 m, spacings between 0.1 and 100 m, and frequencies between 1 Hz and 1 MHz. Therefore, his methodology can be considered appropriate for the solution of Pollazeck's formula. 3.6.2 Wedepohl's Methodology Wedepohl's and Wilcox [71] approximated the earth-return corrections proposed by Pol-laczek to closed-form expressions. The formulae were derived expanding (3.15) into infinite series and taking the most important terms at low frequencies. If R is the radius of the outermost insulation layer of an underground cable (as depicted in Figs. 3.1 and 3.4), the self and mutual impedances of the earth-return path in Cl/m can be obtained from }u)fj, ( f jmR\ 1 4 ) Zself earth-return — i ~ In I — I + - — g 1 1 1 " ' f (3.22) jwp f f ^md\ 1 2 / 1 •^mutual earth-return = | _ ( —2— ) 2 ~ J ™ " f (3.23) The formula for the self impedance yields very accurate results at frequencies for which | m i l | < 0.25 [71]; / i is the magnetic permeability of the earth-return path in H / m , which can be assumed equal to the magnetic permeability in free space (p 0), 7 = 0.577215665 3.6. Earth-Return Impedances 43 (Euler's constant), and m is the reciprocal of the complex penetration depth of the earth-return path; m = y/(jujp)/p in m - 1 , where p is the earth-return resistivity in Q-m, and h is the depth in m at which the cable is buried. This formula is also very accurate if the depth at which the cable is buried is close to 1 m [71]. The formula for the mutual impedance yields very accurate results at frequencies for which |md| < 0.25 [71]; d is the distance between single-circuits a and b for Z a D (see Fig. 3.1 on page 28), the distance between single-circuits a and c for Z a c , and the distance between single-circuits b and c for Z D C ; £ is the sum of the depths of the single-circuits a and b for Z a b , the sum of the depths of the single-circuits a and c for Z a c , and the sum of the depths of the single-circuits b and c for Z b c . This formula is also very accurate if the depths at which the cables are buried are close to 1 m [71]. Unlike Wedepohl's formulae, Saad's formulae keep the first Bessel function K 0 ( m d ) given in (3.15). However, the second Bessel function K 0 ( m D ) is eliminated and the integral is approximated by exponential functions. Uribe shows that Ametani's method performs better than Saad's formulae, and Dommel shows that Wedepohl's formulae perform better than Ametani's method if the conditions |mi?| < 0.25 and \md\ < 0.25 are fulfilled. 3.6.3 Comparison of Methods Pollaczek's and Wedepohl's methods were compared using an example taken from [19] (Srivallipuranandan's example). In that example R = 48.4 mm, p = 100 f2-m, d = 0.3m and hi = hk = 0.75 m. At the highest frequency of interest (1 MHz), \mR\ = 0.0136 < 0.25 and \md\ = 0.0843 < 0.25. Therefore, Wedepohl's formulae are valid. Figs. 3.5 and 3.6 depict the graphs obtained using 100 points per decade. As shown in the plots, both methods produce roughly the same results for self and mutual earth-return impedances. „ The maximum and minimum differences between Wedepohl's method and Pollaczek's method are: self resistance 10.26 percent and 1.32 percent, self reactance 25.01 percent and 10.91 percent, mutual resistance 10.95 percent and 1.32 percent, and mutual reactance 3.6'. Earth-Return Impedances 44 Self Earth-Return Resistance 1 1 1 1 1 1 • 1 i Pol laczek W e d e p o h l li - ' ' / ' I il 1 I 1 1 10 10 10' 10 10 f <Hz> Self Earth-Return Reactance -i-i-riT 1 | 1 | r 10 10 0 10' Pol laczek W e d e p o h l 10 10 10" f <Hz> 10 Figure 3.5: Self earth-return impedance. 3.7. Summary 45 Mutual Earth-Return Resistance 1.51 •— i •— i •— i 1— i •— i 10° 10' 102 103 104 105 10' f <Hz> Mutual Earth-Return Reactance 6 E2 N 0 1 1 1 1 • 1 1 - Pol laczek W e d e p o h l / / / . --/ / 1 / ' / / / ^ ^ ^ ^ 1 I 1 I 1 1 10° 101 102 1 0 3 1 0 4 1 0 5 1 0 6 f <Hz> Figure 3.6: Mutual earth-return impedance. 39.8 percent and 12.64 percent. Maximum and minimum errors occur at 1 MHz and 1 Hz, respectively. 3.7 Summary Alternative methodologies for the calculation of the frequency-dependent parameters of power cable arrangements have been presented. The classical method is appropriate to calculate the parameters of concentric cylindrical conductors. However, an analytic formulation of the field equations and their solution is difficult given arbitrarily-shaped conductors. The approximate method allows the application of the classical method to arbitrarily-shaped cables and pipe-type cables with eccentric, sector-shaped conductors. Nevertheless, 3.7. Summary 46 the equivalent tubular conductors may overlap. Moreover, for the case of pipe-type cables enclosing eccentric conductors the classical method typically represents the conductors as infinitely thin filaments [19]. Such an assumption neglects proximity effects and is only valid when the spacing is large as compared to the conductor radius, which seldom occurs in pipe-type cables. The finite element method is general and allows the analysis of arbitrarily-shaped cables without the approximations mentioned above. Nonetheless, this method requires the modeling of conducting as well as non-conducting regions. Such a modeling is particularly problematic given the infinite nature of the earth-return path. Pollaczek's formula along with Uribe's technique for the numerical integration of Pol-laczek's infinite integral are convenient methodologies for the calculation of earth-return impedances. These techniques cover a wide range of geometrical configurations and are appropriate for the representation of the earth-return path in general-purpose cable-parameter programs. The curves obtained are continuous and can be approximated with rational functions [38]. Wedepohl's formulas also correctly represent the behavior of the earth-return path if the conditions |mi?| < 0.25 and \md\ < 0.25 are fulfilled. The curves are continuous and avoid the numerical integration of infinite integrals. Chapter 4 Geometry Discretization: Proposed Methodology 4.1 In t roduct ion In the previous chapters it has been mentioned that the frequency-dependent parameters of arbitrarily-shaped cable geometries can be calculated with numerical procedures such as finite elements and partial subconductor equivalent circuits. To apply such techniques, the first step is to discretize the system under study, taking the cable geometry and the physics of the problem into consideration. This chapter presents the proposed methodology for cable geometry discretization. The basic idea of the proposed technique is to draw the cable geometry, or scan its photograph, and to use this digital image (pixel map) to automatically identify the spatial coordinates of square-shaped (pixel-shaped) partial subconductors into which the different conductors can be subdivided. Image resolution and penetration depth techniques are then used to determine the number of subconductors needed according to the frequency at which the parameters are calculated. This chapter is organized as follows: Section 4.2 presents the analysis of the modeling problem and describes the application of the fundamental steps of digital image processing to the problem of discretizing cable geometries. Section 4.3 explains the methodology for geometry discretization and discusses the discretization for coordinates and colors, the 47 4.2. Analysis of the Modeling Problem 48 calculation of geometric mean distances, the effect of the resolution on area and geometric mean distance errors, and the proposed methodology for error reduction. Section 4.4 describes the developed algorithm divided into four blocks: 1) processing and segmentation of input image, 2) rediscretization for high penetration depths, 3) rediscretization for low penetration depths, and 4) storage of data. Throughout the chapter, a single-phase coaxial cable illustrates both the proposed methodology and the developed algorithm. Finally, Section 4.5 summarizes the content of the chapter. 4.2 Ana lys is of the Mode l i ng P rob lem Since the starting point is to use an image, the following steps have to be considered to process it. The image has to be acquired, digitized, and subjected to modules of preprocessing, segmentation, representation, description, recognition, and interpretation, under the supervision of a knowledge base [24]. The problem-domain in this case is made up of cross-section drawings of power cable arrangements, and the objective is to determine the frequency-dependent series parameters (resistance and inductance) suitable for electromagnetic transient studies. The image has to be acquired and discretized at spatial coordinates and colors to distinguish between conductors and insulation. To acquire and digitize the image two options can be considered. One option is to trace the cable layout with technical drawing programs capable of discretizing and saving the figure on conventional graphic formats. Another option is to scan a scale drawing and save it using the same aforementioned formats. Regardless of the option, the file stores equally spaced samples of a rectangular matrix in which each element is a discrete color (or gray-level) value. If a scale drawing is initially used, preprocessing is not necessary assuming the user specifies a scale model without noise or alphanumeric information. The segmentation step extracts individual conductors from the background and yields the raw partial subconductor data in which space coordinates and discrete colors of con-ductors and insulation are stored. It is assumed that all distinct conductors (cores, sheaths, 4.2. Analysis of the Modeling Problem 49 armors, etc.) are assigned different colors. A regional representation of the segmentation process, instead of a boundary one, is therefore considered appropriate. The representation step transforms the raw data into a form appropriate for computer analysis, and the description phase converts these data into the quantitative information of interest. Since the image has to be analyzed precisely from a mathematical point of view, lossless graphic formats such as bitmap need to be used to represent it. Considering that current desktop computers commonly have GBs of storage facilities in their hard disks, file sizes are not a problem. Conventional graphic formats also enable the easy retrieval of the variables which describe the image. As a result, simple bidimensional array are sufficient. Rows and columns yield the spatial coordinates of the subconductors (pixels), whereas the matrix elements give the colors. It is unlikely that practical power cables will have more than, perhaps, 50 distinct conductors. Therefore, 256 colors (represented in an integer scale) are enough to describe all possible situations. According to the above points, one byte is more than sufficient to store the color of each subconductor, and the spatial coordinates can be converted into actual dimensions since both resolution (pixels per unit length) and scale factor (drawing length/actual length) are known. The recognition step deals with the process of labeling the different regions. It is necessary to identify conductors (cores, sheaths, armors, etc.), semi-conductive layers, in-sulation, and ground. Once the algorithm recognizes the number of connected components (number of regions with different colors), the process of identifying which color corresponds to what conductor, insulator, and/or ground relies on the user. The interpretation step assigns a meaning to the different regions of the image. For instance, the dc resistance per subconductor and Geometric Mean Distances (GMDs) are calculated in this step. Also, penetration depths are compared with subconductor dimen-sions in order to determine if the resolution used is fine enough to accurately represent skin effects; otherwise the resolution is increased automatically by upsampling pixels. 4.3. Spatial Discretization 50 Normally, a drawing is used to define the cable geometry, but a photograph can be used as well. In that last case, it is necessary to acquire and discretize the image by scanning a picture of the cable system cross-section and preprocess it in order to remove noise, extract conductors and insulation from the background, and, perhaps, smooth the edges and change the colors. Preprocessing does not, however, have to be automatic and can be carried out by the user with the help of one of the many graphics software packages available in the market [2]. 4.3 Spat ia l Discret izat ion Let us make the assumptions indicated in Section 3.2 (page 26) and suppose that the continuous image of a 10-kV power cable cross section, illustrated in Fig. 4.1, is discretized both on spatial coordinates and colors. Let us also assume that the samples are equally spaced, and that the colors of the partial subconductors or pixels are denoted by integer quantities. The discrete image is then represented by a two dimensional matrix with integer elements that can be retrieved from a bitmapped file. The lower the resolution or degree of discernible detail, the more the so-called checkerboard effect becomes visible, and the worse the approximation becomes, as depicted in Fig. 4.2. Each conductor is represented by a different color or gray-level. Insulation is denoted by another distinct color as well (e.g.: white). Therefore, each one of the connected components or regions is identified by subconductors of the same color. These results are summarized through a histogram in which the number of occurrences per color is specified. Fig. 4.3 illustrates, for simplicity, the case of lowest resolution. The coordinate axes are defined following the traditional conventions of Digital Image Processing: rows represent x-coordinates (vertical position) and columns y-coordinates (horizontal position). Assigning labels 0, 1 and 2 to the pixels of insulation, core and sheath, respectively, the image is represented with the integer matrix shown in Fig. 4.4. The integer coordinates of the subconductors are transformed into actual values as 4.3. Spatial Discretization 51 follows x-coordinate rn — 0.5 (4.1) sf res y-coordinate cn — 0.5 (4.2) sf res where rn is the row number of subconductor, cn is the column number of subconductor, sf is the scale factor of the image (drawing length/actual length), res is the resolution of the image in pel/in. To analyze the discretization error, one can compare the actual (analytic) dc resistances and Geometric Mean Distances (GMDs) of the original conductors with those obtained from the dc resistances and the G M D s of the system of discretized subconductors. The discrete (approximate) resistances of the conductors are obtained by performing the par-allel of the subconductors, while the discrete G M D s of the conductors are calculated by bundling the subconductors. core: a =34060 S/mm r. = 10.16 mm sheath:a = 4800 S/mm o=0 roul1=24.38 mm Win= 44.20 mm r i , ,= 40.13 mm (xr = 1 for all conductors and insulators Figure 4.1: Cross section of single-core coaxial cable. 4.3. Spatial Discretization 52 Figure 4.2: Effect of changing the spatial resolution. The dc resistances per unit length and the G M D s of the subconductors are calculated as [27],[34],[48],[72] Ri = ^ T (4-3) OiAi In (GMD i f c ) = J J In (r) dA{dAk (4.4) M Ak where Ri is the dc resistance per unit length in Q/xn of the subconductor i, o~i is the conductivity in S/m of the subconductor i, 4.3. Spatial Discretization 53 Origin 1 2 3 4 5 6 7 8 9 10 11 12 13 I core: 36 subconductors sheath: 16 subconductors Figure 4.3: Coordinate axes for discrete image. 1 2 3 4 5 6 7 8 9 10 11 12 13 1 0 0 0 0 2 0 0 0 2 0 0 0 0 2 0 0 2 0 0 0 0 0 0 0 2 0 0 3 0 2 0 0 0 0 0 0 0 0 0 2 0 4 0 0 0 0 1 1 1 1 1 0 0 0 0 5 2 0 0 1 1 1 1 1 1 1 0 0 2 6 0 0 0 1 1 0 0 0 1 1 0 0 0 7 0 0 0 1 1 0 0 0 1 1 0 0 0 8 0 0 0 1 1 0 0 0 1 1 0 0 0 9 2 0 0 1 1 1 1 1 1 1 0 0 2 10 0 0 0 0 1 1 1 1 1 0 0 0 0 11 0 2 0 0 0 0 0 0 0 0 0 2 0 12 0 0 2 0 0 0 0 0 0 0 2 0 0 13 0 0 0 0 2 0 0 0 2 0 0 0 0 Figure 4.4: Integer matrix for discrete image. 4.3. Spatial Discretization 54 Ai and Ak are the areas in m 2 of the subconductors i and k, respectively, r is the distance in m between an arbitrary point in subconductor i and an arbitrary point in subconductor k. The solution of the integral for square subconductors yields the following results: For the self G M D s GMD,- = 0.447054 (4.5) For the mutual G M D s of neighboring (adjacent) subconductors touching one another (with the same y- or x-coordinate) G M D j f e = 1.00655dtf. For the mutual G M D s of the remaining subconductors i j(x- x'f - 6(x - x'fiv - y'Y + (y- y'Y ln(GMD i f c ) = -24 \ X2 x2 2/2 y'2 J x\ y\ y/Hix - x'Y + {y- y'Y) -( x - x ' f t a n - 1 (1LZJL\ + \x — x'y + ( y - y ' ) 2 t a n - 1 ^ \y - y where is the distance between subconductors i and k, £i is the side of subconductor i, £k is the side of subconductor k, x is the x-coordinate of an arbitrary point in subconductor i, y is the y-coordinate of an arbitrary point in subconductor i, x\ is the upper side x-coordinate of subconductor i, X'i is the lower side x-coordinate of subconductor i, yi is the left-hand side y-coordinate of subconductor i, y2 is the right-hand side y-coordinate of subconductor i, (x - x')(y - y') 25 12 (4.6) (4.7) 4.3. Spatial Discretization 55 x' is the x-coordinate of an arbitrary point in subconductor k, y' is the y-coordinate of an arbitrary point in subconductor k, x[ is the upper side x-coordinate of subconductor k, x'2 is the lower side x-coordinate of subconductor k, y[ is the left-hand side y-coordinate of subconductor k, y'2 is the right-hand side y-coordinate of subconductor k. The mutual G M D s can also be approximated as the distances among subconductors, incurring in a maximum error of 0.655 percent in the calculations (when neighboring subconductors with the same y- or x-coordinate touch one another [7]). The discrete (approximate) values for the dc resistance per conductor, RdC, and the G M D per conductor are given by Rdc = Ri/m (4.8) G M D 1 i=m k=m GMD™ Yl Yl GMDik (4.9) i = l fc=l where m is the number of equal subconductors into which a conductor is subdivided[26]. The actual (analytic) values for tubular conductors can be calculated from Rdc = . 2 l 2-T (4.10) <*nrLt - rin) GMD = r o u t ( r - ^ ] ( r ° u ' i ; " ^ exp { ~ T \ \ ) (4.11) where r o u t and r i n are the inner and outer radii of the tube, a is the conductivity of the conductor, exp (a;) is the exponential function ex, e is the Neper's number. If the conductor has an arbitrarily-shaped geometry, its area can be divided into simpler-shapes [27], e.g., triangles, to obtain its actual G M D from (4.9) as well. Formulas for the 4.3. Spatial Discretization 56 self and mutual G M D s of triangles are proposed in [67], and produce negligible errors in (4.9) if the boundaries of the object are represented with short segments. Figs. 4.5 and 4.6 depict the resistance and G M D errors obtained for the coaxial cable geometry (Fig. 4.1, page 51) when its 200-pel1/in image was rediscretized from 1.8 to 50 pel/in with 0.1-pel/in steps and the nearest neighbor interpolation method [41]. The actual resolution of the sheath is 199.7 pel/in since it is a 3.32-inx3.32-in object represented with a 663x663 matrix. The actual resolution of the core is 200 pel/in since it is a 1.92-inxl.92-in object represented with a 384x384 matrix. Setting one pixel as the minimum limit to represent the thickness of conductor, the initial resolutions for core and sheath were 1.8 pel/in and 12.5 pel/in, as the thicknesses of core and sheath are 0.56 in and 0.08 in, respectively. As shown in the plots, the discretization errors tend to decrease as the resolution grows, but local minima and maxima appear. It was also found that the number of subconductors increases noticeably as the resolution rises. For the core, the number of subconductor varies from 8 (at 1.8 pel/in) to 5,894 (at 50 pel/in). For the sheath, it changes from 114 (at 12.5 pel/in) to 2,011 subconductors (at 50 pel/in). At the maximum resolution (200 pel/in), core and sheath are represented with 95,386 and 32,276 subconductors, respectively. The area errors can be compensated for by dividing the actual area among the subcon-ductors, as proposed in [69], [70]. Similarly, one can think of compensating for the G M D errors by assigning corrected self G M D s to the subconductors. The corrected self G M D s can be obtained from the actual G M D per conductor finding the value of G M D ; in (4.9). However, the inductance depends more on the distances among subconductors (external flux) than on the self G M D s (internal flux). Therefore, such a correction is not effective. In spite of the above-mentioned difficulty, one can analyze the behavior of the G M D error functions to obtain a better estimate for the resolution. Thus, the present thesis proposes that the mutual distances (and therefore the GMDs) can be corrected by an iterative procedure. For instance, keeping the number of subconductors constant, the ' pe^pixel 4.3. Spatial Discretization 57 Figure 4.5: Discretization errors for core. 4.3. Spatial Discretization 58 Figure 4.6: Discretization errors for sheath. 4.3. Spatial Discretization 59 resolution can be varied in discrete steps until the G M D error changes its sign. When that event occurs, linear interpolation can be used to calculate a resolution closer to the crossing through zero of the error function. Fig. 4.7 illustrates the corrected G M D errors after applying the proposed technique. As shown the errors decrease noticeably; for core and sheath, the maximum G M D errors change from —15.6 percent (Fig. 4.5) to —0.029 percent and from —2 percent (Fig. 4.6) to —0.00164 percent, respectively. GMD error for core A V o l_ L_ 0) Q 0 r--0.01 -0.02 -0.03I x 10 20 25 30 35 resolution <pel/in> GMD error for sheath resolution <pel/in> Figure 4.7: G M D errors after corrections. 50 On the other hand, the penetration depth decreases as the frequency steps up, thus increasing the resolution required to have sizes of subconductors with the same order of magnitude as the penetration depth. To reduce memory requirements at high frequencies, another proposal of this thesis is to consider only the boundaries (edges) of the conductors 4.4. Developed Algorithm 60 in such cases. After discretizing the cable geometry with a high resolution, the edges of the conductors can be detected using the magnitude of the gradient associated with the function "pixel value". For example, if the gradient at one pixel is greater than a given threshold, such a pixel belongs to an edge [35]. The edges of the conductors have fewer subconductors, and their resolutions and areas can be corrected in a similar manner as it was described in the above paragraph. 4.4 Developed A lgo r i t hm Once the cable geometry is acquired (either by drawing or scanning) and discretized (by saving its image as a bitmapped file), a rediscretization (RD) algorithm transforms the bitmapped image into the data required by the partial subconductors program. The R D algorithm can be divided into four blocks (Fig. 4.8): processing and segmen-tation of input image, rediscretization for high penetration depths, rediscretization for low penetration depths, and storage of data. As depicted, such blocks are followed sequentially. RD algorithm ("Start) Processing and segmentation of input image Rediscretization for high penetration depths Rediscretization for low penetration depths 1 Storage of data C U D Figure 4.8: Developed algorithm for geometry discretization. The prototype of the R D algorithm was implemented in M A T L A B ® [40] using its 4.4. Developed Algorithm 61 Image Processing Toolbox [41]. 4.4.1 Processing and Segmentation of Input Image This block reads the input data and creates independent images for the different regions (conductors and insulation) of the cable geometry. It can be subdivided into two subblocks (Fig. 4.9): reading of input data and recognition and extraction of conductors. R e a d i n g of input data (input image , d imens ions , etc.) Recogn i t ion and extract ion of conduc tors Figure 4.9: Block of processing and segmentation of input image. The subblock of reading input data reads the bitmapped file where the image sampling and color quantization (color selection) are stored. Data such as input resolution, con-ductor dimensions, scale factor, actual areas of conductors, actual G M D s of conductors, conductivities, permeabilities, frequencies for which the parameters are to be calculated, resolution step, minimum resolution for object representation, and tolerance for G M D error are also entered in this subblock. The subblock of recognition and extraction of conductors splits the cable geometry into conductors, insulation and background, and counts the number of different regions by looking for distinct colors (or grey-levels) in the histogram of the input image. Since the histogram stores the number of occurrences (pixels) per color, the number of different regions is equal to the number of its non-zero entries. Then, by seeking the coordinates associated with pixels of the same color, the input image is segmented into as many images as different regions are present. For example, one image for core (neither sheath nor insulation included), one image for sheath (neither core nor insulation included), and one image for insulation (neither core nor sheath included). 4.4. Developed Algorithm 62 Next, this subblock of the program displays the regions sequentially (using the same colors as the input image) and asks the user to identify each one as conductor or insulation. If the region is a conductor, its label (e.g., core or sheath) is entered by the user, and its actual resolution is calculated as conductor size <pel> divided both by conductor dimen-sion <in> and scale factor. The conductor sizes (width and height) are obtained from the lowest and highest y- and x-coordinates of the subconductors, respectively. Finally, each image of conductor is transformed into a binary image (one bit per pixel) to speed up its processing and downsampling in the next blocks. In the image of the core, for example, pixels are "on" (equal to one) if they belong to the core, otherwise they are "off" (equal to zero). Fig. 4.10 depicts the process of segmenting the coaxial cable of Fig. 4.1 (page 51) into three images. As shown, the objects are not moved from their original positions so that all the conductors use the same coordinate axes. Coaxial cable core sheath insulation + background x x x Figure 4.10: Segmentation of coaxial cable. 4.4. Developed Algorithm 63 4.4.2 Rediscretization for High Penetration Depths This block creates the low resolution images used by the partial subconductors program when the penetration depth is greater than the thickness of conductor. It can be subdi-vided into two subblocks (Fig. 4.11): rediscretization according to thickness of conductor and iteration for G M D error minimization. The subblock of rediscretization according to thickness of conductor can be split into two modules (Fig. 4.12): calculation of initial res-olution per conductor and conductor rediscretization. The subblock of iteration for G M D error minimization can be depicted with three main modules and three iteration loops (Fig. 4.13). The main modules are: G M D calculation, resolution interpolation and con-ductor rediscretization. The iteration loops are: sign change detection (innermost loop), change of number of subconductors (inner loop), and tolerance check (outer loop). Rediscre t iza t ion accord ing to th ickness of conduc tor Iteration for G M D error min imizat ion Figure 4.11: Block of rediscretization for high penetration depths. Calcu la t ion of initial resolut ion per conduc tor r Conduc to r Red iscre t iza t ion Figure 4.12: Subblock of rediscretization according to thickness of conductor. The module of calculation of initial resolution per conductor determines the minimum resolutions that should be used to represent the thickness of each conductor. To obtain such 4.4. Developed Algorithm 64 Resol_c = Initial Resolution Resolution = Initial Resolution *\ Resol c =Resolution r-N Yes d = Change of number of sufjconductors Conductor rediscretization Resolution = Resolution +AResolution inner loopj No Resolution^ Max. Resolution^ !Yes GMD calculation GMD mismatch > 0 1 No i Yes AR = -AResolution AR = AResolution Resol c =Resol c + AR GMD calculation innermost loop: L N o - ^ Change of sign in GMD I Yes Resolution interpolation GMD calculation -No-GMD mismatch < tolerance and Resol c > Initial Resolution outer loopl Resol_c = Resol. for min. GMD mismatch Yes Figure 4.13: Subblock of iteration for G M D error minimization. 4.4. Developed Algorithm 65 values, two criteria are proposed. For rectangular conductors, the thickness of conductor should be represented with at least one pixel. Therefore, the minimum resolution in pel/in can be calculated as 1/thickness. For circular and curved conductors, the thickness of conductor can be represented with the hypotenuse of one pixel. Therefore, the minimum resolution is calculated as v/2/thickness. Of course, the factors of the numerators (1 and \/2) can be changed if higher initial resolutions are desired. For example, one may want to use at least two pixels to represent the thickness of a conductor. In such a case, the factors of the numerators would be 2 and 2\/2 for rectangular and curved conductors, respectively. A factor equal to \/3 or 2 may be applied to curved conductors as well. The module of conductor rediscretization downsamples the binary image using a re-sizing factor rf. After obtaining the initial resolution, rf will be the ratio between the initial resolution and the actual resolution. For example, with rf = 0.05, a new image with a number of rows and columns 20 times smaller is created using nearest neighbor interpolation. If the calculated values of initial resolution are lower than the minimum resolution for object representation specified by the user, the latter is used as value of initial resolution. If the value specified by the user is smaller (e.g., zero), the program uses the calculated values to rediscretize the conductor. In any case, the image obtained is displayed indi-cating its number of subconductors and the user is offered the option of changing the resolution. Thus, new values can be tested until the results look adequate (e.g., object without unexpected asymmetries) and the number of subconductors becomes acceptable (i.e., not too high for low frequencies). The subblock of iteration for G M D error minimization calculates the corrected resolu-tion (ResoLc in Fig. 4.13) that satisfies the constraint imposed by the tolerance for G M D error. Using the initial resolution chosen by the subblock of rediscretization according to thickness of conductor, the iterative process starts calculating the discrete G M D of con-ductor. To estimate it with (4.9), the module of G M D calculation computes first the G M D of subconductor and all the mutual distances among subconductors. Then, it calculates the G M D mismatch as (discrete G M D - actual G M D ) x 100%/actual G M D . 4.4. Developed Algorithm 66 Next, keeping a constant number of subconductors, the innermost loop of iteration is initiated. If the G M D mismatch is positive, then the discrete G M D is greater than the actual G M D and the resolution (res) is increased to decrease both distances and discrete G M D . Conversely, if the G M D mismatch is negative, then the discrete G M D is smaller than the actual G M D and res is decreased to increase both distances and discrete G M D . The changes of resolution are performed with the specified resolution step (typically 0.1 pel/in) until a change of sign in the G M D mismatch is detected. When the change occurs, the module of resolution interpolation calculates a better estimate of the resolution for zero G M D error by interpolating between the two resolution steps bounding the change of sign. Then, the module of G M D calculation is called again to determine the G M D for interpolated resolution and the new G M D mismatch. Later on, the G M D mismatch is compared with the specified tolerance for G M D error, and the interpolated resolution is compared with the initial resolution chosen by the subblock of rediscretization according to thickness of conductor. If both restrictions are satisfied, the program proceeds to the block of rediscretization for low penetration depths. Otherwise the program remains in the outer loop, uses the inner loop to change the number of subconductors, and goes back to the innermost loop of sign change detection to continue the iterative process. Although not depicted in Fig. 4.13, the execution is transferred to the outer and inner loops and the number of subconductors is increased when unlikely events such as reaching a resolution lower than zero or not finding a zero crossing of G M D mismatch occur. If the maximum resolution is reached without satisfying the specified tolerance, the iterative process finishes and the program chooses and displays the resolution for minimum G M D error thus found. For example, specifying minimum resolutions of 10 pel/in and 15 pel/in for core and sheath, respectively, and a tolerance for G M D error equal to 0.1 • 1 0 - 2 percent, the pro-gram converges to corrected resolutions equal to 11.489776 pel/in and 15.088372 pel/in, respectively. In that case, the number of subconductors of the core is 312, the num-ber of subconductors of the sheath is 185, and the G M D errors of core and sheath are —0.700152 • 10~3 percent and —0.453696 • 10~3 percent. If the tolerance is increased ei-4.4. Developed Algorithm 67 ther to 0.1 • 10 _ 1 percent or 0.1 percent, the corrected resolutions for core and sheath are 10.027566 pel/in and 15.088372 pel/in, and the G M D error of core is -0.197687 • 10"2 percent. If the tolerance is reduced to 0.1 • 10~3 percent, the corrected resolutions for core and sheath are 21.396805 pel/in and 21.403158 pel/in, the number of subconductors of the core is 1,091, the number of subconductors of the sheath is 385, and the G M D errors of core and sheath are -0.678624 • 10~4 percent and -0.664774 • 1 0 - 4 percent, respectively. 4.4.3 Rediscretization for Low Penetration Depths This block creates the higher resolution images used by the partial subconductors program when the penetration depth is smaller than the thickness of conductor. Its subblocks and modules are similar to the ones of the block of rediscretization for high penetration depths except for three main differences: instead of the thickness of conductor, the penetration depth is used to determine the initial resolution; the conductor images are rediscretized at each frequency for which the parameters are to be calculated; and the conductors, when the penetration depth is low, are represented by their edges. The block can be subdivided into two subblocks as well (Fig. 4.14): rediscretization according to penetration depth and iteration for G M D error minimization. However, the subblock of rediscretization according to penetration depth displays now three modules inside a counter of frequencies: calculation of initial resolution per conductor, conductor rediscretization and edge detection. The subblock of iteration for G M D error minimiza-tion is the same explained in Section 4.4.2 (depicted with three main modules and three iteration loops in Fig. 4.13). The additional module for edge detection is also shown in Fig. 4.14 and can be subdivided into three submodules: edge determination, edge count and labeling, and edge selection. The module of calculation of initial resolution per conductor determines the minimum resolutions that should be used to have sizes of subconductors with the same order of magnitude as the penetration depth. The criteria are similar to the ones proposed in Subsection 4.4.2, but, instead of the thickness of conductor, the penetration depth is used to calculate the resolution. For rectangular conductors, therefore, the minimum resolution 4.4. Developed Algorithm 68 , ;—No—<(6 j < thickness of conductor^ i = 0 i > nf i - i + 1 nf = number of frequencies Yes Rediscretization according to penetration depth Calculation of initial resolution per conductor Conductor rediscretization 'number of subconductors > n r or -No-Yes Edge determination 1 Edge count and labeling Edge detection Edge selection -4-Iteration for GMD error minimization Figure 4.14: Block of rediscretization for low penetration depths. 4.4. Developed Algorithm 69 in pel/in can be calculated as 1/5, while for circular and curved conductors the minimum resolution in pel/in can be calculated as \/2/5. The penetration depth in the conductor, 5, can be obtained from 5 = ^l / (7r/a/ i ) (4.12) where a and a- are the conductivity and permeability of conductor and / is the frequency in Hz. The minimum resolution is then rounded to start with an initial resolution accurate to one decimal place, i.e., if the minimum resolution is 18.193426 pel/in, the initial resolution is given a value equal to 18.2 pel/in. As previously mentioned in Subsection 4.4.2, different factors for the numerators of the minimum resolution formulae can be used as well. Next, the module of conductor rediscretization downsamples the binary image using the resizing factor. If the number of subconductors associated with the initial resolution is too high (e.g., greater than n m a x = 5,000 subconductors), or if the penetration depth is below a specified threshold (e.g., <5m[n = five percent of the thickness of conductor), the execution is transferred to the module of edge detection, which determines the boundaries of the conductors. Then, the subblock of iteration for G M D error minimization proceeds with the process of conductor rediscretization and resolution correction in a similar manner as it was explained in Section 4.4.2. The restrictions in this case are tolerance for G M D error and minimum resolution according to penetration depth. The module of edge detection creates binary images for the boundaries of the conduc-tors following three steps: First, the submodule of edge determination finds the edges of the conductors using the magnitude of the gradient associated with the function "pixel value". For instance, if a particular point in space has a magnitude higher than a given threshold, such a point belongs to the edge of an object. Once the edges of the conductors have been located, this submodule creates a new binary image in which the pixels are given values equal to one if they belong to an edge. Second, the submodule of edge count and labeling enumerates the existing boundaries by analyzing the adjacency of pixels. "Four-adjacent pixels or four-connected pixels share 4.4. Developed Algorithm 70 an edge (side). Eight-adjacent pixels or eight-connected pixels share an edge (side) or a corner" [35]. If the pixels are eight-connected and have the same values, then the pixels form an eight connected component [24]. In addition, if an eight-connected component is made up of pixels equal to one, then such an eight-connected component represents an edge of conductor. The connected components with pixels equal to one are subsequently labeled to create a new indexed image [41] in which each edge of conductor has a different and consecutive number. The maximum label (number) thus found yields the number of edges. Finally, the module of edge selection offers the user the option of choosing among the detected boundaries. The number of edges of a conductor is typically one (non-hollow, solid conductor) or two (hollow, solid conductor). If there are two edges, the outer edge is usually labeled with one since it is the first boundary found when labeling the connected components. The inner edge is labeled with two since its connected component is found afterward. In such a case, this submodule offers the user the option of choosing between inner edge (number two) and outer edge (number one). For example, to calculate the high frequency loop impedance of a coaxial cable, a user should choose, due to skin and proximity effects, outer edge for core and inner edge for sheath. If the cable is a stranded conductor, more edges will appear since there are holes in its cross section. Under that situation, a user can select among the displayed boundaries the ones to be considered in the calculations. A simple alternative for the user is to choose the outer edge assuming that the conductor is solid, and to assign it the actual area and G M D of the skin layer enclosing the stranded conductor. Tables 4.1 and 4.2 list the resolutions obtained from the rediscretization of the coaxial cable depicted in Fig. 4.1 (page 51) considering a tolerance for G M D error of 0.1 • 10~2 percent and frequencies up to 100 kHz. As shown, the number of subconductors is consid-erably reduced at low frequencies after correcting resolutions and G M D s with the proposed technique.2 For frequencies above 40 kHz, the number of subconductors is too high and 2 The values of G M D error indicated in Table 4.1 for 7.595144-pel/in and 9.703081-pel/in resolutions can be noticed as two small peaks in the graph of corrected G M D error of the core (see Fig. 4.7 on page 59). 4.4. Developed Algorithm 71 it is therefore convenient to switch to a boundary representation of the core. Table 4.3 tabulates the results obtained when its skin layer (outer edge) is given a thickness equal to the penetration depth. As can be observed, the number of subconductors decreases noticeably. For example, at 100 kHz the core is represented with 714 subconductors if a boundary representation is used, whereas at the same frequency it is represented with 41,270 subconductors if a regional representation is employed. Table 4.1: Resolution as a function of frequency for core of coaxial cable FVequency Minimum Corrected Number of GMD <5/thickness resolution resolution Subconductors error <Hz> <pel/in> <pel/in> <%> <%> 1E-06 2.6 7.595144 140 -0.810600E-03 100.00 0.1 2.6 7.595144 140 -0.810600E-03 100.00 1.0 2.6 7.595144 140 -0.810600E-03 100.00 10.0 2.6 7.595144 140 -0.810600E-03 100.00 50.0 3.0 7.595144 140 -0.810600E-03 85.74 60.0 3.3 7.595144 140 -0.810600E-03 78.27 100.0 4.2 7.595144 140 -0.810600E-03 60.63 400.0 8.4 9.703081 216 -0.314144E-03 30.31 700.0 11.1 11.489776 312 -0.700152E-03 22.92 1E+03 13.2 13.283105 422 -0.799748E-03 19.17 4E+03 26.4 26.469469 1,677 -0.303169E-03 9.59 7E+03 34.9 34.998834 2,922 -0.943219E-05 7.25 1E+04 41.7 42.198574 4,242 -0.791129E-05 6.06 4E+04 83.4 - 16,502 - 3.03 7E+04 110.3 - 28,966 - 2.29 1E+05 131.8 - 41,270 - 1.92 Table 4.2: Resolution as a function of frequency for sheath of coaxial cable Frequency Minimum Corrected Number of GMD <5/thickness resolution resolution Subconductors error <Hz> <pel/in> <pel/in> <%> <%> 1E-06 17.7 17.742175 264 -0.774059E-03 100.00 0.1 17.7 17.742175 264 -0.774059E-03 100.00 1.0 17.7 17.742175 264 -0.774059E-03 100.00 10.0 17.7 17.742175 264 -0.774059E-03 100.00 50.0 17.7 17.742175 264 -0.774059E-03 100.00 60.0 17.7 17.742175 264 -0.774059E-03 100.00 100.0 17.7 17.742175 264 -0.774059E-03 100.00 400.0 17.7 17.742175 264 -0.774059E-03 100.00 700.0 17.7 17.742175 264 -0.774059E-03 100.00 1E+03 17.7 17.742175 264 -0.774059E-03 100.00 4E+03 17.7 17.742175 264 -0.774059E-03 100.00 7E+03 17.7 17.742175 264 -0.774059E-03 100.00 1E+04 17.7 17.742175 264 -0.774059E-03 100.00 4E+04 31.3 31.477566 811 -0.175931E-03 56.53 7E+04 41.4 41.665525 1,417 -0.130221E-03 42.73 1E+05 49.5 49.817961 2,011 -0.592936E-04 35.75 4.5. Summary 72 Table 4.3: Resolution as a function of frequency for skin layer of core (coaxial cable). Frequency Minimum Corrected Number of GMD resolution resolution Subconductors error <Hz> <pel/in> <pel/in> <%> 4E+04 83.4 83.629113 451 -0.294933E-04 7E+04 110.3 110.630477 600 -0.173060E-04 1E+05 131.8 131.913958 714 -0.689775E-05 4.4.4 Storage of Data The block of storage of data stores the results obtained in a binary file. Data such as number of conductors, frequencies to be studied, conductivities, permeabilities, resolutions, and integer coordinates of subconductors per resolution and per frequency are stored. If the same resolution is used by different frequencies, repetition of information is avoided. 4.5 Summary This chapter has presented the proposed methodology for the discretization of power cable geometries. Standard graphic formats such as bitmap are used to set up subconductor coordinates of cable configurations, given either their manufacturer drawings or their pictures of actual cross sections. Image resolution, penetration depth, edge detection, as well as area and G M D error minimization techniques are then used to reduce the dimension of the problem and improve the accuracy of the results. At high penetration depths the number of subconductors decreases considerably when the resolution is obtained from the lowest value between thickness of conductor and pene-tration depth. In addition, at low penetration depths, the representation of the conductors by their skin layers reduces the dimension of the problem substantially. The discretization error also diminishes noticeably after using the G M D error minimization technique. The proposed methodology is a novel approach which overcomes the difficulties asso-ciated with the discretization of arbitrarily-shaped power cables. Chapter 5 Calculat ion of Frequency-Dependent Parameters of Power Cables: Proposed Methodology 5.1 In t roduct ion As far as frequency-dependent parameter calculation of power cables is concerned, very useful work has been carried out for the analysis of cylindrically-shaped arrangements [3], [71]. Difficulties arise, nonetheless, when arbitrarily-shaped configurations without closed-form solutions are studied. Typical cases are oval- and sector-shaped cables, cables with concentric neutral conductors, square busbars, and power rails of transit systems, where numerical techniques must be applied to determine the parameters. To include skin and proximity effects into the analysis, one of the alternatives of solution is to represent the cable system with equivalent coupled circuits of partial subconductors [6],[9],[7],[25],[26],[37],[67],[70],[72],[78]. In Chapter 4, it was proposed that the cable geom-etry discretization for arbitrarily-shaped cases can be automatically obtained from digital images [24]. Memory requirements and solution time, however, increase noticeably as the frequency and the number of subconductors rise. The present chapter discusses the details of the proposed Partial Subconductor Equiva-lent Circuit (PSEC) method and its implementation in a general-purpose cable-parameter program. To overcome the memory problems, an algorithm to partition and reduce the 73 5.2. P S E C Method 74 impedance matrices is presented. The proposed algorithm minimizes the interaction be-tween hard disk and C P U and allows the P S E C program to calculate the frequency-dependent parameters of power cable arrangements in reasonable time. This chapter is organized as follows. Section 5.2 presents the P S E C method, describing the derivation of the equations for multi-conductor systems and the modification and reduction procedures of the impedance matrix for a single-phase coaxial cable. Section 5.3 explains the proposed partition methodology, dividing it into three partition stages and one subpartition stage and using memory calculations to exemplify it. Finally, Section 5.4 summarizes the content of the chapter. 5.2 P S E C M e t h o d 5.2.1 M e t h o d F o r m u l a t i o n Let us suppose that an arbitrarily-shaped cable geometry with m conductors is subdivided into n subconductors or filaments, and that a fictitious, lossless and circularly-shaped ring encloses the system under study (Fig. 5.1). r Figure 5.1: Subconductors enclosed by ring. Let us also consider that the subconductors are sufficiently small and that there is mutual coupling among subconductors. Under these conditions, constant current density can be assumed in each subconductor with different subconductors having different current 5.2. P S E C Method 75 densities. For subconductors i and k, for example, with a ring-return path r, the self impedance in fl/m of the loop "subconductor i/ring return" is Z ; i = Ru + jcuLu, and the mutual impedance in Q/m between the loops "conductor z/ring return" and "conductor A;/ring return" is Z;k =]cjLifc. Resorting to the definition of G M D [27], [48], the self and mutual inductances La and Lik in H / m can be calculated as \x ( GMDir2 \ ^ ^ H G M D ^ M D J ( 5 1 ) _ _ (JL . (GMDjrGMDfcr\ U k ~ 2^ i n { G M D r G M D t , ) l 5 ^ j where \i is the permeability of the dielectric between subconductors (equal to 4TT • 10~7 H / m for non-magnetic materials). Assuming arbitrary shapes, distinct positions, and different areas for subconductors i and A;, it is possible to prove that G M D i r = G M D f c r = G M D r = a, where a is the radius of the ring. Therefore, (5.1) and (5.2) become W-Gsk) (5-4) As indicated in Section 4.3, if subconductors i and k are square-shaped filaments, G M D j and G M D f c will both be equal to 0.44705^, where £ is the side of the subconductors. It was also mentioned that in such a case the G M D ^ can be approximated as the distance between subconductors i and k, incurring in a maximum error of 0.655 percent in the calculations. Of course, if more precision is desired, the mutual GMDs of square subconductors can be calculated with (4.6) and (4.7). As proposed in [69] and [70], the resistance per unit length of each subconductor, Ru, can be calculated from Ru = mcRdc (5.5) 5.2. P S E C Method 76 where mc is the number of equal parallel subconductors into which a conductor is subdivided, Rdc is the actual dc resistance of the conductor. To simplify the explanation, let us suppose that the number of conductors, m, is equal to two, e.g., core conductor and conducting sheath, and that the far-end terminals of the conductors are short-circuited. The impedance matrix will thus relate the voltage drops across the conductors with the currents flowing through themselves by dVc/dx _ dVs/dx where x is the direction of propagation, V c and V s are the voltage drops per unit length in V / m across core and sheath, respec-tively, I c and I s are the currents flowing through core and sheath, respectively, Z c c is the self impedance in fi/m of the loop "core/ring return", Z c s is the mutual impedance in Q / m between the loops "core/ring return" and "sheath/ring return", Z s s is the self impedance in Q / m of the loop "sheath/ring return". Assuming, for example, that each conductor is arbitrarily subdivided into three parallel subconductors (mc — 3), the system can be rewritten as follows " dVc/dx ' dVc/dx dVc/dx dVs/dx dVjdx _ dVjdx _ where I c is equal to I C l + I C 2 + I C 3 , I s is equal to I S l + I S 2 + I S 3 , Z i ; and Z i k V i = 1, 2, . . . , 6 and V k = 1,2, . . . , 6 are the self and mutual impedances of the subconductors with ring return. Z C c ^cs Ic Zcs ^ss Is (5.6) " z „ Zl2 Zl3 Z14 Zl5 Zl6 lei Z12 Z22 Z23 Z24 Z25 Z26 lc2 Zl3 Z23 Z33 Z34 Z35 Z36 c^3 z 1 4 Z 2 4 Z34 Z44 Z45 Z46 Isi Zl5 Z25 Z35 Z45 Z55 Z56 is2 Zl6 Z26 Z36 Z46 Z56 Z66 (5.7) 5.2. P S E C Method 77 Since the subconductors of each conductor are in parallel, the cable parameters can be obtained by bundling the subconductors in a similar manner as it is done for overhead transmission lines in [29]. 5.2.2 Equations Modification All voltage drops and current flows per subconductor, except the ones of the first sub-conductor of each conductor in (5.7), need to be moved to the bottom of the vectors of voltages and currents, as shown in (5.8). To do so, rows and columns are interchanged in (5.7) as follows: row 4 and column 4 are moved to row 2 and column 2, row 2 and column 2 are moved to row 3 and column 3, and row 3 and column 3 are moved to row 4 and column 4, respectively. The system then becomes (5.8) where A X 1 = Z n , A12 = Z14, and the other elements of [A] can be obtained from the row and column operations described above. Replacing I C l with I c and I S l with I s in (5.8), errors equal to A i i ( I C 2 + I C 3 ) and Ai2(I S 2 + I s 3 ) are introduced into the equations. However, they can be compensated for if column 1 is subtracted both from column 3 and column 4, and column 2 is subtracted both from column 5 and column 6, respectively. To keep the symmetry, the same arithmetic operations applied to the columns are also applied to the rows, i.e., row 1 is subtracted both from row 3 and row 4, and row 2 is subtracted both from row 5 and row 6, respectively. The system thus becomes dVJdx " " A n A12 A i a A 1 4 A15 A i e dVs/dx A12 A 2 2 A 2 3 A 2 4 A 2 5 A 2 6 dVJdx Aj.3 A 2 3 A 3 3 A 3 4 A 3 5 A 3 6 dVjdx A14 A 2 4 A 3 4 A44 A 4 5 A 4 6 dVs/dx A 1 5 A 2 5 A 3 5 A 4 5 A55 A 5 6 dVJdx . A 1 6 A 2 6 A 3 6 A 4 6 A 5 6 A 6 6 •*-C3 . S3 . " dWjdx ' " B 1 X B12 B13 B 1 4 B 1 5 Bi6 dVs/dx B12 B 2 2 B23 B 2 4 B25 B 2 e 0 B13 B23 B33 B34 B35 B 3 e 0 B14 B24 B34 B44 B45 B46 0 B15 B25 B35 B45 B 5 5 B56 0 . B i 6 B 2 6 B36 B46 B56 B66 "Ic Is Ic 2 Ics I S 2 (5.9) 5.2. P S E C Method 78 where all the equations at the bottom are equal to zero and subsets of equations have been defined for the step of equations reduction. For example, B n = A n , B i 2 = A 1 2 , and B 1 3 = A i 3 - A n . Similarly, the other elements of [B] can be determined with row and column operations. 5.2.3 Equations Reduction Specifying subvectors of voltages and currents and submatrices of impedances for the subsets of equations indicated in (5.9), the system can be rewritten in the following way [V] [Zi] [z 2] [0] [z 3] [z 4] [I] [K] (5.10) where - [ V ] = [Zx]-[I] + [Z 2 ] . [K] [0] = [Z3] • [I] + [Z4] • [K] Applying Kron's reduction, (5.11) becomes - [ V ] = [Z r e d ] • [I] where [Z r e d ] - [Zx] - [Z2] [ Z 4 ] _ 1 [Z3] (5.11) (5.12) (5.13) (5.14) and [Z r ed] stores the frequency-dependent parameters of the cable system under study ( Z c c , Z c s , and Z s s in (5.6)). Using partial Gaussian triangularization [7], [19], the off-diagonal upper elements of the submatrix [Z4] and all the elements of the submatrix [Z2] can be transformed into zeroes. Thus, the matrix [Z r ea] becomes equal to the submatrix matrix [Zx] and the cable parameters with ring return can be calculated without inverting the submatrix [Z 4]. 5.3. Partition Methodology 79 Next, the fictitious ring-return path is eliminated introducing its zero current condition into the equations. This is done in the foregoing example replacing I s with - I c in (5.6). Hence, the second column and the second row can be subtracted from the first column and the first row, respectively, to convert (5.13) into -dVcs/8x = Z c c ' • I c (5.15) where Z c c ' = Z c c — 2 Z C S + Z s s is the loop impedance between core and sheath, V c s = V c - V s is the voltage across the loop core-sheath, I c is the current flowing through core conductor and returning through conducting sheath. Finally, the frequency dependent parameters Rcc' and Lcc' of the loop core-sheath are calculated from Rcc' = Re{Z c c'} and Lcc' = Im{Z c c'}/w. 5.3 Pa r t i t i on Methodo logy Given the symmetry of the impedance matrix, it is sufficient to create its elements con-sidering that the matrix is triangular, i.e., it is sufficient to calculate, operate and store a number of elements, ne, equal to n(n + l)/2. For example, if a given cable system is subdivided into 5,000 subconductors (n = 5,000), a computer requires, assuming double precision per complex element (16 bytes), a physical memory of 190.8 M B to construct the impedance matrix without using virtual memory. In such a case ne is close to 12.5 million complex elements and the number of complex operations needed to triangularize the impedance matrix is approximately 20.83 billion, i.e., n 3 /6 after exploiting the symmetry. Common desktop PCs carry relatively inexpensive hard disks of high capacity so that data storage is not a problem. Nor is speed a problem since CPUs becomes much faster every year. Large amounts of physical memory can be included as well, but physical memory is much more expensive. To allow the P S E C program to find a solution in reasonable time when large cases 5.3. Partition Methodology 80 are studied, we propose to partition the system into smaller subsystems that can be tem-porarily stored on a hard disk. The subsystems can be analyzed sequentially, thus avoiding the use of virtual memory when the physical memory of the computer is lower than that required by the impedance matrix. The proposed partition methodology proceeds in three stages: partition for matrix construction, partition for matrix modification, and partition for matrix reduction. Addi-tionally, each partition stage can be subdivided into subpartitions for storage and retrieval of elements. If M is the physical memory of the computer minus the memory used by the operating system, auxiliary arrays and variables, and N is the memory needed by the impedance matrix ( M and iV are both in bytes), partition of the impedance matrix is required when N is greater than M. The prototype of the P S E C program was created using DIGITAL Visual Fortran Ver-sion 5.0 [16] and Microsoft Visual Studio 97 [43]. 5.3.1 Partition for Matrix Construction Once the cable geometry has been discretized and the subconductor coordinates have been obtained, the program can begin the stage of partition for matrix construction. The number of partitions for this stage, np\, is equal to ceil(N/M), where ceil, like in MATLAB,® is the function rounding the ratio N/M to the nearest greater integer, and N is equal to 16n e bytes. The maximum number of elements per partition, n e p i , is given by floor(M/16), where floor is the function rounding the ratio M/16 to the nearest lower integer. If the ratio N/M is not an integer number, the number of elements of the last partition, nepif, is lower than n e p l . In such a case, n e p i / is equal to {n — (n p l — l )n e p l } . The algorithm of this stage is as follows: 1. Determine M. 2. Read n, all coordinates, self G M D s and resistances per subconductor. 5.3. Partition Methodology 81 3. Calculate ne, N, n p l and nep\. 4. Initialize counters and indices, i.e., ip = 1, istart — 1> Und — nep\. 5. Calculate elements of the impedance matrix for partition ip from index istart to index lend-6. Store elements of the impedance matrix for partition ip from index i s t a r t to index lend-7. Increase counter of partitions, i.e., ip = ip + 1. 8. Assign following lower index, i.e., istart = iend + 1-9. Assign following upper index, i.e., iend — ip * nev\. 10. Ask if ien<i > ne. If that is the case, assign correct upper index, i.e., ien(i — ne. 11. Ask if ip < np\. If that is the case, go to step 5. 12. End. Using the example of the foregoing section where n = 6 subconductors and N = 336 bytes (21 elements, 16 bytes each), and assuming an arbitrary value of physical memory, e.g., M = 288 bytes, the results are: np\ = 2 partitions, nepi = 18 elements, and nepif = 3 elements. When ip = 1, the elements are calculated and stored in a working vector using the following order: A n , A i 2 , A 2 2 , A i 3 , A 2 3 , A 3 3 , . . . , up to A 3 6 . Similarly, when ip = 2, the remaining elements (A 4 6 , A 5 6 and A66) are processed. 5.3.2 Partition for Matrix Modification Once all the values of the impedance matrix have been calculated and stored, they can be retrieved from hard disk to proceed to the stage of equations modification. The first rows of the matrix (the ones associated with the first subconductor of each conductor) have to 5.3. Partition Methodology 82 be subtracted from all the other rows, so it is convenient to keep the first rows in physical memory as the process is carried out. Assuming partitions with an integer number of rows, the number of partitions for this stage, nP2, is equal to ceil{(n — m)/rp}; the maximum number of rows per partition, rp, is equal to floor{(M — M/>)/(16n)}; and the memory used by the first rows of the impedance matrix, Mfr, is equal to 16mn bytes. The maximum number of elements per partition, nep2, is equal to rpn. Similar results can be obtained using np2 equal to ceil{(N — Mfr)/(M — Mfr)}, N equal to 16n2 bytes, nep2 equal to floor{(M — M/ r )/16}, and rp equal to floor{nep2/n}. Note that a factor equal to 8(m 2 —m) bytes can be subtracted from Mjr if the symmetry is taken into account. However, it is not worth storing and manipulating the first rows of the impedance matrix in a unidimensional array since the savings obtained are negligible. For example, if m is equal to 50 conductors and n is equal to 100,000 partial subconductors, then Mjr is equal to 76.3 M B and 16n is equal to 1.5 M B , where 16n is the memory required by one row. In such a case, the factor 8(m 2 — m) is equal to 0.02 M B , representing just a saving of 0.03 percent. As recent PCs can carry physical memory greater than 128 M B , it is unlikely to run out of the resources needed to complete the stage of matrix modification for keeping the first rows of a large system in memory. For example, when N is equal to 128, 256, or 512 M B , the case with 50 conductors and 100,000 partial subconductors can still be solved with 3,029, 855 or 351 partitions, respectively. The rows per partition are 33, 117 or 285, respectively. The algorithm of this stage is as follows: 1. Retrieve elements of first m rows of the impedance matrix. 2. Calculate rp. 3. Initialize counters and indices, i.e., ip = 1, i s t a r t = m + 1, ie„d = istart + rp — I, ir istart j k — !• 5.3. Partition Methodology 83 4. Retrieve elements of the impedance matrix for partition ip from row i3tart to row iend. 5. Subtract row k from row ir. 6. Update available elements of column ir using elements of row ir, i.e., A i j r = A j r ; V i = 1, . . . , m and V i = istart, • • • , iend-7. Subtract A; r k from diagonal element A i P i P . 8. Increase counter of rows, i.e., ir = ir + 1. 9. Verify if row ir is still associated with row k. If not, increase counter of first rows, i.e., k = k + 1. 10. Ask if ir < iend. If that is the case, then go to step 5. 11. Store elements of the impedance matrix for partition ip from row i s t a r to row iend. 12. Increase counter of partitions, i.e., ip = ip + 1. 13. Assign following lower index, i.e., istart = iend + 1-14. Assign following upper index, i.e., iend = istart + rp — l. 15. Ask if iend > n. If that is the case, assign correct upper index, i.e., iend = n. 16. Ask if istart < n - If that is the case, then go to step 4. 17. Store elements of first m rows of impedance matrix. 18. End. The example with n = 6 subconductors is also used to explain the partition for matrix modification. Since the matrix is processed by strips in this stage, N is equal to 576 bytes (36 elements, 16 bytes each). Assuming the previous value of physical memory ( M = 288 bytes), the results are: n p 2 = 4 partitions and rp = 1 row. At the beginning, the elements of the first two rows are retrieved. Second, when ip = 1, the elements of the row 3 are retrieved to perform the operation row 3 minus row 1. Then 5.3. Partition Methodology 84 the available elements of column 3 are updated using the elements of row 3, i.e., A i 3 = A 3 1 and A 2 3 = A 3 2 . Next the diagonal element is updated, i.e., A 3 3 = A 3 3 — A 3 1 , and finally the elements of the row 3 are stored. Similar processes are carried out for ip = 2, ip = 3 and ip = 4 using rows 4, 5 and 6, respectively. At the end, the elements of the first two rows are stored. 5.3.3 Partition for Matrix Reduction Once the impedance matrix has been modified and stored, it can be retrieved from hard disk to proceed to the stage of equations reduction (bundling). Since the idea is to trans-form the off-diagonal upper elements of the right (n — m) columns into zeroes, it is conve-nient to keep their multiplication factors in physical memory and to retrieve the remaining elements of the matrix from hard disk as they are needed. First, the elements of the nih column (the last on the right) are retrieved so that the multiplication factors associated with the first (n — 1) rows can be calculated. Since the upper (n — 1) elements of the nih column will become zeroes, (n — 1) additions can be avoided. Then, the program proceeds with the retrieval of the first (ne — n) elements of the matrix. As M < N, such elements are divided into chunks of data using a maximum number of elements per partition, nep3, equal to floor{(ne — n)/n P 3} and a number of partitions, n p 3 , equal to ceil{16(ne — n)/(M — 16n)}. At the beginning, the first n e p 3 elements are retrieved, operated and stored using the corresponding multiplication factors. Next, the second group with n e p 3 elements, and so forth until completing the nP3 partitions. The multiplication factors are different for each row and are modified while traversing the first (n — 1) rows. Second, the elements of the (n — l)th column are retrieved to obtain the new multi-plication factors. Note that the nlh column is filled with zeroes, so it is only required to calculate (n — 1) factors and to retrieve and operate the first {ne — n — (n — 1)} el-ements. As a result, the partitions can have one more element without exceeding the available physical memory, and the number of partitions decreases to a value equal to ceil{16(ne - n — (n - 1)) / (M - 16(n - 1))}. For the (n - 2)th column two more elements 5.3. Partition Methodology : 85 are allowed in each partition, for the (n-3)*'1 column three more elements are allowed, and so forth until the sum of all the elements and multiplication factors fits in memory. When that is the case, the Gaussian elimination continues without partitioning the matrix. The algorithm of this stage can be summarized as follows: 1. Calculate ne, nP3 and nep3. 2. Initialize index of columns and counter of partitions, i.e., i = n, ip = 1. 3. Initialize indices of the impedance matrix, i.e., istart = 1, iend = tiep3-4. Retrieve elements of the ith column of the impedance matrix. 5. Calculate the multiplication factors F i up to Fj . 6. Retrieve elements of the impedance matrix for partition ip from element i s t a r t to element iend-7. Operate elements of the impedance matrix for partition ip from element i s t a r t to element iend using the multiplication factors associated with the rows of the elements and the corresponding elements of the ith column. 8. Store elements of the impedance matrix for partition ip from element istart to element 9. Increase counter of partitions, i.e., ip — ip + 1. 10. Assign next lower index, i.e., istart — iend + 1-11. Assign next upper index, i.e., iend = istart + neP3 - 1-12. Ask if iend > ne — i. If that is the case, assign correct upper index, i.e., iend — ne — i. 13. Ask if istart < ne — i. If that is the case, then go to step 6. 14. Update index of columns, i.e., i = i — 1. 15. Calculate new number of elements, i.e., ne = i(i + l)/2. 5.3. Partition Methodology 86 16. Check physical memory constraint, i.e., ask if 16ne > M. If that is the case, then cal-culate the new number of partitions, i.e., nP3 = ceil{16(ne — i)/(M — 16i)}, calculate the new maximum number of elements per partitions, i.e., n e p 3 = floor {(ne — i)/nP3], and go to step 3. 17. Retrieve elements of the impedance matrix from element 1 to element ne. 18. Continue Gaussian elimination without partitioning impedance matrix. 19. Store elements of the impedance matrix from element 1 to element m(m + l)/2. 20. End. The partition process for matrix reduction can also be visualized for the example of six subconductors (7Y = 336 bytes). Assuming that M is equal to 160 bytes, the results are: n p 3 = 4 partitions and n e p 3 = 4 elements. Initially, the elements of the 6th column are retrieved to create the array of multiplica-tion factors. Then, groups of four elements are sequentially retrieved and operated. The first group consists of the elements B u , B i 2 , B 2 2 and B13, the second group of the ele-ments B 2 3 , B 3 3 , B 1 4 and B 2 4 , and so forth until reaching the fourth group and the remain-ing elements ( B 3 5 , B 4 5 and B 5 5 ) . For example, B X 1 will become equal to B n -I- F i B i 6 , B i 2 will become equal to B x 2 + F j B ^ , B 2 2 will become equal to B 2 2 + F 2 B 2 6 , B i 3 will become equal to B i 3 + F i B 3 6 , and so forth until using up all the partitions, where Fx = — B 1 6 / B 6 6 is the multiplication factor associated with the elements of the first row and F 2 = — B 2 6 / B 6 6 is the multiplication factor associated with the elements of the sec-ond row. Likewise, F 3 , F 4 and F 5 can be calculated. F 6 will store the ratio — 1 / B 6 6 , which is the value used to avoid more division operations when calculating the factors F i up to F5. Once the calculations of each partition have been carried out, the corresponding results are stored. Later on, the elements of the 5th column are retrieved to create a new array of multipli-cation factors. Similarly, groups of five elements are sequentially retrieved and operated, 5.3. Partition Methodology 87 but the number of partitions is reduced to three. For example, the first group consists of the elements B n , B i 2 , B 2 2 , B 1 3 and B 2 3 . Finally, when the 4th column is reached, the 10 remaining elements fit in memory and no more partitions are required to continue the partial Gaussian elimination process. 5.3.4 Subpartition for Storage and Retrieval of Elements As indicated, the methodology requires storage and retrieval of the chunk of data associated with each partition if M < N. The ideal situation would be to store and retrieve the partition at once so that only one hard disk access per partition is needed. Using single access, however, the time values increase noticeably as the amount of data to be written or retrieved grows. Table 5.1 lists the timings obtained for different record sizes using an A M D - K 6 , ™ 4 3 3 -MHz computer with 60 M B of physical memory. As shown, it is generally faster to write a long record of binary data to a direct access file than to write it to a sequential access file. Also, the process slows down when the record is too long. Therefore, one can conclude that writing 50 M B of data with five 10-MB records is faster (approximately 0.17 sx5 = 0.85 s) than writing 50 M B with only one record (160.49 s). Table 5.1: Hard disk time versus record length. Record Length <MB> Writing Time <s> Sequential Access Direct Access 1.0 0.54 0.05 2.0 0.50 0.50 5.0 0.38 0.22 10.0 0.66 0.17 15.0 2.20 1.70 20.0 2.92 2.31 30.0 96.45 3.35 40.0 246.07 4.56 50.0 369.43 160.49 The delays also vary due to factors such as speed of computer and type of hard disk, and exceed acceptable limits if each record has to be distributed among many physical locations on the disk. Given any arbitrary hardware platform, however, delays out of reasonable limits can be avoided if the P S E C program determines an appropriate size of 5.4. Summary 88 subpartition (record) before starting the construction of the impedance matrix. If a higher speed is needed, the storage and retrieval of data can be carried out incorporating assembly language routines into the code. 5.4 Summary A methodology to calculate frequency-dependent parameters of arbitrarily-shaped power cable arrangements based on the partial subconductor equivalent circuit method has been presented. The proposed methodology partitions the partial subconductor impedance matrices when their sizes exceed the available physical memory, thus allowing the P S E C program to estimate the parameters in low R A M computers. The algorithm also constructs, modifies and reduces the impedance matrix exploiting its symmetry, thereby reducing the number of operations for the stage of matrix reduction from n 3 to n 3 / 6 [7], [19]. The reduction of the number of arithmetic operations in that stage is particularly important given the fact that the triangularization process can take more than 96 percent of the total simulation time when partitions are not required, and up to more than 99 percent of the total simulation time when the number of partitions for the stage of matrix reduction is high. For a large number of subconductors, therefore, it is important to consider the trade-off between physical memory and simulation time. If the R A M is low, for example, the number of partitions for matrix reduction will be high and the simulation time will increase considerably. However, the proposed solution method will be able to determine the cable parameters. The developed algorithm can be used to calculate the series parameters of oval- and sector-shaped cables, L-shaped cables, cables with concentric neutral conductors, square busbars, and power rails of transit systems. Chapter 6 Case Studies 6.1 In t roduct ion In the preceding chapters we have explained methodologies for the calculation of the frequency-dependent parameters of power cables. The objectives of this chapter are to study typical cable arrangements with the developed techniques and to compare the results with those of the alternative methodologies. The cases under study are a coaxial cable, a buried cable system and a sector-shaped cable. First, Section 6.2 analyzes the single-phase coaxial cable. The results of this case allow validation of the proposed algorithms (RD and PSEC) since concentric coaxial cables can also be studied with the analytic method. The cable parameters are obtained disabling and enabling the proposed G M D error reduction technique, using fixed resolutions and res-olutions adapted to the thickness and penetration depth of the conductors, and employing regional and boundary representations of the cables at low penetration depths. Second, Section 6.3 studies a system of buried cables. This section also uses the coaxial cables of Section 6.2 and compares the results obtained with Wedepohl's formulae for the earth-return path with those obtained from the numerical integration of Pollaczek's semi-infinite integrals. Next, Section 6.4 studies a three-phase sector-shaped cable and compares the results obtained with the P S E C method with those obtained from analytical solutions, the Finite 89 6.2. Coaxial Cable 90 Element (FE) method and approximate formulae. The results of this case are evaluated with fixed resolutions and without the G M D error reduction technique. Finally, Section 6.5 summarizes the content of the chapter. 6.2 Coax ia l Cable Fig. 4.1 on page 51 depicts the geometry of a coaxial cable. Its parameters, which can be obtained analytically through Bessel functions, are used to validate the proposed method-ologies. The area of the core conductor is 1,543.63 m m 2 and the area of the conducting sheath is 525.35 mm 2 . The sheath is considered the return path and the ground is not included in the calculations. Three different simulations were performed for a frequency range of 0-100 kHz: • Simulation a: regional representation of conductors, unchanging resolution for all frequencies, and G M D error minimization option of the R D program disabled. • Simulation b: regional representation of conductors, unchanging resolution for all frequencies, and G M D error minimization option of the R D program enabled. • Simulation c: boundary representation of conductors at low penetration depths, variable resolution (adapted to the thickness of conductor and the penetration depth per frequency), and G M D error minimization option of the R D program enabled. In all the cases mentioned above, the rediscretization was carried out with the developed R D program and the parameters were calculated with the developed P S E C program. The analytic solution was obtained using the EMTP-support routine T U B E [19], [20]. 6.2.1 Simulation a Applying segmentation, the cable geometry can be divided into two images, one for the core conductor and the other for the conducting sheath. Using square-shaped subcon-ductors, looking for the minimum area errors, and setting a limit of approximately 5,000 6.2. Coaxial Cable 91 subconductors, subconductor densities (resolutions) of 42.9 subconductors/in and 29.9 sub-conductors/in are obtained for core and sheath, respectively. These resolutions yield area discretization errors of 1.51 percent and 2.20 percent and G M D errors of -0.69 percent and —1.17 percent for core and sheath, respectively. With the above resolutions, the core is subdivided into 4,337 subconductors, the sheath is subdivided into 712 subconductors, and the cable geometry is thus represented with a total of 5,049 partial subconductors. As indicated in Chapter 4, the subconductors must have sizes with the same order of magnitude as the penetration depth (S) to take skin effects into consideration [19]. Given the existing circular symmetry, the hypotenuse of every square subconductor should be at least equal to <5, where 5 can be obtained from (4.12), page 69, as a function of the frequency and the conductivity and permeability of the conductor. From the chosen resolutions, the lengths of the hypotenuses of the square subconductors of core and sheath are 0.84 mm and 1.2 mm, respectively, so the cut-off frequencies [37] associated with the internal resistances should be close to 10.61 kHz for the core and 36.56 kHz for the sheath. Fig. 6.1 compares the loop parameters obtained with the two methods using a scale proportional to y/J along the abscissae. Four points per decade (e.g., 10 kHz, 40 kHz, 70kHz, and 100 kHz) are considered up to a frequency of 100 kHz. As shown in the plots, both models predict the right behavior of the frequency-dependent parameters. The inductance decreases since less magnetic flux linkages exists inside the conductors. The loop inductance, therefore, tends to a constant value given by zero internal inductances Lcore-out and LSheath-in, and the frequency-independent induc-tance Lcore_sheathjinsulation [19], which, in this case, represents the external inductance due to the magnetic flux outside the conductors and inside the insulation. The loop resistance increases since currents tend to flow through the outermost and innermost layers of core and sheath, respectively. As a result, less cross-sectional area is available for currents to flow through. 6.2. Coaxial Cable 92 loop resistance T 1 1 r — PSEC 0 1.6 6.4 14.4 25.6 40.0 57.6 78.4 102.4 f <kHz> loop inductance 0 1.6 6.4 14.4 25.6 40.0 57.6 78.4 102.4 f <kHz> Figure 6.1: Resistance and inductance of coaxial cable (Simulation a). 6.2. Coaxial Cable 93 The error in the loop resistance is positive (calculated resistance greater than actual resistance) at 10 KHz, and becomes negative (smaller calculated resistance) at 40 KHz. The errors are 0.12 percent and -1.92 percent at 10 kHz and 40 kHz, frequencies which bound the 10.61-kHz and 36.56-kHz values predicted for the cut-off frequencies of core and sheath. The maximum and minimum differences between the values of the P S E C method and those of the analytical method are —6.49 percent and 0 percent on the resistance curve, occurring at 100 kHz and 0.1 Hz, and —1.11 percent and —0.78 percent on the inductance curve, taking place at 10 kHz and 0.1 Hz, respectively. The average simulation time was 59.2 min per frequency using an Intel P e n t i u m ™ I I I , 733-MHz computer with 256 M B of physical memory, three partitions for matrix modifi-cation, and 10-MB subpartitions (records). Setting limits of 192 M B for physical memory and 1 M B for each subpartition, the simulation time was l h 38 r a 32 s per frequency, with two partitions for matrix construction, three partitions for matrix modification, and a total of 222 partitions for matrix reduction (74 for multiplication factors and 148 for the rest of the elements). In addition, a study with 128 M B of physical memory and 1-MB subpartitions resulted in a simulation time of 8 h10m48 s per frequency, with two partitions for matrix construction, four partitions for matrix modification, and a total of 2,979 partitions for matrix reduction (993 for multiplication factors and 1,986 for the rest of the elements). 6.2.2 Simulation b Fig. 6.2 depicts the results of the proposed method after using the G M D error minimization technique described in Section 4.3 (page 56) with 35.9 pel/in and 49.7 pel/in as minimum resolutions for core and sheath, respectively. The corrected resolutions were 35.907699 pel/in (3,077 subconductors) and 49.720173 pel/in (2,013 subconductors), resulting in G M D errors of -0.549855 • 1 0 - 4 percent and -0.650632 • 10 - 4 percent for core and sheath, respectively. The system is thus represented with a total of 5,090 partial subconductors. As shown the differences between the inductance values of the P S E C method and those 6.2. Coaxial Cable 94 1.6 1.2 r-loop resistance 25.6 40.0 f <kHz> 102.4 loop inductance - - analytic — PSEC 0 1.6 6.4 14.4 25.6 40.0 57.6 78.4 102.4 f <kHz> Figure 6.2: Resistance and inductance of coaxial cable (Simulation b). 6.2. Coaxial Cable 95 of the analytic method are not noticeable. The maximum and minimum inductance errors are -0.065 percent and -0.002 percent, occurring at 10 kHz and 400 Hz, respectively. However, the cut-off effect appears on the resistance curve at a lower frequency (i.e., the resistance error becomes negative between 7 kHz and 10 kHz) since the resolution used to discretize the core in this case was lower than that used to discretize it in Simulation a. The maximum and minimum resistance errors are —7.99 percent and 0 percent, taking place at 100 kHz and 0.1 Hz, respectively. The average simulation time (PSEC program) was 58.2 min per frequency using the 733-MHz computer described above. 6.2.3 Simulation c Fig. 6.3 illustrates the results obtained after using the images and corrected resolutions listed in Tables 4.1, 4.2 and 4.3 (pages 71 and 72). As explained in Section 4.4, each reso-lution is calculated according to the thickness of conductor and the penetration depth per frequency, and the core conductor is represented with its skin layer at very low penetration depths. As shown in the plots, the cut-off effect does not appear on the resistance curve above frequencies of 40 kHz when the boundary representation is used. The maximum and minimum differences between the values of the P S E C method and those of the analytical method are 1.91 percent and 0 percent on the resistance curve, occurring at 70 kHz and 0.1 Hz, and —0.94 percent and —0.0023 percent on the inductance curve, taking place at 40 kHz and 60 Hz, respectively. The differences between the values of the P S E C method and those of the analytic method are not noticeable below frequencies of 10 kHz. The saving on number of sub-conductors is also substantial. For example, no more than 686 subconductors are used at frequencies below 1 kHz and no more than 2,725 subconductors are needed at 100 kHz. As compared to Simulation b, Simulation c yields greater inductance errors (negative errors) above frequencies of 40 kHz. This is due to the fact that the flux linking the inner 6.2. Coaxial Cable 96 loop resistance i 1 1 1 1 r — PSEC I I I I I I I I I 0 1.6 6.4 14.4 25.6 40.0 57.6 78.4 102.4 f <kHz> loop inductance 0 1.6 6.4 14.4 25.6 40.0 57.6 78.4 102.4 f <kHz> Figure 6.3: Resistance and inductance of coaxial cable (Simulation c). 6.2. Coaxial Cable 97 subconductors (the ones at penetration depths of 25 and 35, beyond the skin layer of 15 used to represent the core conductor) have been neglected. However, the inductance error is still low and did not exceed —0.94 percent. Above 40 kHz, the resistance is also higher than that of the analytic solution since the skin layer does not represent the total area through which the currents flow (wave penetration and current flows beyond one 5 are not considered). In spite of that, the resistance error was also small (below 1.91 percent). For example, the differences are barely noticeable if a logarithmic scale is used to plot the results (Fig. 6.4). loop resistance 2 I . — I • — I ' — I ' — f <Hz> loop inductance f <Hz> Figure 6.4: Resistance and inductance of coaxial cable (Simulation c, log scale). To analyze the relationship between error and thickness of skin layer, the thickness of the skin layer of the core conductor was modified at one frequency (40 kHz). Layers with thicknesses of 35 (subdivided into three sublayers of one 6 each), layers with thicknesses of 25 (subdivided into two sublayers of one 5 each), and layers with thicknesses of 1.55 were tested. 6.2. Coaxial Cable 98 As expected, the loop inductance errors improved as the layer became thicker (—0.1431 percent for 35 thickness, -0.3265 percent for 25 thickness, and -0.4899 percent for 1.56 thickness). However, the cut-off effect associated with the loop resistance showed up again (the resistance errors were —5.3847 percent, —8.3865 percent and —10.3755 percent, respectively). Alternative solutions to minimize the errors at low penetration depths are discussed in the section of future work of the chapter of conclusions (Chapter 7). Section 6.3 shows the curves of the impedances with return through ring Z c c , Z c s and Z s s . Additional results showing the effect of the resolution in the calculations can be found in [58]. Table 6.1 tabulates the computer timings per frequency obtained with the 733-MHz computer. Table 6.1: Simulation time versus number of subconductors (coaxial cable, Simulation c). Frequency Number of Simulation Time <Hz> Subconductors 1E-06 404 2.74 s 0.1 404 2.14 s 1.0 404 2.15 s 10.0 404 2.14 s 50.0 404 2.31 s 60.0 404 3.13 s 100.0 404 4.12 s 400.0 576 5.44 s 700.0 576 5.49 s 1E+03 686 10.16 s 4E+03 1,941 3.93 min 7E+03 3,186 17.13 min 1E+04 4,506 44.85 min 4E+04 1,262 1.1 min 7E+04 2,017 4.42 min 1E+05 2,725 10.99 min 6.2.4 Summary of Results Tables 6.2 and 6.3 summarize the results of the three simulations. For Simulation b (Fig. 6.2) the resistance does not change noticeably as compared to that of Simulation a (Fig. 6.1), but the inductance for Simulation b greatly improves as compared to that of Simulation a. The improvement in the inductance results is due to the effectiveness of the G M D discretization error reduction algorithm. The results of 6.2. Coaxial Cable 99 Table 6.2: Loop resistance as a function of frequency for coaxial cable. frequency analytic Simulation a Simulation b Simulation c resistance resistance error resistance error resistance error <Hz> <n/km> <fi/km> <%> <0./km> <%> <n/km> <%> 1E-06 0.415578 0.415578 -4.812577E-07 0.415578 -4.812577E-07 0.415578 -4.812577E-07 0.1 0.415578 0.415578 0.000000 0.415578 0.000000 0.415578 0.000000 1.0 0.415579 0.415579 2.406284E-07 0.415579 -2.165656E-06 0.415579 1.227205E-05 10.0 0.415655 0.415655 -4.835738E-05 0.415655 -1.012859E-04 0.415661 1.279666E-03 50.0 0.417405 0.417401 -1.023706E-03 0.417396 -2.132939E-03 0.417517 2.665204E-02 60.0 0.418143 0.418138 -1.352646E-03 0.418131 -2.854523E-03 0.418296 3.646931E-02 100.0 0.421810 0.421800 -2.324745E-03 0.421788 -5.288166E-03 0.422153 8.145473E-02 400.0 0.445352 0.445385 7.483068E-03 0.445371 4.240690E-03 0.445484 2.961163E-02 700.0 0.459583 0.459670 1.894023E-02 0.459644 1.340477E-02 0.459761 3.874951E-02 1E+03 0.471146 0.471294 3.137947E-02 0.471252 2.254780E-02 0.471956 0.171977 4E+03 0.544636 0.545434 0.146522 0.545071 7.985726E-02 0.547557 0.536337 7E+03 0.597182 0.598187 0.168380 0.597326 2.427302E-02 0.601428 0.711133 1E+04 0.644366 0.645129 0.118511 0.643610 -0.117262 0.648481 0.638708 4E+04 1.104503 1.083251 -1.924147 1.072261 -2.919138 1.120123 1.414186 7E+04 1.514878 1.451906 -4.156925 1.431088 -5.531122 1.543875 1.914172 1E+05 1.834888 1.715868 -6.486460 1.688326 -7.987478 1.866049 1.698276 Table 6.3: oop inductance as a function of frequency for coaxial cable frequency <Hz> analytic inductance </nH/km> Simulation a Simulation b Simulation c inductance <fiH/km> error <%> inductance <^H/km> error <%> inductance <f*H/km> error <%> 1E-06 139.743472 138.652430 -0.780746 139.683630 -4.282275E-02 139.915350 0.122995 0.1 139.743466 138.652430 -0.780742 139.683620 -4.282562E-02 139.915350 0.123000 1.0 139.742869 138.651830 -0.780748 139.683030 -4.282079E-02 139.914670 0.122941 10.0 139.683311 138.592590 -0.780853 139.624070 -4.241094E-02 139.847390 0.117465 50.0 138.329990 137.246210 -0.783474 138.283680 -3.347792E-02 138.369240 2.837418E-02 60.0 137.762621 136.681510 -0.784764 137.721320 -2.997983E-02 137.759480 -2.280009E-03 100.0 134.981303 133.911030 -0.792905 134.960650 -1.530064E-02 134.790790 -0.141140 400.0(*) 120.616647 119.551060 -0.883449 120.613950 -2.236010E-03 120.955040 0.280553 700.0 116.334045 115.261700 -0.921781 116.324160 -8.497083E-03 116.614410 0.241000 1E+03 114.175803 113.100060 -0.942181 114.162890 -1.130975E-02 113.900090 -0.241481 4E+03 108.604180 107.487040 -1.028634 108.556290 -4.409591E-02 108.312930 -0.268176 7E+03 107.224307 106.088880 -1.058927 107.162440 -5.769867E-02 107.026680 -0.184312 1E+04 106.508679 105.362480 -1.076156 106.439960 -6.451963E-02 106.351630 -0.147452 4E+04 104.190030 103.038420 -1.105298 104.150730 -3.771954E-02 103.210810 -0.939840 7E+04 103.187786 102.065440 -1.087673 103.201600 1.338724E-02 102.384540 -0.778431 1E+05 102.596750 101.501180 -1.067841 102.659170 6.084013E-02 101.928110 -0.651717 (*) Simulation c used an 11.489776-pel/in resolution at 400 Hz (instead of a 9.703081-pel/in resolution) since the inductance error decreased from 0.93 percent to 0.28 percent. The skin effect is more developed at 400 Hz (the penetration depth is 30.31 percent at that frequency), and the jagged shape of the core conductor for a 9.703081-pel/in resolution might have altered current distributions and proximity effects a bit more. 6.3. Buried Cable System 100 Simulation c (adaptive resolution scheme, Fig. 6.3) are the best for the resistance (note that the cut-off effect does not appear on the curve) and show the effectiveness of the boundary representation of the conductors at high frequencies. The inductance for Simulation c (maximum error of —0.94 percent) worsens as compared to that of Simulation b (maximum error of —0.065 percent) due to the fact that the inner subconductors of the core conductor are neglected at high frequencies when the boundary representation is used. Nonetheless, the inductance error in Simulation c is still low. As shown in Table 6.1, the timings are reduced noticeably when the adaptive resolution scheme is used. For 16 frequencies, Simulation c was completed in l h 23 m 6 s , whereas Simulation a and Simulation b were accomplished in 15.94 h and 15.51 h, respectively. 6.3 Bu r i ed Cab le System Fig. 6.5 depicts the geometry of a three-phase system of buried cables. Each phase con-sists of a single-core coaxial cable similar to the one illustrated in Fig. 4.1 (page 51) and whose loop parameters were calculated with the P S E C program in Section 6.2. The laying conditions (spacing of 2 diameters between cables and depth of burial of 1 m) are typical of 10-kV systems and were taken from [5]. The value of 100 Jl-m for the earth resistivity is also typical in this type of study. The equations system (3.8) on page 33 describes the voltages of the system as a function of currents and cable parameters. Due to symmetry Z a b — Zbc. To calculate the parameters of the conductors, Section 6.2 assumed current return through a ring of radius R, where R is the outer radius of the outermost insulation. 6.3. Buried Cable System 101 air / / / / / A earth h=lm a b c R=AA2 rnmL.^0^ \ d~ 17.68 cm d ' d ' Figure 6.5: Geometry of buried cable system. Hence the parameters of the conductors obtained from the P S E C program include the impedance Z s h e ath/earth-insuiation) a n d the impedances needed to complete the equations are the mutual earth-return impedances Z a b , Z a c and the self earth-return impedance1 Z s e i f earth-return, which must be added to the impedances with ring return Z c c , Z c s and Z s s . The self and mutual earth-return impedances can be calculated with Pollaczek's and Wedepohl's formulas on the assumption that the cables are filaments inside the infinite earth-return path, as explained in Section 3.6, page 37. Wedepohl's formulae are valid in this case since \mR\ — 0.004 < 0.25 and | m d a c | = 0.03 < 0.25 at the highest frequency of interest (100 kHz). Figs. 6.6 and 6.7 depict the results obtained for the impedances Z c c and Z c s using the frequencies listed in Tables 6.2 and 6.3 (from 1 Hz to 100 kHz) and the results of Simulation c. The results for the impedances with earth return Z s s , Z a b and Z a c as well as the curves for the impedances with return through ring are illustrated in Figs. 6.8, 6.9, 6.10, 6.11, 6.12, and 6.13. The impedance Z s e i f earth-return is two orders of magnitude greater than the impedances of the conductors. Therefore the impedances with return through earth Z c c , Z c s and Z s s 'Zseif earth-retum is part of Z 2 2 , and Z 2 2 is part of the impedances with return through earth Z c c , Z c s and Z s s (see equation of Z 2 2 on page 30 and equation (3.8) on page 33). 6.3. Buried Cable System 102 self resistance of core with return through earth 1201 . — . — • — ' — I f <Hz> self inductance of core with return through earth 31 1 • •— • •— I • •— f <Hz> Figure 6.6: Self core impedance Z c c of buried cable system. 6.3. Buried Cable System mutual resistance between core and sheath with return through earth 120 A E XL O V CO O or 100 80 40 20 P S E C + Po l laczek P S E C + W e d e p o h l 10° 101 102 103 104 10 s f <Hz> mutual inductance between core and sheath with return through earth E 2 5 F 1 2 v CO q 1.5 P S E C + Po l laczek P S E C + W e d e p o h l 10 10 10 10 f <Hz> 10 10 Figure 6.7: Mutual core-sheath impedance Z c s of buried cable system. 6.3. Buried Cable System 104 display curves similar to those of the self earth-return impedance. As expected, the resistance increases and the inductance decreases as the frequency grows. The differences in the resistance are not noticeable. The differences in the in-ductances, however, are present in the entire range of frequencies. For example, the in-ductances corrected with Wedepohl's model are always greater than those corrected with Pollaczek's model. In spite of that, the results from Wedepohl's method can be considered very accurate given the simplicity of the formulae. The differences are due to the fact that Wedepohl's formulae only consider the first terms of the infinite series. Tables 6.4 tabulates the maximum and minimum differences between the impedances corrected with Wedepohl's formulae and those corrected with Pollaczek's formulae. For resistance and inductances the differences are maximum at 100 kHz and minimum at 1 Hz. For the resistance the differences vary between 0.00322 percent and 3.56 percent, whereas for the inductance the differences change from 10.14 percent to 29.77 percent. Table 6.4: Differences in percentage terms between corrected impedances (Wedepohl's Impedance Error real part <%> Error imaginary part <%> max min max min Z c c 3.35 6.41E-02 18.31 10.14 Z Ss 3.41 3.22E-03 19.86 10.76 ZCs 3.41 1.32 19.86 10.75 Z a b 3.48 1.32 25.57 12.06 Zac 3.56 1.32 29.77 12.82 It is also important to mention that the P S E C program incurred in an maximum error of -2.4 percent (at 10 kHz Lss = 12.3 ^ H / k m with the P S E C program while Lss = 12.6 £ iH/km with the analytic solution) when calculating the self inductance of sheath with return through ring (see Fig. 6.12). This effect is caused by the jagged shape of the sheath and the development of the skin effect in the core conductor. Those phenomena, in turn, increase proximity effects between conductors and decrease the inductance of the sheath. Nonetheless, the self inductance of sheath with return through earth was not affected since the self earth-return impedance is much greater than the impedance of the 6.3. Buried Cable System 105 conducting sheath (as previously mentioned, the earth-return impedance is two orders of magnitude greater). For a very low earth resistivity the solution of such a problem is as follows: the error is drastically reduced at 40 kHz when the conducting sheath is modeled with 811 subconduc-tors (see Fig. 6.12). Therefore, the error will also be reduced at frequencies between 100 Hz and 10 kHz, inclusive, using 811 subconductors instead of 264 subconductors (equivalent to change the resolution of the conducting sheath from 17.7 pel/in to 26.4 pel/in). The simulation time will not increase substantially. The most critical case is that of 10 kHz where 5,053 subconductors should be used (approximately a 1-h simulation with the 733-MHz computer). After those changes, the total simulation time for 16 frequencies should be around 2 h. self resistance of sheath with return through earth 120 100 80 J 60 O V 40 co CO CC 20 P S E C + Po l laczek P S E C + W e d e p o h l 0 10' 10 10 10 10' f <Hz> self inductance of sheath with return through earth P S E C + Po l laczek P S E C * W e d e p o h l 10 10 10 10 f <Hz> 10 10" Figure 6.8: Self sheath impedance Z s s of buried cable system. 6.3. Buried Cable System 106 mutual resistance between sheath a and sheath b with return through earth 1201 • — >— . — • — ' II 10° 101 102 1 0 3 1 0" 105 f <Hz> mutual inductance between sheath a and sheath b with return through earth 2.51 • — • — • I • — 1 I P S E C + Pol laczek 10° 101 102 103 10" 105 f <Hz> Figure 6.9: Mutual impedance Z a b = Z b c of buried cable system. 6.3. Buried Cable System 107 mutual resistance between sheath a and sheath c with return through earth 1201 . i — • ' — • — I 10° 101 102 1 0 3 1 0" 105 f <Hz> mutual inductance between sheath a and sheath c with return through earth 2.51 • — i — • . . . ..| .—•— •— •—•— i •—•— MI P S E C + Pol laczek 10° 101 102 103 10* 105 f <Hz> Figure 6.10: Mutual impedance Z a c of buried cable system. 6.3. Buried Cable System 108 self resistance of core with return through ring 3.51 , — 1 — 1 — 1 — f <Hz> self inductance of core with return through ring f <Hz> Figure 6.11: Self core impedance with return through ring (coaxial cable). 6.3. Buried Cable System self resistance of sheath with return through ring 1.21 1 • • -< f <Hz> Figure 6.12: Self sheath impedance with return through ring (coaxial cable). 6.3. Buried Cable System 110 mutual resistance between core and sheath with return through ring 1.51 , i • • 1 ' P S E C f <Hz> mutual inductance between core and sheath with return through ring 14.51 , .— i • •— i • • — ' I • I f <Hz> Figure 6.13: Mutual core-sheath impedance with return through ring (coaxial cable). 6.4. Sector-Shaped Cable 111 6.4 Sector-Shaped Cab le Fig. 6.14 shows the geometry of a sector-shaped cable studied with the Finite Element Method in [75] and with approximate formulae in [4]. Its parameters, which can not be obtained through closed-form formulae, are calculated with the P S E C algorithm and compared with those obtained with the other two methods. The area of each sector-shaped conductor of the core is 299.76 mm 2 and the area of the conducting sheath is 326.73 mm 2 . The sheath is considered the return path and the ground is not included in the calculations. Also, Z n — Z22 = Z 3 3 and Z 1 2 = Z 1 3 = Z 2 3 due to symmetry. H=1 for all conductors and insulators 27 mm 2 5 m m 19 mm 4.255 mm core b =58000 S/mm a= 0 sheath a =1100 S/mm Figure 6.14: Geometry of sector-shaped cable. After applying segmentation, the cable geometry is divided into three images, one for conductor 1, one for conductor 2, and the other for the conducting sheath (conductor 3 is not needed because of the symmetry). Using square-shaped subconductors, resolutions of 65.5 subconductors/in and 46 subconductors/in result in area discretization errors of 0.02 percent and 0.01 percent for the conductors and the conducting sheath, respectively. With the above resolutions, each sector-shaped conductor is subdivided into 1,994 subconductors, the sheath is subdivided into 1,055 subconductors, and the cable geometry is thus represented with a total 5,043 partial subconductors. Since the lengths of the hypotenuses of the square subconductors of core and sheath are 6.4. Sector-Shaped Cable 112 0.55 mm and 0.78 mm, respectively, resistance values smaller than the actual resistances should start to appear after 14.44 kHz and 416.34 kHz on the core and sheath curves, respectively. Figs. 6.15 and 6.16 depict the results obtained for the loop parameters using a scale proportional to v7 along the abscissae and frequencies up to 600 kHz. The results of the P S E C method are plotted with 10 points per decade whereas the results of the approximate method and the F E method are only plotted for the following six points: 6 Hz, 60 Hz, 600 Hz, 6 kHz, 60 kHz and 600 kHz. As shown, the P S E C method and the F E method display similar curves, and a cut-off effect is hardly noticed around 490 kHz. This result can be explained analyzing the conductivities of the conductors. The conductivity of the sheath is much smaller than the conductivity of the core, so the internal resistance of the sheath is always much greater than the internal resistance of the core. As a result, the cut-off frequency of the loop resistances is given by that of the internal resistance of the sheath. It can also be noticed that the approximate method yields very good estimates of self resistance and self inductance. For those parameters, substantial differences between the results of the approximate method and those of the numerical methods are only visible above 60 kHz, where the approximate method gives much greater values of self resistance, and below 60 Hz, where it produces much smaller values of self inductance. For the mutual parameters, the approximate method yields results too different from those of the numerical methods. Tables 6.5, 6.6, 6.7 and 6.8 summarize and compare the results of the three methods. Table 6.5: Sel resistance as a function of frequency for sector-shaped cable frequency Approximate [75] F E [75] PSEC resistance dif. PSEC-Approx. resistance dif. PSEC-FE resistance <Hz> <n/km> <%> <n/km> <%> <n/km> 6.0 2.83996 1.55 2.84006 1.55 2.8840129 60.0 2.84334 1.66 2.84987 1.43 2.8906402 600.0 2.93723 2.11 2.96756 1.06 2.9991337 6E+03 3.52328 6.29 3.52505 6.24 3.7448973 6E+04 5.89100 -7.89 5.15051 6.01 5.4602898 6E+05 18.8742 -27.03 15.9166 -7.12 14.858436 6.4. Sector-Shaped Cable Figure 6.15: Self resistance and self inductance of sector-shaped cable. 6.4. Sector-Shaped Cable mutual resistance E o v 160 250 f <kHz> mutual inductance 640 90 160 250 f <kHz> 640 Figure 6.16: Mutual resistance and mutual inductance of sector-shaped cable. 6.4. Sector-Shaped Cable 115 Tab e 6.6: Se f inductance as a function of frequency for sector-shaped cable frequency <Hz> Approximate [75] F E [75] PSEC inductance </zH/km> dif. PSEC-Approx. <%> inductance </iH/km> dif. PSEC-FE <%> inductance </tH/km> 6.0 196.081 17.72 232.020 -0.52 230.82391 60.0 191.259 17.56 220.673 1.89 224.84486 600.0 164.736 5.98 156.797 11.34 174.58279 6E+03 125.467 7.93 120.400 12.48 135.42119 6E+04 108.540 6.01 105.090 9.49 115.06827 6E+05 100.191 9.26 99.0031 10.58 109.47334 Table 6.7: Mutual resistance as a function of frequency for sector-shaped cable frequency <Hz> Approximate [75] F E [75] PSEC resistance <f!/km> dif. PSEC-Approx. <%> resistance <f)/km> dif. PSEC-FE <%> resistance <f2/km> 6.0 2.78239 1.58 2.78246 1.58 2.8264260 60.0 2.78105 1.64 2.78317 1.57 2.8267877 600.0 2.75504 2.25 2.78308 1.22 2.8170241 6E+03 2.58640 7.07 2.69462 2.77 2.7693209 6E+04 2.22645 37.85 2.88824 6.26 3.0690638 6E+05 6.42831 41.56 8.82314 3.13 9.0996821 Table 6.8: Mutual inductance as a function of frequency for sector-shaped cable frequency <Hz> Approximate [75] F E [75] PSEC inductance </iH/km> dif. PSEC-Approx. <%> inductance </iH/km> dif. PSEC-FE <%> inductance </iH/km> 6.0 25.0563 51.21 40.4454 -6.75 37.887812 60.0 27.3958 34.40 40.1914 -9.15 36.821023 600.0 36.4771 3.26 35.7225 5.44 37.665900 6E+03 47.2753 -13.02 38.3245 9.14 41.827974 6E+04 52.1332 -22.90 40.6921 4.25 42.420645 6E+05 51.1033 -29.09 38.0609 4.01 39.587825 For the self resistance the maximum and minimum differences between the values of the P S E C method and those of the F E method are —7.12 percent and 1.06 percent, occurring at 600 kHz and 600 Hz, respectively. Between the values of the P S E C method and those of the approximate method, the maximum and minimum differences are —27.03 percent and 1.55 percent, taking place at 600 kHz and 6 Hz, respectively. For the self inductance the maximum and minimum differences between the results of the P S E C method and those of the F E method are 12.48 percent and —0.52 percent, occurring at 6 kHz and 6 Hz, respectively. Between the results of the P S E C and those of the approximate method, the maximum and minimum differences are 17.72 percent and 5.98 percent, taking place at 6 Hz and 600 Hz, respectively. 6.5. Summary 116 For the mutual parameters the maximum and minimum differences between the values of the P S E C method and those of the F E method are 6.26 percent and 1.22 percent on the resistance curve, occurring at 60 kHz and 600 Hz, and -9.15 percent and 4.01 percent on the inductance curve, taking place at 60 Hz and 600 kHz, respectively. The average simulation time was 58.32 min per frequency with three partitions for matrix modification. 6.5 Summary Typical power cable arrangements have been studied with the developed techniques. The results have been compared with those of the analytic, F E and approximate methods and are in agreement with those of the first two methodologies. The cable geometries were represented with standard bitmapped files and were redis-cretized with the developed R D algorithm. The results improved noticeably when the G M D error minimization technique was used. The simulation timings also decreased sub-stantially when the resolution of the image adapted itself to the thickness of conductor and the penetration depth per frequency. Furthermore, the representation of the conductors as edges contributed to a considerable reduction of the timings at frequencies with low penetration depths. In the case of underground cables, the developed programs can be used in connection with Pollaczek's and Wedepohl's formulae. As opposed to the F E method, the proposed techniques do not require an explicit modeling (subdivision) of the earth-return path. The earth-return path can be incorporated later on once the parameters of the conductors has been calculated by the P S E C program. This is possible since the penetration depth of the earth-return path is generally much greater than the cross-section of the conductor. Therefore, the conductors can be considered infinitely thin filaments as compared to the earth-return path. The proposed methodologies could also be applied to tunnel-installed underground 6.5. Summary cables as will be explained in the section of future research of the chapter of conclusions (Chapter 7). Chapter 7 Conclusions and Recommendations for Future Research 7.1 Conclusions A new methodology for the calculation of the frequency-dependent parameters of arbitrar-ily-shaped power cable arrangements has been presented. The proposed algorithm uses digital images to discretize the cable geometry and the partial subconductor equivalent circuit method to estimate the cable parameters. The digital images are rediscretized according to the thicknesses of the conductors and the penetration depths per frequency, and the discretization errors are corrected with an it-erative procedure. To reduce the dimension of the problem, the conductors are represented with their edges at very low penetration depths. The partition of the partial subconductor impedance matrices has also been proposed. Thus, the P S E C program can estimate the parameters when the sizes of the matrices exceed the available physical memory, even in low R A M computers. Coaxial cables, buried cables, and sector-shaped cables have been studied with the developed program, and it has been shown that the results obtained are in good agreement with those of classical methods and finite element methods. The developed algorithms can also be used to calculate the frequency-dependent pa-118 7.1. Conclusions 119 rameters of oval-cables, L-shaped cables, cables with concentric neutral conductors, square busbars, and power rails of traction systems. The main contributions of the present thesis work can be summarized as follows: 1. It was proved that it is feasible to use digital images for an accurate calculation of the frequency-dependent parameters of arbitrarily-shaped power cables. These digital images are obtained from off-the-shelf computer graphics formats for human vision such as bitmapped files, which are available to most computer users. 2. A technique to reduce G M D discretization errors has been developed. Given a user-defined tolerance, an iteration procedure and a linear interpolation scheme are in-troduced to correct the resolution of the input images. The proposed technique allows a simple rectangular model, specifically a bitmapped file, to represent curved conductors without incurring in high discretization errors. An essential requirement to use the proposed technique is the knowledge of the actual dc G M D s of the conductors, which for simple geometries are given by formulae, and which for unusual shapes can be obtained from the subdivision of the cross sections of the conductors into simpler shapes, e.g., triangles. The knowledge of the actual dc resistances of the conductors is also paramount, but the resistance corrections are simply performed by the division of the actual areas among all the subconductors. The proposed technique was developed by drawing an analogy with the methods used to solve non-linear equations, in which a residual is approximated to values lower than a specified tolerance by consecutive iterations (e.g., load flow equations). 3. A n adaptive resolution scheme has been proposed and implemented. At low fre-quencies the proposed scheme assigns low resolutions to the conductors, whereas at higher frequencies it assigns higher resolutions. The proposed adaptive resolution scheme reduces the requirements of simulation time and memory noticeably. For very low penetration depths (typically, below five percent of the thickness of the conductor), the adaptive scheme switches to a boundary representation. This 7.2. Recommendations for Future Research 120 boundary representation only subdivides the skin layers associated with the edges of the conductors, and reduces the resistance errors substantially. Simple criteria to choose the resolutions and to estimate the cut-off frequencies have been proposed as well. 4. The partition methodology is another contribution of the present thesis work. The main advantage of the proposed method is that the maximum number of subcon-ductors is not limited by the amount of physical memory of the computer, but by the amount of storage (hard disk) capacity, which is usually very high. However, a disadvantage of the method is that the simulation time increases noticeably when the number of partitions for matrix reduction is high. In spite of that, the developed program is still capable of calculating the reduced impedance matrix. 7.2 Recommendat ions for Future Research 1. For low penetration depths the resistance and inductance errors should diminish, and for intermediate penetration depths the number of subconductors should be lower. To achieve such reductions the following method can be explored: The R D algorithm could extract, one by one, the outer concentric layers of subcon-ductors from a conductor to create a new image with the coordinates of the extracted subconductors. For instance, it could be possible for the program to extract layers until reaching the penetration depth 5, and then bundle the remaining inner subcon-ductors with formula (4.9) on page 55. Next, the equivalent conductor of the bundled subconductors could be placed at the geometric center, and the G M D s among the outer subconductors and the equivalent subconductor could be approximated by the mutual distances. Thus, the program could be modified such that when it switches to a boundary representation (e.g., at a low penetration depth), there will be magnetic flux and current through an inner equivalent subconductor, and the inductance and resistance of the conductor will slightly increase and decrease, respectively. Different schemes 7.2. Recommendations for Future Research 121 of subdivision could be analyzed as well. For instance, the program could choose the resolution so that the layers have thicknesses of half a penetration depth, or could extract the concentric layers located at two or three penetration depths before performing the bundling of the inner subconductors. The program could also bundle neighboring subconductors located in the same layer if more memory savings were desired. 2. The introduction of a second interpolation into the G M D error minimization tech-nique could further improve the accuracy of the results. After finding the interval where the G M D error function crosses through zero, the R D program could reduce the resolution step (e.g., from 0.1 pel/in to 0.01 pel/in) to obtain a better estimate for the resolution. As a result of such an error reduction, it would also be possible to correct the self G M D s of the subconductors. To do so, the R D program could calculate the value of G M D j using (4.9) along with the actual G M D of the conductor and the products of the corrected mutual G M D s of the subconductors. 3. The parameters of tunnel-installed underground cables could be calculated as follows: If a fictitious circular ring encloses the tunnel-installed cables, the R D and P S E C programs will be able to determine the parameters, but one more conductor will be present in the reduced impedance matrix. This additional conductor is the partial earth-return path enclosed by the ring, and it should be bundled with the rest of the earth-return path. In this case Pollaczek's self earth-return impedance can be calculated on the assump-tion that the radius of the outermost insulation (or inner radius of the semi-infinite earth-return path) is the radius of the fictitious ring. Therefore, it is possible to add the self earth-return impedance to all the ring-return impedances of the reduced ma-trix, and then to eliminate the row and column of the partial earth-return path with a Kron's reduction. Since the voltage of the partial earth-return path with respect to the semi-infinite earth-return path is equal to zero, this elimination is similar to that performed for ground wires of overhead transmission lines. 7.2. Recommendations for Future Research 122 4. The discrete model of the cable geometry could be obtained from photographs. This would be particularly useful for the analysis of cases with unusual shapes or without blueprints. 5. Mathematical information on the boundaries of the conductors and, therefore, on the value of their actual areas and GMDs, should be directly obtainable either from vector graphics files, which are traditionally used to store line art and C A D infor-mation, or from metafiles, which contain either bitmapped or vector graphics data [45]. 6. A support routine for the calculation of the dc G M D s of cables with unusual shapes should be developed. The basic idea is to subdivide the cross sections into triangles and then to bundle those triangles using (4.9). This routine would also be useful to calculate the actual G M D of the partial earth-return path associated with tunnel-installed underground cables. The sides of the triangles delimiting the boundaries of the conductors should be sufficiently short, but the mesh do not have to have a complicated structure since it is only required for dc calculations. The formulae for the self and mutual G M D s of triangles are given in [67], and the information on the boundaries of the conductors could be obtained from the vector graphics files used to represent the cable geometry. 7. The partition algorithm for matrix reduction might become faster if the P S E C pro-gram were to load two or more columns of factors every time it has to access the hard disk. This would imply more partitions for the rest of the elements, but the overall performance of the partition algorithm might improve. 8. For the case of pipe-type cables, the accuracy of the P S E C program should be eval-uated. The reason is that the subconductor method is based on the assumption that all the solution region has the same permeability, and the pipe of pipe-type cables is a magnetic conductor. For example, previous investigations [7], [78] have proposed that the subconductor method should simply assign the proper value of pr (relative permeability) to the subconductors located in the pipe. 7.2. Recommendations for Future Research 123 9. The developed programs should include routines for the calculation of the cable ca-pacitances. A methodology that can be explored is the advanced method of subareas [77]. This method subdivides the boundaries of the conductors into segments, and then bundles the potential coefficients of the segments to obtain the self and mutual potential coefficients of the conductors. 10. Analytic methods such as those in [22], [33] could be incorporated into the subroutine T U B E [19], [20]. In this way, T U B E would be able to analyze pipe-type cables as well. 11. The solution method could be hybridized [67], i.e., the P S E C method could be used for low frequencies, whereas a surface integral equation analysis could be used for high frequencies. 12. The currents in the subconductors could be calculated after obtaining the reduced impedance matrix. To do so (5.12) should be used. This would allow the graphic visualization of the skin and proximity effects in the cable system under study. List of Symbols Alphanumeric Symbols boldface vector quantity or complex quantity A , A magnetic vector potential, W b / m Aj area of subconductor i, m 2 A;; self impedance after row and column operations for equations mod-ification (Chapter 5), f2/m Aik mutual impedance after row and column operations for equations modification (Chapter 5), Q/m Ak area of subconductor k, m 2 An value of field distribution A at node n, W b / m a radius of ring-return path, m a, b, c subscripts of phases a, b, and c, respectively B , B magnetic flux density, T = W b / m 2 B\ capacitive susceptance of insulation layer between core conductor and conducting sheath, 15/m B2 capacitive susceptance of insulation layer between conducting sheath and earth, 15/m Bcc self shunt susceptance of core conductor, 15/m Bcs mutual shunt susceptance between core conductor and conducting sheath, 15/m Ber(v) Bessel function Bei(v) Bessel function Bii self impedance after row and column operations for equations mod-ification, 0,/m Bjt mutual impedance after row and column operations for equations modification, fi/m Bss self shunt susceptance of conducting sheath, 15/m Ci shunt capacitance of tubular insulation, F / m ceil(x) function rounding x to the nearest greater integer cn column number of subconductor cosech(x) hyperbolic cosecant of x coth(x) hyperbolic cotangent of £ D electric flux density, C / m 2 D distance between cable i and image of cable k, m d distance between cables i and k, m 124 List of Symbols 125 dik E , E Es, Es e exp(x) F ; / floor(x) G M D G M D j G M D , , G M D i r G M D f e G M D f e r G M D r H , H h hi hk I, / Io II h Ic I •••core Em Is Isheath i lend Ip ir ^start J , J Jsk J. J, j K distance between subconductors i and k, m electric field intensity, V / m source electric field, V / m Neper's number exponential function e x multiplication factor for Gaussian elimination frequency, H z function rounding x to the nearest lower integer geometric mean distance, m self G M D of subconductor i, m mutual G M D between subconductors i and k, m mutual G M D between subconductor i and ring-return path, m self G M D of subconductor k, m mutual G M D between subconductor k and ring-return path, m self G M D of ring-return path, m magnetic field, A / m average depth at which two cables are buried, m; depth at which cable is buried, m depth at which cable i is buried, m depth at which cable k is buried, m current, A modified Bessel function with complex argument, first kind, order 0 current flowing through loop "core/sheath", A; modified Bessel function with complex argument, first kind, order 1 current flowing through loop "sheath/ground", A current flowing through core conductor, A current flowing through subconductor i of core conductor, A current flowing through core conductor, A imaginary part of complex number current flowing through conducting sheath, A current flowing through conducting sheath, A current flowing through subconductor i of conducting sheath, A subscript of matrix or vector element; index of columns for matrix reduction stopping index counter of partitions counter of rows starting index volume current density, A / m 2 source current density of the kth conductor, A / m 2 current density component along z direction, A / m 2 imaginary number number of conductors List of Symbols 126 K 0 modified Bessel function with complex argument, second kind, order 0 K i modified Bessel function with complex argument, second kind, order 1 k subscript of matrix or vector element; counter of first rows for matrix modification Lcc' loop inductance between core and sheath, H / m Lcore-out inductance of Z c o r e _ o u t , H / m Lcore/sheath—insulation inductance of Z c o r e/ sheath—insulation! H / m La self inductance of loop "subconductor i/ring return", H / m Lik mutual inductance between loops "conductor z/ring return" and "conductor /c/ring return", H / m LD dc inductance of conductor, H / m Lsheath-in inductance of Z s h e a t h - i n , H / m Lss inductance of Z s s , H / m £ sum of the depths of burial of two cables (Chapter 3), m; side of subconductor (Chapter 5), m £ c perimeter of cross section of conductor, m £i side of subconductor i, m 4 side of subconductor k, m M physical memory of computer minus memory used by operating system, auxiliary arrays, and variables, bytes Mfr memory used by first rows of impedance matrix when performing matrix modification, bytes m reciprocal of complex penetration depth, m - 1 m number of conductors rnc number of equal parallel subconductors into which a conductor is subdivided N number of samples per lobe for Pollaczek's integrand (Chapter 3); memory needed by impedance matrix (Chapter 5), bytes; number of nodes with unknown values (Appendix D) NB number of boundary nodes with known values NT total number of nodes n number of subconductors ne number of elements of symmetrical impedance matrix after exploit-ing symmetry nep\ maximum number of elements per partition for matrix construction nepif number of elements of last partition for matrix construction nep2 maximum number of elements per partition for matrix modification n e p 3 maximum number of elements per partition for matrix reduction nf number of frequencies n m a x maximum number of subconductors np\ number of partitions for matrix construction nP2 number of partitions for matrix modification n P 3 number of partitions for matrix reduction List of Symbols 127 P P P c pel ft R R(A) RJ Rdc Me Resol _c Ri Rn R0 f rn res rf Tin rn To 1out rv fsh—in Tsh—oul s Sck SR sf t thickness U R U Z U ^ v,v V 2 Pollaczek's integral complex penetration depth, m complex penetration depth of conductor, m pixel inside radius of insulation, m outer radius of outermost insulation of cable, m residual of two-dimensional diffusion equation, A / m 2 loop resistance between core and sheath, Q / m dc resistance of conductor, Q/m real part of complex number corrected resolution, pel/in dc resistance of subconductor i, Q,/m self resistance of loop "subconductor i/ring return" (equal to dc internal resistance of subconductor z), Q / m dc resistance of conductor, fi/m ring-return path (Chapter 5); radial direction (Appendix A) radius of boundary, m radius of core conductor, m inner radius of earth-return path (equal to outer radius of outermost insulation of cable), m resolution of image, pel/in resizing factor outside radius of insulation, m inner radius of conductor, m row number of subconductor external radius of conductor, m outer radius of conductor, m maximum number of rows per partition for matrix modification inner radius of conducting sheath, m outer radius of conducting sheath, m area of conductor, m 2 cross sectional area of kth conductor, m 2 solution region scale factor of image (drawing length/actual length) time, s thickness of conductor, m or in unit vector in r direction unit vector in z direction unit vector in <f> direction voltage or electric scalar potential, V voltage between core and sheath, V voltage between sheath and ground, V voltage drop across core with respect to ground (Chapter 3), V; voltage drop across core with respect to ring (Chapter 5), V List of Symbols 128 V c V c V , V s h e a t h X Xi x\ X2 I 2 X. y y' yi Vi yi y'2 z J12 J13 J22 J23 voltage drop across core with respect to ground, V voltage drop across loop "core/sheath", V voltage drop across sheath with respect to ground (Chapter 3), V; voltage drop across sheath with respect to ring (Chapter 5), V voltage drop across sheath with respect to ground, V direction of propagation (Chapter 3); horizontal distance between cables i and k (Chapter 3), m x-coordinate of an arbitrary point in subconductor i (Chapter 4), m x-coordinate of an arbitrary point in subconductor k, m upper side x-coordinate of subconductor i, m upper side x-coordinate of subconductor k, m lower side x-coordinate of subconductor i, m lower side x-coordinate of subconductor k, m y-coordinate of an arbitrary point in subconductor i, m y-coordinate of an arbitrary point in subconductor k, m left-hand side y-coordinate of subconductor i, m left-hand side y-coordinate of subconductor k, m right-hand side y-coordinate of subconductor i, m right-hand side y-coordinate of subconductor k, m complex impedance per unit length, 0,/m self impedance of loop "core/sheath" (Chapter 3), fi/m; self impedance of loop "subconductor 1/ring-return path" (Chapter 5), Q/m; self impedance of loop "sector-shaped conductor 1/sheath" (Chap-ter 6), Q/m mutual impedance between loop "core/sheath" and loop "sheath/ground" (Chapter 3), Q/m; mutual impedance between loop "subconductor 1/ring-return path" and loop "subconductor 2/ring-return path" (Chapter 5), f2/m; mutual impedance between loop "sector-shaped conductor 1/sheath" and loop "sector-shaped conductor 2/sheath" (Chapter 6), 0,/m mutual impedance between loop "subconductor 1/ring-return path" and loop "subconductor 3/ring-return path" (Chapter 5), 0,/m; mutual impedance between loop "sector-shaped conductor 1/sheath" and loop "sector-shaped conductor 3/sheath" (Chapter 6), fi/m self impedance of loop "sheath/ground" (Chapter 3), f2/m; self impedance of loop "subconductor 2/ring-return path" (Chapter 5), fi/m; self impedance of loop "sector-shaped conductor 2/sheath" (Chap-ter 6), fl/m mutual impedance between loop "sheath/armor" and loop "ar-mor/ground" (Chapter 3), fi/m; mutual impedance between loop "subconductor 2/ring-return path" and loop "subconductor 3/ring-return path" (Chapter 5), fi/m; List of Symbols 129 -*33 z be ^core-out Zcore/sheath-Z c s Zdc Zhf Zi Z;i Zik Zinsulation Z mutual insulation Jmutual earth—return -'pipe—in Zself earth Z return sheath/earth—insulation Jsheath—in Jsheath—mutual Jsheath—out mutual impedance between loop "sector-shaped conductor 2/sheath" and loop "sector-shaped conductor 3/sheath" (Chapter 6), Cl/m self impedance of loop "armor/ground" (Chapter 3), Cl/m; self impedance of loop "subconductor 3/ring-return path" (Chapter 5), Cl/m; self impedance of loop "sector-shaped conductor 3/sheath" (Chap-ter 6), Cl/m mutual impedance between sheath a and sheath b with return through earth, Cl/m mutual impedance between sheath a and sheath c with return through earth, Cl/m mutual impedance between sheath b and sheath c with return through earth, Cl/m frequency-dependent impedance of conductor, Cl/m self impedance of core with return through earth (Chapter 3), or with return through ring (Chapter 5), Cl/m loop impedance between core and sheath, Cl/m internal impedance of core conductor with return through outer conducting sheath, Cl/m impedance of insulation between core and sheath, Cl/m mutual impedances between core and sheath with return through earth (Chapter 3), or with return through ring (Chapter 5). Cl/m impedance (resistance) of conductor at dc, Cl/m high frequency impedance of conductor, Cl/m internal impedance, Cl/m self impedance of loop "subconductor i/ring return", Cl/m mutual impedance between loops "conductor i/ring return" and "conductor A;/ring return", Cl/m impedance of insulation between inside eccentric conductor and pipe, Cl/m mutual impedance between inside eccentric conductors (with return through pipe), Cl/m mutual impedance of earth-return path (between outermost loop "sheath/earth" of one single-circuit and outermost loop "sheath/earth" of another single-circuit), Cl/m internal impedance of pipe with return through inside eccentric conductor, Cl/m self impedance of earth-return path, Cl/m impedance of insulation between sheath and earth-return path, Cl/m internal impedance of conducting sheath with return through inner core conductor, Cl/m mutual impedance of conducting sheath (between loop "core/sheath" and loop "sheath/earth" of the same single-circuit), Cl/m internal impedance of conducting sheath with return through outer List of Symbols 1 3 0 earth-return path, 0,/m Z s s self impedance of sheath with return through earth (Chapter 3), or with return through ring (Chapter 5), Cl/m Ztube-in internal impedance of tubular conductor with return path inside the tube, Q/m Ztube-mutuai mutual impedance of tubular sheath (between the inside loop and the outside loop of the conducting sheath), Q/m •• Ztube-out internal impedance of tubular conductor with return path outside the tube, fi/m z direction of propagation (Appendix A); Greek Symbols 7 Euler's constant A R resolution step, pel/in AResolution resolution step, pel/in 5 penetration depth, m <5m;n minimum penetration depth, m or percentage of conductor thick-ness e permittivity, F / m e0 permittivity of vacuum, F / m er relative permittivity er relative error p. permeability, H / m fj,r relative permeability Ho permeability of vacuum, H / m p resistivity, • m a conductivity, 13/m or S/m r distance between an arbitrary point in subconductor i and an ar-bitrary point in subconductor k, m 4> azimuthal direction <Pi base function of A for node i to angular frequency, rad/s Acronyms A C S R Aluminum Conductors, Steel Reinforced A I E E American Institute of Electrical Engineers ac alternating current C A D Computer-Aided Design C I N V E S T A V "Centro de Investigation y Estudios Avanzados del Instituto Politec-nico National", (Center of Research and Advanced Studies of the National Polytechnic Institute) CONICIT "Consejo National para el Desarrollo de las Investigaciones Cienti-ficas y Tecnologicas", (National Council for the Development of Scientific and Technological Research) C P U Central Processing Unit List of Symbols 131 C T E "Departamento de Conversion y Transporte de Energia" (Depart-ment of Energy Conversion and Delivery) dc direct current D D P "Direction de Desarrollo Profesoral" (Faculty Development Office) E M T P Electromagnetic Transients Program [17, 19] E P R Ethylene-Propylene Rubber EPRI Electric Power Research Institute F E Finite Element G M D Geometric Mean Distance H P F F High-Pressure Fluid-Filled I E E E Institute of Electrical and Electronics Engineers P C Personal Computer P E Polyethylene P E E C Partial Element Equivalent Circuit P S E C Partial Subconductor Equivalent Circuit P T Pipe-Type R A M Random Access Memory R D Rediscretization S C L F Self-Contained Liquid-Filled S F 6 Sulfur Hexafluoride T E M Transverse Electromagnetic U B C University of British Columbia USB "Universidad Simon Bolivar", (Simon Bolivar University) VLSI Very Large Scale Integrated X L P E Crosslinked Polyethylene Bibl iography [1] A I E E Committee Report, "A-C resistance of pipe-cable systems with segmental conductors," AIEE Transactions, vol. 71, pt. Ill, pp. 393-414, January 1952. [2] Alden, B., "Graphics software," The Institute, a news supplement to IEEE, column: Traveling the Information Highway, February 2000. [3] Ametani, A. , "A general formulation of impedance and admittance of cables," IEEE Transactions on Power Apparatus and Systems, vol. PAS-99, no. 3, pp. 902-910, May/June 1980. [4] Ametani, A . and K. Fuse, "Approximate method for calculating the impedances of multiconductors with cross section of arbitrary shapes," Electrical Engineering in Japan, vol. 112, no. 2, pp. 117-123, 1992. [5] Anders, J . A. , Rating of Electric Power Cables: Ampacity Computations for Trans-mission, Distribution, and Industrial Applications, I E E E Press Power Engineering Series/McGraw-Hill, New York, NY, 1997. [6] Antonini, G . , A . Orlandi, and C. R. Paul, "Internal impedance of conductors of rectangular cross section," IEEE Transactions on Microwave Theory and Techniques, vol. 47, no. 7, pp. 979-985, 1999. [7] Arizon, P. de and H. W. Dommel, "Computation of cable impedances based on subdivision of conductors," IEEE Transactions on Power Delivery, vol. PWRD-2, pp. 21-27, January 1987. [8] Bianchi, G . and G . Luoni, "Induced currents and losses in single-core submarine cables," IEEE Transactions on Power Apparatus and Systems, vol. PAS-95, no. 1, pp. 49-58, January/February 1976. [9] Cao, M . , P. P. Biringer, and V . Bulat, "Eddy current losses and unbalance in electrode arms," IEEE Transactions on Magnetics, vol. 24, no. 6, pp. 2670-2672, 1988. [10] Carson, J . R., "Wave propagation in overhead wires with ground return," Bell System Technical Journal, vol. 5, pp. 539-554, 1926. [11] Castellanos, F . , J . R. Marti, and F. Marcano, "Phase-domain multiphase transmission line models," Electrical Power & Energy Systems, vol. 19, no. 4, pp. 241-248, 1997. 132 B I B L I O G R A P H Y 133 Cheng, D. K. , Field and Wave Electromagnetics, Addison-Wesley Publishing Com-pany Inc., second edition, 1989. Comellini, E . , A . Invernizzi, and G . Manzoni, "A computer program for determining electrical resistance and reactance of any transmission line," IEEE Transactions on Power Apparatus and Systems, vol. PAS-92, no. 1, pp. 308-314, January/February 1973. Cristina, S. and M . Feliziani, "A finite element technique for multiconductor cable parameters calculation," IEEE Transactions on Magnetics, vol. 25, no. 4, pp. 2986-2988, 1989. Deri, A . , G . Tevan, A . Semlyen, and A. Castanheira, "The complex ground return plane - a simplified model for homogeneous and multi-layer earth return," IEEE Transactions on Power Apparatus and Systems, vol. PAS-100, no. 8, pp. 3686-3693, August 1981. Digital Equipment Corporation, DIGITAL Fortran, Language Reference Manual, Digital Equipment Corporation, Maynard, M A , April 1997. Dommel, H. W. , "Digital computer solution of electromagnetic transients in single-and multi-phase networks," IEEE Transactions on Power Apparatus and Systems, vol. PAS-88, no. 4, pp. 388-399, April 1969. Dommel, H. W., "Overhead line parameters from handbook formulas and computer programs," IEEE Transactions on Power Apparatus and Systems, vol. PAS-104, no. 2, pp. 366-372, February 1985. Dommel, H. W., EMTP Theory Book, Microtran Power System Analysis Corporation, Vancouver, B C , second edition, May 1992. Latest update: April 1996. Dommel, I. I., "Program documentation for subroutine T U B E , " U B C Internal hand-written report, May 1978. EPRI, Underground Transmission Systems Reference Book, Electric Power Research Institute, TR-101670, 1992. Ferkal, K . , M . Poloujadoff, and E . Dorison, "Proximity effect and eddy current losses in insulated cables," IEEE Transactions on Power Delivery, vol. 11, no. 3, pp. 1171— 1178, 1996. Fink, D. G . and H. W Beaty, Standard Handbook for Electrical Engineers, McGraw-Hill, Inc., New York, NY, 14th edition, September 1999. Gonzalez, R. and R. Woods, Digital Image Processing, Addison-Wesley Publishing Company Inc., 1993. Graneau, P., "Alternating and transient conduction currents in straight conductors of any cross-section," International Journal Electronics, pp. 41-59, 1965. B I B L I O G R A P H Y 134 [26] Graneau, P., Underground Power Transmission, John Wiley & Sons, 1979. [27] Grover, F . W., Inductance Calculations: Working Formulas and Tables, D. Van Nostrand Company, Inc, New York, NY, second edition, May 1947. [28] Gustavsen, B., J . Sletbak, and T. Henriksen, "Calculation of electromagnetic tran-sients in transmission cables and lines taking frequency dependent effects accurately into account," IEEE Transactions on Power Delivery, vol. 10, no. 2, pp. 1076-1084, 1995. [29] Hesse, M . H . , "Electromagnetic and electrostatic transmission-line parameters by digital computer," IEEE Transactions on Power Apparatus and Systems, vol. 82, pp. 282-291, June 1963. [30] Hill, R. J . , S. Brillante, and P. J. Leonard, "Railway track transmission line parame-ters from finite element field modelling: series impedance," IEE Proceedings-Electric Power Applications, vol. 146, no. 6, pp. 647-660, 1999. [31] Huang, C . - C , "Computation of resistance and inductance matrices in a symmet-ric structure," IEEE Transactions on Components, Packaging, and Manufacturing Technology, vol. 18, part A, no. 3, pp. 674-676, 1995. [32] Hwang, L . - T . and I. Turlik, "A review of the skin effect as applied to thin film interconnections," IEEE Transactions on Components, Hybrids, and Manufacturing Technology, vol. 15, no. 1, pp. 43-55, 1992. [33] Kane, M . , A . Ahmad, and P. Auriol, "Multiwire shielded cable parameter computa-tion," IEEE Transactions on Magnetics, vol. 31, no. 3, pp. 1646-1649, 1995. [34] Kimbark, E . W. , Electrical Transmission of Power and Signals, John Wiley & Sons, Inc., New York, NY, sixth edition, November 1962. [35] Laplante, P. A. (Editor), Electrical Engineering Dictionary, C R C Press L L C , 2000. [36] Lewis, W. A. and P. D. Tuttle, "The resistance and reactance of aluminum conductors, steel reinforced," AIEE Transactions, vol. 78, pt. Ill, pp. 1189-1215, February 1959. [37] Lucas, R. and S. Talukdar, "Advances in finite element techniques for calculating cable resistances and inductances," IEEE Transactions on Power Apparatus and Systems, vol. PAS-97, no. 3, pp. 875-883, May/June 1978. [38] Marti, J. R., "Accurate modelling of frequency-dependent transmission lines in elec-tromagnetic transient simulations," IEEE Transactions on Power Apparatus and Systems, vol. PAS-101, no. 1, pp. 147-157, January 1982. [39] Marti, L . , "Simulation of transients in underground cables with frequency-dependent modal transformation matrices," IEEE Transactions on Power Delivery, vol. 3, no. 3, pp. 1099-1110, July 1988. B I B L I O G R A P H Y 135 MathWorks Inc., The, Using Matlab, The MathWorks, Inc, Natick, M A , January 1999, fourth printing, revised for M A T L A B 5.3 (Release 11). MathWorks Inc., The, Image Processing Toolbox User's Guide, The MathWorks, Inc, Natick, M A , January 1999, revised for Version 2.2 (Release 11) (Online only). Meyerhoff, L . and G . S. Eager, Jr., "A-C resistance of segmental cables in steel pipe," AIEE Transactions, vol. 68, pt. II, pp. 816-834, 1949. Microsoft Corporation, Developer Studio Environment User's Guide, Microsoft Cor-poration, Redmond, WA, 1994-1997. Microtran Power System Analysis Corporation, MicroTran Reference Manual, Mi -crotran Power System Analysis Corporation, Vancouver, B C , 1997. Murray, J. D., "Graphics file formats F A Q (part 1 of 4): General graphics for-mat questions," http://www.faqs.org/faqs/graphics/fileformats-faq/partl, January 20 1997, last update September 25, 1998. Nguyen, T . T . , "Earth-return path impedances of underground cables. Part 1: Nu-merical integration of infinite integrals," IEE Proceedings-Generation, Transmission and Distribution, vol. 145, no. 6, pp. 621-626, 1998. Noda, T . , N. Nagaoka, and A. Ametani, "Phase domain modeling of frequency-dependent transmission lines by means of an A R M A model," IEEE Transaction on Power Delivery, vol. 11, no. 1, pp. 401-411, January 1996. Oeding, P. and K. Feser, "Geometric mean distances of rectangular conductors (in german)," ETZ-A, vol. 86, no. 16, pp. 525-533, 1965. Ollivier-Gooch, C , " G R U M M P - Generation and refinement of unstructured, mixed-element meshes in parallel," http://tetra.mech.ubc.ca/GRUMMP/, September 02 1999. Owen, S., "Meshing software survey," http://www.andrew.cmu.edu/user/sowen/ softsurv.html, September 1998. Papagiannis, G . K. , D. G. Triantafyllidis, and D. P. Labridis, "A one-step finite element formulation for the modeling of single and double-circuit transmission lines," IEEE Transactions on Power Systems, vol. 15, no. 1, pp. 33-38, February 2000. Parker, S. P., McGraw-Hill Dictionary of Scientific and Technical Terms, McGraw-Hill, Inc., New York, NY, fifth edition, 1994. Parsi-Feraidoonian, R. and H. W. Dommel, "Impedances of distribution cables with semi-conductive layers," in Canadian Electrical Association, Analytical Techniques of Insulation Co-ordination Committee, Power System Planning and Operation Section, Enginnering and Operation Division, Vancouver, B C , April 1992, pp. 1-11. B I B L I O G R A P H Y 136 [54] Pollaczek, F . , "Uber das feld einer unendlich langen wechsel stromdurchflossenen einfachleitung ("On the field produced by an infinitely long wire carrying alternating current", in german)," E.N.T., vol. 3, no. 9, pp. 339-360, 1926. [55] Ramo, S., J . R. Whinnery, and T. Van Duzer, Fields and Waves in Communication Electronics, John Wiley & Sons, Inc., New York, N Y , third edition, 1994. [56] Restle, P., A. Ruehli, and S. G. Walker, "Dealing with inductance in high-speed chip design," Proceedings 1999 Design Automation Conference, pp. 904-909, 1999. [57] Rivas, R. A. , Overhead Transmission Lines and Underground Cables, Section 9 of The Handbook for Electric Power Calculations, Wayne Beaty (Editor), McGraw-Hill, Inc., New York, NY, 3rd edition, 2000. [58] Rivas, R. A . and J. R. Marti, "Calculation of frequency-dependent parameters of power cable arrangements using pixel-shaped conductor subdivisions," in Proceed-ings of International Conference on Power Systems Transients (IPST'99), Budapest, Hungary, June 1999, pp. 335-340. [59] Rodolakis, A . J . , "Point-and-click cable ampacity studies," IEEE Computer Appli-cations in Power, vol. 11, no. 2, pp. 53-56, 1998. [60] Saad, O., G . Gaba, and M . Giroux, "A closed-form approximation for ground return impedance of underground cables," IEEE Transactions on Power Delivery, vol. 11, no. 3, pp. 1536-1545, 1996. [61] Satsios, K. J . , D. P. Labridis, and P. S. Dokopoulos, "Finite element computation of field and eddy currents of a system consisting of a power transmission line above conductors buried in nonhomogeneous earth," IEEE Transactions on Power Delivery, vol. 13, no. 3, pp. 876-882, 1998. [62] Satsios, K. J . , D. P. Labridis, and P. S. Dokopoulos, "Inductive interference caused to telecommunication cables by nearby A C electric traction lines, measurements and fern calculations," IEEE Transactions on Power Delivery, vol. 14, no. 2, pp. 588-594, 1999. [63] Schelkunoff, S. A. , "The electromagnetic theory of coaxial transmission lines and cylindrical shields," Bell System Technical Journal, vol. 13, pp. 532-579, 1934. [64] Silvester, P., "The accurate calculation of skin effect in conductors of complicated shape," IEEE Transactions on Power Apparatus and Systems, vol. PAS-87, no. 3, pp. 735-742, March 1968. [65] Smith, D. R. and J. V . Barger, "Impedance and circulating current calculations for the U D multi-wire concentric neutral circuits," IEEE Transactions on Power Apparatus and Systems, vol. PAS-91, no. 3, pp. 992-1006, May-June 1972. B I B L I O G R A P H Y 137 [66] Triantafyllidis, D. G . , G . K. Papagiannis, and D. P. Labridis, "Calculation of over-head transmission line impedances a finite element approach," IEEE Transactions on Power Delivery, vol. 14, no. 1, pp. 287-293, 1999. [67] Tsuk, M . J. and J. A. Kong, "A hybrid method for the calculation of the resistance and inductance of transmission lines with arbitrary cross sections," IEEE Transactions on Microwave Theory and Techniques, vol. 39, no. 8, pp. 1338-1347, 1991. [68] Uribe, F . A . , J . L . Naredo, P. Moreno, and Guardado L . , "Calculating earth impedances for underground transmission cables," in Proceedings of International Conference on Power Systems Transients (IPST'01), Rio de Janeiro, Brazil, June 2001. [69] Wang, Y . - J . , "A precise method for the impedance calculation of a power rail tak-ing into account the skin effect and complex geometry," European Transactions on Electrical Power, vol. 10, no. 1, pp. 19-27, January/February 2000. [70] Wang, Y . - J . and J . -H. Wang, "Modeling of frequency-dependent impedance of the third rail used in traction power systems," IEEE transactions on Power Delivery, vol. 15, no. 2, pp. 750-755, April 2000. [71] Wedepohl, L . M . and D. J. Wilcox, "Transient analysis of underground power-transmission systems - system-model and wave-propagation characteristics," Proceed-ings of the Institution of Electrical Engineers, vol. 120, no. 2, pp. 253-260, February 1973. [72] Weeks, W . T . , L . L . Wu, M . F. McAllister, and A. Singh, "Resistive and inductive skin effect in rectangular conductors," IBM Journal of Research and Development, vol. 23, no. 6, pp. 652-660, November 1979. [73] Wiseman, R. J . , "A-C resistance of large size conductors in steel pipe or conduit," AIEE Transactions, vol. 67, pt. II, pp. 1745-1758, 1948. [74] Xie, X . and J .L . Prince, "Frequency response characteristics of reference plane effec-tive inductance and resistance," IEEE Transactions on Advanced Packaging, vol. 22, no. 2, pp. 221-229, May 1999. [75] Yin, Y . , Calculation of Frequency-Dependent Parameters of Underground Power Ca-bles with Finite Element Method, Ph.D. thesis, University of British Columbia, Van-couver, B C , September 1990. [76] Yin , Y . and H. W. Dommel, "Calculation of frequency-dependent impedances of underground power cables with finite element method," IEEE Transactions on Mag-netics, vol. 25, no. 4, pp. 3025-3027, July 1989. [77] Zhou, D. H . , Numerical Methods for Frequency-Dependent Line Parameters with Applications to Microstrip Lines and Pipe-Type Cables, Ph.D. thesis, University of British Columbia, Vancouver, B C , 1993. B I B L I O G R A P H Y 138 [78] Zhou, D. H . and J. R. Marti, "Skin effect calculations in pipe-type cables using a linear current subconductor technique," IEEE Transactions on Power Delivery, vol. 9, no. 1, pp. 598-604, 1994. Append ix A Derivat ion of the Frequency-Dependent Parameters of a Cy l ind r i ca l Conductor To determine the frequency-dependent parameters (resistance and inductance) of one cylin-drical, solid conductor, a Bessel's differential equation can be derived using Maxwell's equations. Such a Bessel's equation will describe either current density or impedance, and its solution will be given by Bessel functions. The theory employed was taken from [55]. To start, let us consider Faraday's law and the generalized Ampere's law, both in differ-ential form and with sinusoidal time variation V x E = - — = -jtvaH (A.l) V x H = J + — = J + jueE (A.2) Applying the curl to (A.l) , the following equation is obtained V x V x E = -juuV x H (A.3) Including (A.2) in (A.3) yields V x V x E = -JLUU-{J + JUJCE) (A.4) Assuming the presence of a good conductor, the expression J = crE (Ohm's law) can be included in (A.4) 139 140 V x V x E = -ju!u(o + j we)E (A.5) Since a » we in good conductors, displacement currents can be neglected, and (A.5) becomes V x V x E = -jtuno-E (A.6) Applying Ohm's law again, (A.6) can be written as V x V x J = -JLjuoJ (A.7) For vector identity V x V x J = V ( V • J) — V 2 J , so (A.7) becomes V ( V • J) - V 2 J = -juju-aJ (A.8) The assumption of Ohm's law implies the absence of charge density (the net charge is zero in a good conductor). Therefore V - D = p = 0 ^ V - E = 0 ^ V - J = 0 (A.9) and (A.8) becomes V 2 J = jujuo-J (A.10) Using the formula of the Laplacian operator in cylindrical coordinates, and assuming only a current density component along the z direction (propagation direction), the following equation can be written 7 2 T 1 d (dJ,& \ ldJ,A , d2Jz r or \ or J r or or* Matching (A.10) and (A.11), a Bessel's differential equation is finally obtained 5# + - ^ + T 2 J 2 = 0 (A.12) or* r or 141 where T 2 = -jufia (A.13) Since 6 = * I —-—, T can also be written as (A.14) where 5 is the "skin depth" or "depth of penetration", which is the depth at which fields and currents have been decreased by a factor equal to e _ 1 (0.368), or 63.2 percent with respect to the values at the surface. The two independent solutions of the Bessel's differential equation are given by Jz = AJ0(Tr) + BH0M{Tr) (A.15) For a solid wire, the case of r = 0 must be included in the solution; since BH0{x\Tr) tends to infinity when Tr tends to zero, the constant B must be equal to zero. Therefore Jz = AJ0(Tr) (A.16) The constant A can be found since J z = oE0 at r„, where r0 is the external radius of the conductor. Therefore ' • = ' * $ £ j < A- 1 7> Using (A.14), (A.17) can be written as follows Jz = oE0 (A.18) J0{r*Y r0) Dividing J0 into real and imaginary parts, the following expressions can be defined Ber{v) — real part of J0(j ^v) (A.19) 142 Bei(v) = imag part of Ja(j ?v) (A.20) Jo(j 2v) = Ber(v) + jBei(v) (A.21) j l2J'0{j 2 « ) = Ber'(v) + jBei'(v) and the current density can be written as Jz = aE0 Ber(%)+jBei{%) Ber(^r0)+jBei(^r0) (A.22) (A.23) where Ber(v) and Bei(v) are Bessel functions that can be calculated with polynomial approximations. To find the internal impedance, let us recall that Ez{rQ) E0 Zi I I (A.24) where the current I can be obtained either by integrating J or by applying Ampere's law; from Ampere's law 2Trr0Htj)\r=ro = I Hj, as a function of E can be obtained from (A.l); therefore (A.25) r 1 r d_ dr a dd, a dz Er rE, E, - -ju/j.(Hrur + HiMi, + Hzuz) (A.26) Since E = Ezuz and there are only variations of Ez with respect to r 143 From (A. 17) Deriving (A.28) Matching (A.27) and (A.29) Using (A. 13), (A.30) becomes Including (A.31) in (A.25) ^ = junH* (A.27) a Jo{Tr0) dEz = E0J'0{Tr)T Or J0{Tr0) (A.29) * = ^ £ (A.30) ju>u-J0{Tr0) H*--^j7TrJ (A-31) / = _ ^ ^ 1 ( A , 2 ) T J0(Tr0) Including (A.32) in (A.24), an expression for the internal impedance is finally obtained T J o ^ T r o ) (A.33) 2nr0a J'0(Tr0) Using (A.21), (A.22), and (A.14), (A.33) becomes Y 2 1 Ber(^r0)+jBei(^r0) Zi = j v \ j "/ A •»/ (A.34) 6 2ixr0a Ber'{^r0)+jBei'(^r0) 144 Defining the quantities Rs and q Rs = ~~~r (A.35) oo ~ r 0 (A.36) the internal impedance (A.34) can be written as = • Rs Ber(q)+jBei(q) V27rr0Ber'(q)+jBei'(q) V ' Using modified Bessel functions with complex argument, the expression for the internal impedance becomes z_ = VIRa iojqVJ) , A 38^ V2nr0 h {qy/j) where UQy/j) = Ber{q) + jBei{q) (A.39) VJh (qy/1) = Ber'(l) + JBei!(q) (A.40) and I0 and Ii are modified Bessel functions of the 1 s t kind, order 0 and 1, respectively. The expression (A.38) is equal to the expression (3.7) on page 32. To prove it, let us recall that Rs = l/(<7<5), a = 1/p, 8 = VWi^i^) a n d m = \Jiw\xo. Introducing the above relations into (A.38) we obtain Since q = y/2T0/8 and 8 = Vfilmi f n e argument qV] can be written as m.r0 and (A.41) becomes 145 pm I0(mr0) 2nr0 J i (mr0) (A.42) which is (3.7) on page 32 when r0 = r c o r e . Using M A T L A B ® [40], a program was written to evaluate the internal impedance given by (A.38). The same case studied in [55], page 183, was used. The data of the conductor is as follows: 1-mm diameter Copper wire, a — 5.80 • 107 S/m @ 300K. The dc resistance and the dc inductance were calculated as RQ = l/(airrl) and L0 = p,0/(8ir), respectively. As can be seen, the same curves indicated in [55] (page 185) were obtained. S o l i d - w i r e sk in effect quanti t ies c o m p a r e d with d c v a l u e s 3.5 2.5 1.5 L/Lo 2 3 4 5 ro/del ta , Rat io of rad ius to depth of penetrat ion Figure A . l : Solid-wire skin effect quantities compared with dc values. Append ix B Wedepohl 's Formulae for the Frequency-Dependent Parameters of Tubular Conductors Wedepohl and Wilcox [71] propose the use of formulae based on hyperbolic functions for the calculation of the frequency-dependent parameters in systems with tubular conductors (e.g., the system depicted in Fig. 3.1, page 28). The formula for Z C O re-out was derived com-paring the series expansion of the classical formula (equation (3.7) on page 32) with that of a preliminary formula based on a hyperbolic cotangent. The coefficients of the series expansion of the preliminary formula were then adjusted to predict parameters accurate to two decimal places. The other formulae were derived assuming that the conducting sheath is usually much thinner than its average radius. As a result of such an assumption, a Bessel's second order differential equation becomes a second order differential equation with constant coefficients. The approximate formulae for the other impedances come from the solution of such a differential equation with constant coefficients. Zcore-out = pcZTcore coth(0.777m c o r e r c o r e ) + 0-356af">" in Q/m. This formula yields a max-imum error of four percent in the resistive part (it occurs when | m r c o r e | = 5) and a maximum error of five percent in the reactive component (it occurs when | m r c o r e | = 3.5 [71]). For other values of | m r c o r e | , the formula is very accurate and avoids the evaluation of Bessel functions; rcore is the radius of the core conductor, pcore is the resistivity of the core conductor in fi-m, and m c o r e is the reciprocal of the complex penetration depth of the core; m c o r e — \J(j ujpcore)/' pcore in m _ 1 , where pcore is the magnetic permeability of the core in H / m , \icore = u-Tcore • //„ with p,rcore ^ 1 if the material of the core conductor is magnetic, and w is the angular frequency in rad/s. Z Sheath-in = 2Jshu . { m s h c o t h ( m s h A s / t ) - 1 } in Q/m. This formula yields a good'accuracy if the condition r ' h ~ o u t ~ r ^ ~ i n < \ is satisfied [71]; rsit-oui is the outer ra-T s h — out~T~Ts}l — i n o dius of the sheath, rsh_in is the inner radius of the sheath, psh is the resistivity of the 146 147 sheath in fl-m, and m s h is the reciprocal of the complex penetration depth of the sheath; m s h = V(iu>Hsh)/Psh in n i _ 1 , where /ish is the magnetic permeability of the sheath in H / m , nSh — HrSh ' Mo with jj,rsh ^ 1 if the material of the sheath is magnetic, and A s / j is the thickness of the sheath, which can be calculated as (rsh-0ut — Tsh-in)-ZSheath-out = ^Jf—{mshcoth(mshA.s/l) + r ^ * } in fi/m. This formula yields a good accuracy if the condition between radii mentioned for Z s n ea th - in is satisfied. Zsheath-mutuai = n ( r , p'.'l+r*\—?> cosech(m s h A s h ) in ft/m. This formula yields a good ac-" V s h —in* ' sh — o u t ) curacy if the condition between radii mentioned for Z s h e ath - in and Z S h e a t h - o u t is satisfied. If the single-circuits have more concentric conductors, e.g., armors, the formulae for the new impedances can be derived by analogy utilizing the right electric properties as well as the appropriate radii. Append ix C Capacit ive Susceptances Associated wi th Underground Tubular Conductors Consider again the underground cable system depicted in Fig. 3.1, page 28. Determine the capacitive susceptances of such a system. 1. Calculate the self and mutual susceptances The method and the formulae proposed in [19] and [71] are also appropriate to calculate the shunt susceptances of this cable system. By assuming that there is no capacitive coupling among the three phases because of shielding effects, the following six nodal equations can be written Ulcere* /9x ' Bcca BCSa 0 0 0 0 " v v corea d l s h e a t h a / ^ B C S a BsSa 0 0 0 0 V s heat h a ^Icoreb / dx 1 0 0 BCcb Bcsb 0 0 v v corej, dlsheathb/dx ~ J 0 0 BCsb BSsb 0 0 V sheathb <9Icorec /dx 0 0 0 0 " BCCc Bcsc v v corec dlsheathjdx _ 0 0 0 0 BCSc BSSc . Vsheath c where a,b,c = subscripts denoting quantities associated with the single-circuits of phases a, b, and c, respectively, Icore = charging current flowing through core conductor, Isheath = charging current flowing through sheath, V C o r e = voltage of core conductor with respect to ground, V s heath = voltage of sheath with respect to ground, Bcc = self shunt susceptance per unit length of core conductor, S c., = mutual shunt susceptance per unit length between core conductor and sheath, 148 149 Bss = self shunt susceptance per unit length of sheath, with Bcc = Bi, Bcs = —Bi, and Bss = B\ + B2, B\ = capacitive susceptance per unit length of insulation layer between core and sheath, B2 = capacitive susceptance per unit length of insulation layer between sheath and earth, Bi = LoCi in cJ/m and C; = 2ire0€ri/\n(ri/qi), where C{ is the shunt capacitance of the tubular insulation in F / m , ^ is the inside radius of the insulation, is the outside radius of the insulation, eri is the relative permittivity of the insulation material, and e0 is the permittivity of free space in F / m . 2. Eliminate the sheaths The Matrix of self and mutual susceptances can be reduced by applying a Kron's reduction. To do so, use the same procedure explained in Subsection 3.3, page 34. Append ix D Fini te Element M e t h o d D . l Ca lcu la t ion Procedure This section summarizes the calculation procedure proposed by Yin in his Ph.D. thesis work [75]. 1. Obtain the basic equations of the solution region: The principal equations of the solution region are a two-dimensional diffusion equation and the equations relating the magnetic vector potential with the current density of each conductor u.-1V2A-JLuaA + Js = 0 in SR (D.l) - j w J aAds + SChJSk = h (A: = 1,2,... ,K) (D.2) sck where // and a are permeability and conductivity, co is angular frequency, A is the longi-tudinal component of the magnetic vector potential (its other components are zero), SR is the solution region, K is the number of conductors, and Sck a n d Jsk are the cross sectional area and the source current density of the kth conductor, respectively. J$ is constant over the cross section of conductor and can be visualized as the current distribution that exists in the absence of time variation, i.e., uniform current distribution near dc. Equations (D.l) and (D.2) are derived from Maxwell's equations in Section D.2. 2. Assign base functions to the field distribution A: The base functions are usually polynomials since they can be differentiated and inte-grated straightforwardly, and their selection depends on the type of shape assigned to the elements. The base functions approximate the field distribution by the following interpo-lation formula A = Y^Ancpn (D.3) n=l 150 D . l . Calculation Procedure 151 where An is the value of the field distribution A at the node n, (pi,if2,--- ,VNT are the base functions of A, and NT is the total number of nodes. NT = N + NB, where N is the number of nodes with unknown values and NB is the number of boundary nodes with known values (A — 0 at the boundary nodes). The field values at arbitrary points in the solution region can be obtained from (D.3) once the node values are known. As cables usually display circular geometries, the solution is obtained more efficiently using curved-sided elements (e.g., cylindrical shells or wedge-shaped subregions). Yin proved this by comparing for the same discretization error the results of meshes of straight-sided elements (e.g., rectangles or triangles) with those of meshes of curved-sided elements. The former yielded more elements, more nodes and smaller span angles. Therefore, their C P U time and storage requirements were higher. A n exception to the above rule is the case of power rails of transit systems, where straight-sided elements represent the geometry more appropriately. Since the interpolation polynomials are simpler in straight-sided elements, Y in uses a non-linear transformation to map straight-sided elements in the local (auxiliary) coordinate system into curved-sided elements in the global (actual) coordinate system. 3. Approximate the basic equations to a residual: Equation (D.l) can be approximated to a residual R(A), and such a residual can be forced to satisfy the integral equation of the Galerkin technique R(A) = n~1V2A-juaA + Js (D.4) JR{A)<pmds = 0 (m = 1,2,... ,N) (D.5) SR Introducing (D.3) and (D.4) into (D.5), applying Green's formula V • {vVu) = Vv • V u + vV2u to the result, and introducing (D.3) into (D.2), expressions (D.l) and (D.2) become - ^ A „ / {^lV(pm-Vipn+$wo-ipm<Pn)ds+ I Js<pmds = 0 (m = 1, 2 , . . . , N) (D.6) n=l J J SR SR - j u) An I / aipnds I + SCk n=l \J J JSk=Ik {k = l,2,...,K) (D.7) Equations (D.6) and (D.7) can be rewritten as a set of (N + K) algebraic complex equations r to] - [F] i [ [Ao] 1 [I E 1 ] . [FH] [S C ] . . [Js] . . [I] + [IE 2] . (D.8) D.2. Derivation of the Basic Formulae 152 where the unknowns are the field values at the nodes of the mesh [An] and the source current densities of the conductors [J s], and the forcing functions are the currents of the conductors [I] and the known field values at the boundary nodes [As]. The nodes are sorted so that the first nodes are those with unknown values and the last nodes are those of known values (boundary nodes). As a result, the boundary nodes can be moved to the right-hand side of the set of equations by matrix operations and the subvectors of equivalent currents [IEJ] and [IE2] can be obtained. The formulae for the calculation of the above-mentioned arrays are in [75]. The matrix [D] is sparse and banded, and can be factorized using Choleski decompo-sition algorithms. Once [D] is factorized the system becomes [Du] [0] [F'] [sy [Au] [J.] [ I E J I (D.9) where the reduced system of the lower part can be solved to obtain [J s]. Subsequently, the upper part is solved to obtain [Au]. Because the impedances can be derived from [Js] and [I], the solution of the upper part can be skipped if the field values are not to be plotted. 4. Calculate the impedances from the field solution: The impedances can be obtained directly from [Js] with the Js method, or from the power losses and the stored magnetic energy with the loss-energy method. The Js method is usually preferred since it avoids the calculation of all the field variables. With the Js method the impedances of the kih column of the impedance matrix can be calculated from Ji ik Si ho~i (i = 1,2,... ,K) (D.10) forcing U — 0 and Ik ^ 0, V i, k = 1, 2, • • • , K; k. To find all the elements of the impedance matrix, (D.8) is solved K times. However, the system is factorized only once since the solution method only modifies the vector of currents [I]. D.2 Der ivat ion of the Basic Formulae Since the magnetic flux density B is a solenoidal field, i.e., V • B = 0, B can be expressed as the curl of a vector function A , where A is the magnetic vector potential. B - V x A From Faraday's law the electric field intensity E can be written as (D. l l ) (D.12) D.2. Derivation of the Basic Formulae 153 Introducing (D. l l ) into (D.12) „ „ 3 ( V x A ) at or v * ( E + - ) = o If a vector field is curl-free, then it can be expressed as the gradient of a scalar Hence E + ^ = - V V (D.15) at where V is the electric scalar potential and (—VV) is the voltage drop per unit length along a conductor; (—VV) is also called the source electric field E s [75]. E and A have only longitudinal components, so they can be expressed as Euz and A u z , respectively. Assuming also sinusoidal time variation (D.15) becomes (E + jojA)% = — V V (D.16) E s = Esuz can be associated with a source current density Js = J u z , which is a constant over the cross section of conductor [75]. Therefore Jsuz = aEsuz = - a V V (D.17) where a is the conductivity. The volume current density vector J has only longitudinal component, so it can also be expressed as J u z = aEuz. Using this relation and (D.17), equation (D.16) becomes - + juA ) u z a Js&z (D.18) or -juoA + Js (D.19) D.2. Derivation of the Basic Formulae 154 The integration of (D.19) over the cross section of the kth conductor gives h= j Jkds = -ju J oAds + SCkJSk (D.20) sck sck which is the equation (D.2) on page 150. From Ampere's law and neglecting displacement currents V x B = |iiJ (D.21) where /v, is the permeability. Including (D. l l ) in (D.21) V x V x A = AiJ (D.22) Since V x V x A = V ( V • A ) - V 2 A , (D.22) becomes V ( V • A ) - V 2 A = HJ (D.23) Making V • A = 0 for Coulomb condition - / x " 1 V 2 A - J (D.24) Substituting A for Auz and J for J u z in (D.24) and introducing (D.19) into the result H~1V2A-juaA + Js = 0 (D.25) which is the two-dimensional diffusion equation (D.l) on page 150. 

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

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

Comment

Related Items