- Email: [email protected]
Wire antenna versus modified transmission line approach to the transient analysis of grounding grid
Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
Contents lists available at ScienceDirect
Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound
Wire antenna versus modified transmission line approach to the transient analysis of grounding grid Damir Cavka a, Basma Harrat b, Dragan Poljak a,n, Bachir Nekhoul b, Kamal Kerroum c, Khalil El Khamlichi Drissi c a
Department of Electronics, University of Split, Croatia LAMEL Laboratory, University of Jijel, Algeria c LASMEA Laboratory, University of Blaise Pascal, France b
a r t i c l e i n f o
abstract
Article history: Received 1 October 2010 Accepted 13 April 2011 Available online 7 June 2011
The paper deals with transient analysis of grounding grids using two different approaches, wire antenna theory and modified transmission line model. The Pocklington integro-differential equations, in frequency domain, arising from the wire antenna theory are numerically handled via the Galerkin– Bubnov variant of indirect Boundary Element Method (GB-IBEM), while the transient response was obtained using inverse Fourier transform. The modified transmission line equations are treated using the finite difference time domain (FDTD) method. Some illustrative numerical results are presented and discussed in the paper. & 2011 Elsevier Ltd. All rights reserved.
1. Introduction Transient analysis of grounding systems is of great importance in lightning protection systems (LPS) design. One of the most important parameters arising from the transient analysis is the transient impedance of the grounding electrode. Grounding systems can be modeled using simple electric circuit methods [1,2], the transmission line based model (TLM) [3–5] or antenna (full-wave) model (AM) [6–9]. While the circuit methods can be considered to be rather oversimplified, the TLM methods have advantage of simplicity and relatively low computational cost. On the other hand, though valid for long horizontal conductors, simplified TL approach is not convenient for vertical and interconnected conductors. Also, the influence of earth–air interface is usually neglected [10]. Within the framework of the TL method, the effect of mutual coupling between different conductors of the grounding system is neglected. In general, TLM based solutions are limited to a certain upper frequency, depending on the electrical properties of the ground and configuration of particular grounding system. On other hand, the rigorous electromagnetic models based on antenna theory are the most accurate. The AM approach is based on solution of the Pocklington’s integro-differential equation for the half space problems. The earth–air interface effect is usually taken into account by the Sommerfeld integrals [11], or certain approximate approaches like modified image theory (MIT) [10,12], or reflection coefficient approximation (RC) [6,7]. The main drawback of the
n
Corresponding author. E-mail addresses: [email protected] (D. Cavka), [email protected] (D. Poljak), [email protected] (B. Nekhoul). 0955-7997/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.enganabound.2011.05.002
Sommerfeld integral formulation is rather long computational time, especially for the evaluation of broadband frequency spectrum. On the other hand, MIT approach accounts only for the electrical properties of the soil, but not the burial depth, while RC approximation produces error within 10% comparing to the rigorous Sommerfeld integral approach [11]. In this work, an assessment of the transient behavior of different grounding grid configurations using both the antenna and modified transmission line approaches has been carried out. This paper can be considered as a sequel of already published papers dealing with comparison of these two methods. In [13] the comparison of direct time domain approach based on antenna theory and TL model is presented for buried cables, while trade-off between the direct TLM approach and indirect frequency domain approach based on antenna theory model for transient analysis of grounding electrodes has been discussed in [14]. Furthermore, comparison of wire antenna and MTL approach to the assessment frequency response of horizontal grounding electrodes has been presented in [15]. It is worth mentioning that research presented in [13–15] is limited to a rather simple geometry of a single grounding electrode, while the research presented in this paper is extended to a complex grounding grid configuration, where mutual coupling between wires has to be taken into account along with the wire junctions effects. The approach presented in this paper, which can be considered as an extension of the work published in [15], within AM is based on the integro-differential equation of Pocklington type, with ground-air interface effects being taken into account trough exact Sommerfeld integral formulation. Contrary to the usual approach featuring the Moment Method [9], in this work the current distribution along the grounding grid is obtained by solving Pocklington integro-differential equation in the frequency domain
1102
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
via the Galerkin–Bubnov indirect Boundary Element Method [16] featuring linear isoparametric elements. Finally, the corresponding transient response is obtained by means of IFFT algorithm. The modified transmission line (MTL) model is based on the corresponding telegrapher’s equations. The current distribution along the grid conductors is assessed by solving partial differential equation for the scalar potential via finite difference technique directly in time domain. This procedure is followed by an integration of the related telegrapher’s equations to obtain the related current distribution.
2. Formulation The physical problem of interest is illustrated in Fig. 1. The square grounding grid is buried below ground and is subjected to a lighting lightning stroke at a certain point. Several grounding grid configurations are chosen for the analysis with dimensions varying from 10 10 to 30 30 m2, with or without additional vertical electrodes. All grids consist of wire conductors with radius a¼5 mm and buried at d ¼1.5 m depth. Fig. 2 shows various grid configurations. Two values of soil conductivity are considered; s1 ¼0.001 S/m (dry soil) and s2 ¼0.01 S/m (wet soil). In both cases relative permittivity is er ¼9. In all cases the current injection point is placed at the center of the grid. The transient analysis of grounding grids is undertaken using the two different approaches, wire antenna theory and modified transmission line model.
equations is derived from Maxwell equations by expressing the electric field in terms of magnetic vector and electric scalar potential, utilizing Lorentz gauge and by satisfying certain continuity conditions for the tangential field components at the electrode surface [6–9]. The outline of the derivation of this integral equation is given in Appendix, for the sake of completeness. The set of Pocklington equations is given by 2R 3 0 ^ 2 0 0 ^0 Cn0 In ðs ÞUsm Usn U½k1 þ rr g0n ðsm ,sn Þds N R 6 7 W X6 0 n 07 ^ ^n 2 Eexc 6 þRU Cn0 Iðsn ÞUsm Usn U½k1 þ rrgin ðsm ,sn Þds 7 sm ðsÞ ¼ CU 4 R 5 ! n¼1 þ C 0 In ðs0 ÞUs^ m U G s ðsm ,s0n Þds0 n
m ¼ 1,2,. . .NW ;
C¼
1 ; j4poeefec
R¼
k20 k21 k20 þk21
where I(s0 ) is the induced current along the line, Eexc s ðsÞ the excitation function, g0(s,s0 ) denotes the lossy medium Green function, while gi(s,s0 ) arises from the image theory. These functions are given by g0 ðs,s0 Þ ¼
ejk1 R0 R0
g1 ðs,sn Þ ¼
ejk1 R1 R1
2.1.1. Set of Pocklington integro-differential equations for arbitrarily shaped wires The currents flowing along the grounding grid configuration are governed by the set of coupled Pocklington integro-differential equations for wires of arbitrary shape. This set of Pocklington
ð2Þ
and R0 and R1 are the distances from the source point and its image to the observation point, respectively. Furthermore, k0 and k1 are propagation constants of air and lossy ground, respectively k20 ¼ o2 m0 e0
ð3Þ
s k21 ¼ o2 m0 eefec ¼ o2 m0 e0 er j
o
2.1. Antenna theory approach
ð1Þ
ð4Þ
where er and s are relative permittivity and conductivity of the ground, respectively, and o is operating frequency. The third term ! G s ðs,s0 Þ contains the Sommerfeld integrals and is constructed from the vector components for horizontal and vertical dipoles [11] ! ! ! 0 H H ! G s ðs,s0 Þ ¼ ðx^ Us^ ÞU GH r U r þ Gf U f þGz U z ! ! 0 þ ðz^ Us^ ÞU GVr U r þGVz U z ð5Þ The vector components are given in the Appendix. In the full wave analysis of grounding systems excited by the current source the left hand side of Eq. (1) vanishes and corresponding Pocklington Eq. (1) simplifies reducing to the homogenous one. Consequently, the excitation is incorporated into formulation through the boundary condition I1 ¼ Ig
Fig. 1. Square grounding grid subjected to a lightning stroke.
ð6Þ
where Ig denotes current generator and I1 current in the injection node. At a junction consisting of two or more segments the continuity properties of the electric field must be satisfied [17], which is ensured by applying the Kirchhoff current law n X
Ik ¼ 0
ð7Þ
k¼1
and continuity equation @I1 @I2 @In ¼ ¼ ¼ 0 0 @s1 at junctuion @s2 at junctuion @s0n at junctuion
Fig. 2. Different grounding grid configurations.
ð8Þ
The condition (8) ensures the discontinuities in charge per unit length to be ruled out in passing from one conductor to another across the junction. At the conductor free ends, the total current vanishes, i.e. is forced to be zero.
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
2.1.2. The evaluation of the input impedance spectrum The input impedance is given by the ratio Zin ¼
Vg Ig
þ k21 Us^ m Us^ 0n
1
Direct calculation of above integral is very time consuming. On the other hand, by carefully choosing an integration path, the computational cost can be appreciably reduced. In the case of horizontal arrangement of wires, like above examples, the best path is vertical, i.e. over z axis. The frequency dependent input impedance is then obtained and multiplied with current spectrum to obtain the frequency response of grounding system. Finally, transient response is calculated by means of the IFFT. 2.1.3. Numerical solution The set of Pocklington integro-differential Eq. (11) is numerically handled by means of the Galerkin–Bubnov variant of Indirect Boundary Element Method (GB-IBEM). The boundary element solution technique used in this work is an extension of the method applied to single wire cases and presented elsewhere, e.g. in [16]. The unknown current Ine ðzÞ along the nth wire segment is expressed by the sum of a finite number of linearly independent basis functions fni, with unknown complex coefficients Ini n X
Ine ðs0 Þ ¼
T
Ini fni ðs0 Þ ¼ ff gn fIgn
ð11Þ
i¼1
Z
Z
Z
Z
Dl m Dl n
fjm ðsm Þfin ðs0n Þg0nm ðsm ,s0n Þds0n dsm
dfjm ðsm Þ dfin ðs0n Þ ginm ðsm ,snn Þds0n dsm dsnn dsm Dlm Dln Z Z n þ RUk21 Ub s m Ub sn fjm ðsm Þfin ðs0n Þginm ðsm ,snn Þds0n dsm
ð9Þ
where Vg and Ig are the values of the voltage and the current at the driving point. After calculating the current distribution, a feeding point voltage is obtained by integrating the normal electric field component from infinity to the electrode surface Z r ! ! Vg ¼ Edl ð10Þ
Z
Z
1103
RU
þ s^ m
Dl m Dl n
Dlm Dln
! fjm ðsm Þfin ðs0n Þ G snm ðsm ,s0n Þds0n dsm
ð15Þ
Implementing isoparametric elements yields following expression for mutual impedance matrix: Z 1Z 1 dsn 0 dsm T fDgj fD0 gi g0nm ðsm ,s0n Þ 0 dz dz ½Zeji ¼ dz dz 1 1 Z 1Z 1 dsn 0 dsm T þ k21 Us^ m Us^ 0n ff gj ff 0 gi g0nm ðsm ,s0n Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm T RU fDgj fD0 gi ginm ðsm ,snn Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm n T þ RUk21 Us^ m Us^ n ff gj ff 0 gi ginm ðsm ,snn Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm T! þ s^ m U ff gj ff 0 gi G snm ðsm ,s0n Þ 0 dz dz ð16Þ dz dz 1 1 Matrices {f} and {f0 } contain the shape functions while {D} and {D0 } contain their directional derivatives. The excitation function in the form of the current source Ig is taken into account as a forced boundary condition at the certain node i of the grounding system Ii ¼ Ig The current related (7) and
ð17Þ treatment of wire junctions is related to Kirchhoff’s law in its integral and differential form, respectively, to continuity of currents and charges at the junction (8).
and the use of isoparametric elements yields Ine ðzÞ ¼
n X
T
Ini fni ðzÞ ¼ ff gn fIgn
2.2. Modified transmission line method approach ð12Þ
i¼1
where n is the number of local nodes on the element. A linear approximation over a boundary element along nth wire is used in this work and the corresponding shape functions are given by f1 ¼
1z 2
f2 ¼
1þz 2
ð13Þ
as this choice was proved to be optimal in modeling various wire structures [16]. Applying the weighted residual approach featuring the Galekin–Bubnov procedure the set of Pocklington equations is transformed into a system of algebraic equations. Performing a certain mathematical manipulations, the following matrix equation is obtained: Nw X Nn X
½Zeji fIn gi ¼ 0,
m ¼ 1,2,. . .,Nw ; j ¼ 1,2,. . .,Nm
ð14Þ
n¼1i¼1
where Nw is the total number of wires, Nm is number of elements on the mth antenna and Nn is number of elements on the nth antenna. ½Zji is the mutual impedance matrix for the jth observation boundary element on the mth antenna and ith source boundary element on the nth antenna, defined as follows: Z Z dfjm ðsm Þ dfin ðs0n Þ ½Zeji ¼ g0nm ðsm ,s0n Þds0n dsm ds0n dsm Dl m Dl n
Neglecting the transverse propagation effects the grounding system is simulated by means of a complex network. Using transmission lines equations a propagation equation for scalar potential for the assessment of both current and potential at an arbitrary point on the network (grounding grid) is deduced. The corresponding transmission line equations for the scalar potential and the current in the frequency domain for onedimensional (1D) case of propagation are given by 8 @I < @U @Z þ RI þ L @t ¼ 0 Z ¼ x or y ð18Þ @I : @Z þGU þ C @U @t ¼ 0 Combining the two telegrapher’s Eq. (18) the current or voltage can be eliminated and the second order partial differential equation for either potential or current is obtained. If the propagation occurs in two-directions; x and y, the two-dimensional (2D) differential equation and it is given by @2 U @2 U @U @2 U 2LC 2 ¼ 0 þ 2 2RGU2ðRC þ LGÞ 2 @t @x @y @t
ð19Þ
where R, L, C and G are per unit length parameters of the buried interconnected conductors. For the grounding grid, the per unit length parameters are calculated taking into account the soil–air interface effects. There are various approaches for the assessment of these parameters, i.e. using the formulas suggested by Sunde [18] or by Liu et al. [5].
1104
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
2.2.1. Finite difference solution of 2D propagation equation for scalar potential The spatial discretization of 2nd order differential operator at certain point (i, j) using the finite difference approximation is shown in Fig. 3. The finite difference approximation of spatial and temporal derivatives at certain point (i, j) are given by @2 U 1 n n ¼ ðU n 2Ui,j þUi1,j Þ @x2 Dx2 i þ 1,j
ð20aÞ
@2 U 1 n n n ¼ Ui,j þ 1 2Ui,j þ Ui,j1 2 2 @y Dy
ð20bÞ
@U 1 ¼ ðU n U n1 Þ @t Dt i,j i,j
ð20cÞ
@2 U 1 n1 n2 ¼ ðU n 2Ui,j þ Ui,j Þ @t 2 Dt2 i,j
ð20dÞ
"
Substituting (20a–d) into (19), the following expression arises: # 2 2 2ðRC þLGÞ 2LC n 2RG Ui,j Dt ðDxÞ2 ðDyÞ2 ðDtÞ2 " # " # " # " # 1 1 1 1 n n n n þ þ þ þ U U U Ui,j1 i þ 1,j i1,j i,j þ 1 ðDxÞ2 ðDxÞ2 ðDyÞ2 ðDyÞ2 ¼
! 2ðRC þ LGÞ 4LC 2LC n2 n1 þ Ui,j Ui,j Dt ðDtÞ2 ðDtÞ2
ð21Þ
The elements of vector [B] are given by ! 2ðRC þ LGÞ 4LC 2LC n2 n1 Bk ¼ Ui,j þ Ui,j Dt ðDtÞ2 ðDtÞ2
2.2.2. Boundary conditions and node injection The solution of the propagation Eq. (19) requires the knowledge of the conditions at the injection point and grid edges, respectively, shown in Fig. 4a. In the case of direct injection of the lightning strike, at an arbitrary point on the border of the grid (Fig. 4) Kirchhoff’s current low has to be imposed knowing that in the current capacity is very low 1 UðZ,tÞ ¼ ð26Þ IðtDtÞ Z ¼ x or y G0 where G0 is the equivalent conductance of corresponding nodes, IðtDtÞ is the transversal current known at (t Dt). The presence of two-media configuration (air and soil) is taken in account through the calculation of the linear parameters of the electrical circuit for the case of the grounding electrode [5,18]. Such a treatment of the non-homogeneous media is identical to the case of transmission line with ground return. At the injection node the value of current is known (lightning strike generator), which allows one to deduce the voltage at the node and to start calculations. The input impedance is defined by relation Zin ¼
which 2 A11 6 ^ 6 6 6 A1k 6 6 ^ 4 AN1
can also be expressed in matrix form 32 3 2 3 A1k A1N U1 B1 6 7 6 7 & ^ & ^ 7 76 ^ 7 6 ^ 7 76 7 6 7 6 Uk 7 ¼ 6 Bk 7 Akk AkN 7 76 7 6 7 7 6 7 6 & ^ & ^ 54 ^ 5 4 ^ 7 5 ANk ANN UN BN
2 ðDxÞ2
2 ðDyÞ2
2RG
2ðRC þ LGÞ 2LC Dt ðDtÞ2
ð22Þ
ð23Þ
while the non-diagonal elements of matrix [A] are, as follows: Akl ¼
Akl ¼
1 ðDxÞ2 1 ðDyÞ2
if 1 is the adjacent node k in x direction
ð24aÞ
if 1 is the adjacent node k in y direction
ð24bÞ
Akl ¼ 0 elsewhere:
Fig. 3. Spatial discretization of the square grid.
Uðk,tÞ Iðk,tÞ
ð27Þ
where k denotes the number of injection node.
where [A] denotes the coefficients matrix, [U] the unknown voltage vector, [B] the entire right-hand side and N the total number of nodes. The diagonal elements of matrix [A] are given by Akk ¼
ð25Þ
ð24cÞ
2.2.3. Currents distribution along the grid At every calculation step, once the transient voltages are computed, the currents induced in interconnected conductors of grounding grid are obtained by numerically integrating the following telegrapher’s equation: @U @I þ RI þ L ¼ 0 @Z @t
Z ¼ x or y
ð28Þ
3. Numerical results In all computational examples the lightning current is expressed by the double exponential function with parameters: I0 ¼1.1043, a ¼0.07924 106 s 1 and b ¼0.07924 106 s 1. Fig. 5 shows the transient voltage at the feeding point calculated by AM and TLM approach, respectively, for all four scenarios for grid configuration and soil conductivity s1 ¼0.001 S/m, while Fig. 6 shows resulting transient impedance of the grounding systems. Comparing the results obtained by two different approaches a very good agreement can be observed for grid types 1 and 2, good
Fig. 4. (a) Conditions at the grid edges; (b) equivalent electrical network of grounding grid.
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
1105
Such a behavior is due to the fact that parts of the grid behave as single antennas and there are many reflections from discontinuities which MTL fails to take into account and consequently ensure accurate results. This effect is more evident in lower conductive soils then in higher ones, what can be seen in Figs. 7 and 8, that show transient voltages and impedances at injection point calculated for ground conductivity of s2 ¼0.01 S/m. The agreement between results obtained via two methods is found to be satisfactory especially for later time instants which correspond to lower frequency part of the spectrum. For very early times the results are very similar, although the values of
Fig. 5. Transient feeding-point voltage for dry soil: (a) linear scale and (b) log scale.
Fig. 7. Transient fed-point voltage for wet soil: (a) linear scale and (b) log scale.
Fig. 6. Transient impedance for dry soil: (a) linear scale and (b) log scale.
agreement for type 3, while major differences occur for type 4 grid. Some differences appear in the very early (10 8–10 7 s) time instants which corresponds to the high frequency content of the input signal spectrum, which cannot be accurately predicted by the MTL method. Consequently, there are some differences in rather early time instants, as it is visible from Fig. 5, i.e. MTL fails to accurately predict the early time behavior of grounding grids. After that very early stage, the same can be drawn for the transient impedance, as well. Greater the grid size, the worse agreement between results is achieved. One would expect that MTL method would work better if the wires are longer, but it is not a case in this particular grid configuration.
Fig. 8. Transient impedance for wet soil: (a) linear scale and (b) log scale.
1106
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
peak transient voltage, for all cases of grid configuration are somewhat higher. Comparing the results for both values of ground conductivities it can be noticed that peak value of voltage is advanced for the case of higher conductivity. Also, the values of the peak voltages are pretty much alike regardless of the grid size. It is well-known that for very early time instants the higher frequency part of the impedance spectrum is important [15]. As the grid density (mesh) remains the same for all grid configurations that part of the frequency spectrum is unchanged regardless of the grid size, as shown in Fig. 9. Depending on the conductivity that part of spectrum starts at different frequencies. In the case of s2 ¼0.01 S/m that frequency is about 3 MHz (Fig. 9), while for s1 ¼0.001 S/m is above 30 MHz (Fig. 10). This is because in high conductivity environment effective length of the grounding wires becomes very short at higher frequencies. Due to different frequency spectrum at lower conductivity up to 30 MHz peak values of voltage are therefore very different and decrease with the grid size. Very early time behavior is identical for all configurations since the content of the spectrum related to higher frequencies is pretty much the same. Note that Fig. 9 is in logarithmic scale while Fig. 10 in linear, just for sake of clarity. Also, it should be noted that this analysis is based only on the results arising from the wire antenna approach, since MTL results are obtained directly in the time domain.
4. Conclusions Transient behavior of various grounding grid configurations has been analyzed using both the wire antenna theory and the modified transmission line method. The formulation based on wire antenna theory is related to the set of corresponding Pocklington integro-differential equation for curved wires, while the modified transmission line (MTL) approach is based on the time domain telegrapher’s equation. The integro-differential relationships arising from the antenna theory are numerically treated via the Galerkin–Bubnov scheme of the indirect Boundary Element Method. The MTL equations are solved using the finite difference method. The analysis is undertaken for different ground conductivities. Generally, the MTL method fails to predict accurate results for the very early time instants of the transient impedance, especially for the lower conductivity of the soil. At later time instants, there is a good agreement between the methods, although differences are higher as the grid size increases for the case of low conductivity soil. For the higher conductivity scenario an agreement between results obtained by different methods is quite satisfactory.
Appendix A. Derivation of the Pocklington equations set The set of Pocklington equations for a configuration of interconnected buried wires can be obtained as an extension of the Pocklington integro-differential equation for a single buried wire of arbitrary shape [11]. An extension to the case of multiple wires is straightforward [19]. The Pocklington equation for a single grounding electrode can be derived by enforcing the continuity conditions for the tangential components of the electric field along the wire surface. For the sake of simplicity, wire is first placed in infinite lossy medium, and then formulation is extended to a half space problem. For the PEC wire the total field composed from the excitation !exc !sct field E and scattered field E vanishes !exc !sct þ E Þ¼0 s^ Uð E
Fig. 9. Frequency spectrum of impedance for s2 ¼ 0.01.
on the wire surface
ðA1Þ
where s^ is the unit vector tangent at the observation point. Starting from Maxwell’s equations and Lorentz gauge the scattered electric field can be expressed in terms of the vector ! potential A ! !sct E ¼ jo A þ
! 1 r ðr A Þ jomeefec
ðA2Þ
The vector potential is defined by the particular integral Z 0 Iðs0 Þg0 ðs,s0 Þs^ ds0 ðA3Þ
! m A ðsÞ ¼ 4p
C
0 where I(s0 ) is the induced current along the line, s^ is the unit 0 vector tangent at the source point and g0(s,s ) is the corresponding Green’s function of the form
g0 ðs,s0 Þ ¼
ejk1 R0 R0
ðA4Þ
where R0 is the distances from the source point to the observation point qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðA5Þ R0 ¼ ðxx0 Þ2 þ ðyy0 Þ2 þ ðzz0 Þ2 þ a2
Fig. 10. Frequency spectrum of impedance for s1 ¼0.001.
Combining Eqs. (A1)–(A3) follows the Pocklington integro-differential equation for the unknown current distribution along the
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
single arbitrary wire antenna insulated in unbounded medium Z 0 Eexc ðsÞ ¼ CU Iðs0 ÞUs^ Us^ U½k21 þ rrg0 ðs,s0 Þds0 ðA6Þ C0
The integral Eq. (A6) can be extended for a case of a wire located near the interface between half-spaces by modifying the kernel to account for the field reflecting from the interface. Although various reflection coefficients approximations exist, the rigorous Sommerfeld integral approach is used in this work for the sake of accuracy. The excitation field component can be written as the sum of !inc the incident (direct) field E and field reflected from the inter!ref face E !exc !inc !ref E ¼ E þE
ðA7Þ
Now Eq. (A1) can be rewritten !inc !ref !sct s^ U E ¼ s^ U E s^ U E
ðA8Þ
!inc where first term s^ U E is given by (A6) while the field !ref is component due to the interface E 2 R 3 0 ^n 2 RU C 0 Iðs ÞUs U½k1 þ rrgi ðs,sn Þds0 !ref 5 ðA9Þ E ðsÞ ¼ CU4 R ! þ C 0 Iðs0 ÞU G s ðs,s0 Þds0 where gi ðs,sn Þ is the Green function arising from the image theory gi ðs,sn Þ ¼
where D1 ðlÞ ¼
D2 ðlÞ ¼
2
g0 þ g1
ðA10Þ
2k21
ðA14aÞ
g1 ðk21 þk20 Þ
ðA14bÞ
and
g0 ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffi l2 k20 ;
g1 ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffi l2 k21
ðA15Þ
Finally, combining Eqs. (A6)–(A9) yields the Pocklington integro-differential equation for the unknown current distribution along the single wire antenna of arbitrary shape buried in a lossy ground 2R 3 0 ^ ^0 Us U½k21 þ rrg0 ðs,s0 Þds0 0 Iðs ÞUs 6 C R 7 6 þRU 0 Iðs0 ÞUs^ Us^ n U½k2 þ rrgi ðs,sn Þds0 7 ðA16Þ Eexc 1 7 C s ðsÞ ¼ CU6 4 R 5 ! þ C 0 Iðs0 ÞUs^ U G s ðs,s0 Þds0 To derive the set of coupled Pocklington integro-differential equations for NW wires of arbitrary shape (1) the influence of each antenna has to be summarized, i.e. one obtains " NW R X 0 ^ 2 0 0 exc ^0 Esm ðsÞ ¼ CU C 0 In ðs ÞUsm Usn U½k1 þ rr g0n ðsm ,sn Þds n
n¼1
þ RU
Z
Cn0
þ Cn0
with R1 is the distance from the image point to the observation point, respectively qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi R1 ¼ ðxx0 Þ2 þ ðyy0 Þ2 þðz þ z0 Þ2 ðA11Þ
2 2 k21 g0 þ k20 g1 g1 ðk21 þk20 Þ
Z
ejk1 R1 R1
1107
n Iðs0n ÞUb s m Ub s n U½k21 þ rrgin ðsm ,snn Þds0
! In ðs’ÞUs^ m U G s ðsm ,s0n Þds0
m ¼ 1,2,. . .NW :
#
ðA17Þ
References
n
and s^ is the unit vector tangent at the source point of the image wire. ! The kernel G s ðs,s0 Þ is a correction term (5) containing the Sommerfeld integrals and is constructed from the following components for horizontal and vertical dipoles [11]: GVr ¼
@2 2 R k V @r@z 0
ðA12aÞ
! @2 2 2 R þk 1 k0 V @z2
GVz ¼
GH r ¼ cos j
ðA12bÞ
@2 2 R k V þ k21 U R @r2 1
!
1 @ 2 R k1 V þ k21 U R GH f ¼ sin f r @r V GH z ¼ j4poeefec cos f Gr
ðA12cÞ
ðA12dÞ ðA12eÞ
The superscript on the G denotes a vertical (V) or horizontal (H) current element and the subscript indicates the cylindrical component of the field vector. The horizontal current element is oriented along x axis. The Sommerfeld integral terms are: Z 1 0 UR ¼ D1 ðlÞeg1 9z þ z 9 J0 ðlrÞldl ðA13aÞ 0
VR ¼
Z
1 0
D2 ðlÞeg1 9z þ z 9 J0 ðlrÞldl 0
ðA13bÞ
[1] Meliopoulos AP, Moharam MG. Transient analysis of grounding systems. IEEE Transactions on Power Apparatus and Systems 1983;PAS-102(2):389–99. [2] Ramamoorty M, Narayanan MMB, Parameswaran S, Mukhedkar D. Transient performance of grounding grids. IEEE Transactions on Power Delivery 1989;4(4): 2053–9. [3] Liu Y, Zitnik M, Thottappillil R. Animproved transmission line model of grounding system. IEEE Transactions on EMC 2001;43(3):348–55. [4] Lorentzou MI, Hatziargyriou ND, Papadias C. Time domain analyisis of grounding electrodes impulse response. IEEE Transactions on Power Delivery 2003;2:517–24. [5] Liu Yaqing, Theethayi Nelson, Thottappillil Rajeev. An engineering model for transient analysis of grounding system under lightning strikes: nonuniform transmission-line approach. IEEE Transactions on Power Delivery 2005;20(2): 722–30. [6] Poljak D, Doric´ V. Wire antenna model for transient analysis of simple grounding systems. Part I: the Vertical grounding electrode. Progress in Electromagnetics Research 2006;64:149–66. [7] Poljak D, Doric´ V. Wire antenna model for transient analysis of simple grounding systems. Part II: the horizontal grounding electrode. Progress in Electromagnetics Research 2006;64:167–89. [8] Grcev L, Dawalibi F. An electromagnetic model for transients in grounding systems. IEEE Transactions on Power Delivery 1990;5(4):1773–81. [9] Grcev L, Heimbach M. Frequency dependent and transient characteristics of substation grounding systems. IEEE Transactions on Power Delivery 1997;12(1): 172–8. [10] Bridges GE. Transient plane wave coupling to bare and insulated cables buried in a lossy half-space. IEEE Transactons on EMC 1995;37(1):62–70. [11] Burke GJ, Miller EK. Modeling antennas near to and penetrating a lossy interface. IEEE Transactions on Antenna and Propagation 1984;AP32(10):1040–9. [12] Poljak D, Gizdic I, Roje V. Plane wave coupling to finite length cables buried in a lossy ground. Engineering Analysis with Boundary Elements 2002;26(1): 803–6. [13] Poljak D, Lucic R, Roje V. Transient analysis of electromagnetic field coupling to buried cables. In: Proceedings of ICEAA 01 international conference on electromagnetics in advanced applications. Torino, Italy; 2001. p. 335–8. [14] Poljak D, Doric V, Lucic R, Antonijevic S, Murko V. Comparison of wire antenna and transmission line model for transient analysis of grounding electrodes. Kanazawa: ICLP; 2006.
1108
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
[15] Poljak D, Doric V, Drissi KE, Kerroum K, Medic I. Comparison of wire antenna and modified transmission line approach to the assessment of frequency response of horizontal grounding electrodes. Engineering Analysis with Boundary Elements 2008;32:676–81. [16] Poljak D, Brebbia CA. Boundary element methods for electrical engineers. Southampton-Boston: WIT Press; 2005.
[17] King RWP, Wu TT. Analysis of crossed wires in a plane-wave field. IEEE Transactions on Electromagnetic Compatibility 1975;17(4):255–65. [18] Sunde ED. Earth conducting effects in transmission systems. New York, N.Y: Dover publications, Inc; 1968. [19] Poljak D. Advanced modeling in computational electromagnetic compatibility. New York: John Wiley and Sons; 2007.
Contents lists available at ScienceDirect
Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound
Wire antenna versus modified transmission line approach to the transient analysis of grounding grid Damir Cavka a, Basma Harrat b, Dragan Poljak a,n, Bachir Nekhoul b, Kamal Kerroum c, Khalil El Khamlichi Drissi c a
Department of Electronics, University of Split, Croatia LAMEL Laboratory, University of Jijel, Algeria c LASMEA Laboratory, University of Blaise Pascal, France b
a r t i c l e i n f o
abstract
Article history: Received 1 October 2010 Accepted 13 April 2011 Available online 7 June 2011
The paper deals with transient analysis of grounding grids using two different approaches, wire antenna theory and modified transmission line model. The Pocklington integro-differential equations, in frequency domain, arising from the wire antenna theory are numerically handled via the Galerkin– Bubnov variant of indirect Boundary Element Method (GB-IBEM), while the transient response was obtained using inverse Fourier transform. The modified transmission line equations are treated using the finite difference time domain (FDTD) method. Some illustrative numerical results are presented and discussed in the paper. & 2011 Elsevier Ltd. All rights reserved.
1. Introduction Transient analysis of grounding systems is of great importance in lightning protection systems (LPS) design. One of the most important parameters arising from the transient analysis is the transient impedance of the grounding electrode. Grounding systems can be modeled using simple electric circuit methods [1,2], the transmission line based model (TLM) [3–5] or antenna (full-wave) model (AM) [6–9]. While the circuit methods can be considered to be rather oversimplified, the TLM methods have advantage of simplicity and relatively low computational cost. On the other hand, though valid for long horizontal conductors, simplified TL approach is not convenient for vertical and interconnected conductors. Also, the influence of earth–air interface is usually neglected [10]. Within the framework of the TL method, the effect of mutual coupling between different conductors of the grounding system is neglected. In general, TLM based solutions are limited to a certain upper frequency, depending on the electrical properties of the ground and configuration of particular grounding system. On other hand, the rigorous electromagnetic models based on antenna theory are the most accurate. The AM approach is based on solution of the Pocklington’s integro-differential equation for the half space problems. The earth–air interface effect is usually taken into account by the Sommerfeld integrals [11], or certain approximate approaches like modified image theory (MIT) [10,12], or reflection coefficient approximation (RC) [6,7]. The main drawback of the
n
Corresponding author. E-mail addresses: [email protected] (D. Cavka), [email protected] (D. Poljak), [email protected] (B. Nekhoul). 0955-7997/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.enganabound.2011.05.002
Sommerfeld integral formulation is rather long computational time, especially for the evaluation of broadband frequency spectrum. On the other hand, MIT approach accounts only for the electrical properties of the soil, but not the burial depth, while RC approximation produces error within 10% comparing to the rigorous Sommerfeld integral approach [11]. In this work, an assessment of the transient behavior of different grounding grid configurations using both the antenna and modified transmission line approaches has been carried out. This paper can be considered as a sequel of already published papers dealing with comparison of these two methods. In [13] the comparison of direct time domain approach based on antenna theory and TL model is presented for buried cables, while trade-off between the direct TLM approach and indirect frequency domain approach based on antenna theory model for transient analysis of grounding electrodes has been discussed in [14]. Furthermore, comparison of wire antenna and MTL approach to the assessment frequency response of horizontal grounding electrodes has been presented in [15]. It is worth mentioning that research presented in [13–15] is limited to a rather simple geometry of a single grounding electrode, while the research presented in this paper is extended to a complex grounding grid configuration, where mutual coupling between wires has to be taken into account along with the wire junctions effects. The approach presented in this paper, which can be considered as an extension of the work published in [15], within AM is based on the integro-differential equation of Pocklington type, with ground-air interface effects being taken into account trough exact Sommerfeld integral formulation. Contrary to the usual approach featuring the Moment Method [9], in this work the current distribution along the grounding grid is obtained by solving Pocklington integro-differential equation in the frequency domain
1102
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
via the Galerkin–Bubnov indirect Boundary Element Method [16] featuring linear isoparametric elements. Finally, the corresponding transient response is obtained by means of IFFT algorithm. The modified transmission line (MTL) model is based on the corresponding telegrapher’s equations. The current distribution along the grid conductors is assessed by solving partial differential equation for the scalar potential via finite difference technique directly in time domain. This procedure is followed by an integration of the related telegrapher’s equations to obtain the related current distribution.
2. Formulation The physical problem of interest is illustrated in Fig. 1. The square grounding grid is buried below ground and is subjected to a lighting lightning stroke at a certain point. Several grounding grid configurations are chosen for the analysis with dimensions varying from 10 10 to 30 30 m2, with or without additional vertical electrodes. All grids consist of wire conductors with radius a¼5 mm and buried at d ¼1.5 m depth. Fig. 2 shows various grid configurations. Two values of soil conductivity are considered; s1 ¼0.001 S/m (dry soil) and s2 ¼0.01 S/m (wet soil). In both cases relative permittivity is er ¼9. In all cases the current injection point is placed at the center of the grid. The transient analysis of grounding grids is undertaken using the two different approaches, wire antenna theory and modified transmission line model.
equations is derived from Maxwell equations by expressing the electric field in terms of magnetic vector and electric scalar potential, utilizing Lorentz gauge and by satisfying certain continuity conditions for the tangential field components at the electrode surface [6–9]. The outline of the derivation of this integral equation is given in Appendix, for the sake of completeness. The set of Pocklington equations is given by 2R 3 0 ^ 2 0 0 ^0 Cn0 In ðs ÞUsm Usn U½k1 þ rr g0n ðsm ,sn Þds N R 6 7 W X6 0 n 07 ^ ^n 2 Eexc 6 þRU Cn0 Iðsn ÞUsm Usn U½k1 þ rrgin ðsm ,sn Þds 7 sm ðsÞ ¼ CU 4 R 5 ! n¼1 þ C 0 In ðs0 ÞUs^ m U G s ðsm ,s0n Þds0 n
m ¼ 1,2,. . .NW ;
C¼
1 ; j4poeefec
R¼
k20 k21 k20 þk21
where I(s0 ) is the induced current along the line, Eexc s ðsÞ the excitation function, g0(s,s0 ) denotes the lossy medium Green function, while gi(s,s0 ) arises from the image theory. These functions are given by g0 ðs,s0 Þ ¼
ejk1 R0 R0
g1 ðs,sn Þ ¼
ejk1 R1 R1
2.1.1. Set of Pocklington integro-differential equations for arbitrarily shaped wires The currents flowing along the grounding grid configuration are governed by the set of coupled Pocklington integro-differential equations for wires of arbitrary shape. This set of Pocklington
ð2Þ
and R0 and R1 are the distances from the source point and its image to the observation point, respectively. Furthermore, k0 and k1 are propagation constants of air and lossy ground, respectively k20 ¼ o2 m0 e0
ð3Þ
s k21 ¼ o2 m0 eefec ¼ o2 m0 e0 er j
o
2.1. Antenna theory approach
ð1Þ
ð4Þ
where er and s are relative permittivity and conductivity of the ground, respectively, and o is operating frequency. The third term ! G s ðs,s0 Þ contains the Sommerfeld integrals and is constructed from the vector components for horizontal and vertical dipoles [11] ! ! ! 0 H H ! G s ðs,s0 Þ ¼ ðx^ Us^ ÞU GH r U r þ Gf U f þGz U z ! ! 0 þ ðz^ Us^ ÞU GVr U r þGVz U z ð5Þ The vector components are given in the Appendix. In the full wave analysis of grounding systems excited by the current source the left hand side of Eq. (1) vanishes and corresponding Pocklington Eq. (1) simplifies reducing to the homogenous one. Consequently, the excitation is incorporated into formulation through the boundary condition I1 ¼ Ig
Fig. 1. Square grounding grid subjected to a lightning stroke.
ð6Þ
where Ig denotes current generator and I1 current in the injection node. At a junction consisting of two or more segments the continuity properties of the electric field must be satisfied [17], which is ensured by applying the Kirchhoff current law n X
Ik ¼ 0
ð7Þ
k¼1
and continuity equation @I1 @I2 @In ¼ ¼ ¼ 0 0 @s1 at junctuion @s2 at junctuion @s0n at junctuion
Fig. 2. Different grounding grid configurations.
ð8Þ
The condition (8) ensures the discontinuities in charge per unit length to be ruled out in passing from one conductor to another across the junction. At the conductor free ends, the total current vanishes, i.e. is forced to be zero.
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
2.1.2. The evaluation of the input impedance spectrum The input impedance is given by the ratio Zin ¼
Vg Ig
þ k21 Us^ m Us^ 0n
1
Direct calculation of above integral is very time consuming. On the other hand, by carefully choosing an integration path, the computational cost can be appreciably reduced. In the case of horizontal arrangement of wires, like above examples, the best path is vertical, i.e. over z axis. The frequency dependent input impedance is then obtained and multiplied with current spectrum to obtain the frequency response of grounding system. Finally, transient response is calculated by means of the IFFT. 2.1.3. Numerical solution The set of Pocklington integro-differential Eq. (11) is numerically handled by means of the Galerkin–Bubnov variant of Indirect Boundary Element Method (GB-IBEM). The boundary element solution technique used in this work is an extension of the method applied to single wire cases and presented elsewhere, e.g. in [16]. The unknown current Ine ðzÞ along the nth wire segment is expressed by the sum of a finite number of linearly independent basis functions fni, with unknown complex coefficients Ini n X
Ine ðs0 Þ ¼
T
Ini fni ðs0 Þ ¼ ff gn fIgn
ð11Þ
i¼1
Z
Z
Z
Z
Dl m Dl n
fjm ðsm Þfin ðs0n Þg0nm ðsm ,s0n Þds0n dsm
dfjm ðsm Þ dfin ðs0n Þ ginm ðsm ,snn Þds0n dsm dsnn dsm Dlm Dln Z Z n þ RUk21 Ub s m Ub sn fjm ðsm Þfin ðs0n Þginm ðsm ,snn Þds0n dsm
ð9Þ
where Vg and Ig are the values of the voltage and the current at the driving point. After calculating the current distribution, a feeding point voltage is obtained by integrating the normal electric field component from infinity to the electrode surface Z r ! ! Vg ¼ Edl ð10Þ
Z
Z
1103
RU
þ s^ m
Dl m Dl n
Dlm Dln
! fjm ðsm Þfin ðs0n Þ G snm ðsm ,s0n Þds0n dsm
ð15Þ
Implementing isoparametric elements yields following expression for mutual impedance matrix: Z 1Z 1 dsn 0 dsm T fDgj fD0 gi g0nm ðsm ,s0n Þ 0 dz dz ½Zeji ¼ dz dz 1 1 Z 1Z 1 dsn 0 dsm T þ k21 Us^ m Us^ 0n ff gj ff 0 gi g0nm ðsm ,s0n Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm T RU fDgj fD0 gi ginm ðsm ,snn Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm n T þ RUk21 Us^ m Us^ n ff gj ff 0 gi ginm ðsm ,snn Þ 0 dz dz dz dz 1 1 Z 1Z 1 dsn 0 dsm T! þ s^ m U ff gj ff 0 gi G snm ðsm ,s0n Þ 0 dz dz ð16Þ dz dz 1 1 Matrices {f} and {f0 } contain the shape functions while {D} and {D0 } contain their directional derivatives. The excitation function in the form of the current source Ig is taken into account as a forced boundary condition at the certain node i of the grounding system Ii ¼ Ig The current related (7) and
ð17Þ treatment of wire junctions is related to Kirchhoff’s law in its integral and differential form, respectively, to continuity of currents and charges at the junction (8).
and the use of isoparametric elements yields Ine ðzÞ ¼
n X
T
Ini fni ðzÞ ¼ ff gn fIgn
2.2. Modified transmission line method approach ð12Þ
i¼1
where n is the number of local nodes on the element. A linear approximation over a boundary element along nth wire is used in this work and the corresponding shape functions are given by f1 ¼
1z 2
f2 ¼
1þz 2
ð13Þ
as this choice was proved to be optimal in modeling various wire structures [16]. Applying the weighted residual approach featuring the Galekin–Bubnov procedure the set of Pocklington equations is transformed into a system of algebraic equations. Performing a certain mathematical manipulations, the following matrix equation is obtained: Nw X Nn X
½Zeji fIn gi ¼ 0,
m ¼ 1,2,. . .,Nw ; j ¼ 1,2,. . .,Nm
ð14Þ
n¼1i¼1
where Nw is the total number of wires, Nm is number of elements on the mth antenna and Nn is number of elements on the nth antenna. ½Zji is the mutual impedance matrix for the jth observation boundary element on the mth antenna and ith source boundary element on the nth antenna, defined as follows: Z Z dfjm ðsm Þ dfin ðs0n Þ ½Zeji ¼ g0nm ðsm ,s0n Þds0n dsm ds0n dsm Dl m Dl n
Neglecting the transverse propagation effects the grounding system is simulated by means of a complex network. Using transmission lines equations a propagation equation for scalar potential for the assessment of both current and potential at an arbitrary point on the network (grounding grid) is deduced. The corresponding transmission line equations for the scalar potential and the current in the frequency domain for onedimensional (1D) case of propagation are given by 8 @I < @U @Z þ RI þ L @t ¼ 0 Z ¼ x or y ð18Þ @I : @Z þGU þ C @U @t ¼ 0 Combining the two telegrapher’s Eq. (18) the current or voltage can be eliminated and the second order partial differential equation for either potential or current is obtained. If the propagation occurs in two-directions; x and y, the two-dimensional (2D) differential equation and it is given by @2 U @2 U @U @2 U 2LC 2 ¼ 0 þ 2 2RGU2ðRC þ LGÞ 2 @t @x @y @t
ð19Þ
where R, L, C and G are per unit length parameters of the buried interconnected conductors. For the grounding grid, the per unit length parameters are calculated taking into account the soil–air interface effects. There are various approaches for the assessment of these parameters, i.e. using the formulas suggested by Sunde [18] or by Liu et al. [5].
1104
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
2.2.1. Finite difference solution of 2D propagation equation for scalar potential The spatial discretization of 2nd order differential operator at certain point (i, j) using the finite difference approximation is shown in Fig. 3. The finite difference approximation of spatial and temporal derivatives at certain point (i, j) are given by @2 U 1 n n ¼ ðU n 2Ui,j þUi1,j Þ @x2 Dx2 i þ 1,j
ð20aÞ
@2 U 1 n n n ¼ Ui,j þ 1 2Ui,j þ Ui,j1 2 2 @y Dy
ð20bÞ
@U 1 ¼ ðU n U n1 Þ @t Dt i,j i,j
ð20cÞ
@2 U 1 n1 n2 ¼ ðU n 2Ui,j þ Ui,j Þ @t 2 Dt2 i,j
ð20dÞ
"
Substituting (20a–d) into (19), the following expression arises: # 2 2 2ðRC þLGÞ 2LC n 2RG Ui,j Dt ðDxÞ2 ðDyÞ2 ðDtÞ2 " # " # " # " # 1 1 1 1 n n n n þ þ þ þ U U U Ui,j1 i þ 1,j i1,j i,j þ 1 ðDxÞ2 ðDxÞ2 ðDyÞ2 ðDyÞ2 ¼
! 2ðRC þ LGÞ 4LC 2LC n2 n1 þ Ui,j Ui,j Dt ðDtÞ2 ðDtÞ2
ð21Þ
The elements of vector [B] are given by ! 2ðRC þ LGÞ 4LC 2LC n2 n1 Bk ¼ Ui,j þ Ui,j Dt ðDtÞ2 ðDtÞ2
2.2.2. Boundary conditions and node injection The solution of the propagation Eq. (19) requires the knowledge of the conditions at the injection point and grid edges, respectively, shown in Fig. 4a. In the case of direct injection of the lightning strike, at an arbitrary point on the border of the grid (Fig. 4) Kirchhoff’s current low has to be imposed knowing that in the current capacity is very low 1 UðZ,tÞ ¼ ð26Þ IðtDtÞ Z ¼ x or y G0 where G0 is the equivalent conductance of corresponding nodes, IðtDtÞ is the transversal current known at (t Dt). The presence of two-media configuration (air and soil) is taken in account through the calculation of the linear parameters of the electrical circuit for the case of the grounding electrode [5,18]. Such a treatment of the non-homogeneous media is identical to the case of transmission line with ground return. At the injection node the value of current is known (lightning strike generator), which allows one to deduce the voltage at the node and to start calculations. The input impedance is defined by relation Zin ¼
which 2 A11 6 ^ 6 6 6 A1k 6 6 ^ 4 AN1
can also be expressed in matrix form 32 3 2 3 A1k A1N U1 B1 6 7 6 7 & ^ & ^ 7 76 ^ 7 6 ^ 7 76 7 6 7 6 Uk 7 ¼ 6 Bk 7 Akk AkN 7 76 7 6 7 7 6 7 6 & ^ & ^ 54 ^ 5 4 ^ 7 5 ANk ANN UN BN
2 ðDxÞ2
2 ðDyÞ2
2RG
2ðRC þ LGÞ 2LC Dt ðDtÞ2
ð22Þ
ð23Þ
while the non-diagonal elements of matrix [A] are, as follows: Akl ¼
Akl ¼
1 ðDxÞ2 1 ðDyÞ2
if 1 is the adjacent node k in x direction
ð24aÞ
if 1 is the adjacent node k in y direction
ð24bÞ
Akl ¼ 0 elsewhere:
Fig. 3. Spatial discretization of the square grid.
Uðk,tÞ Iðk,tÞ
ð27Þ
where k denotes the number of injection node.
where [A] denotes the coefficients matrix, [U] the unknown voltage vector, [B] the entire right-hand side and N the total number of nodes. The diagonal elements of matrix [A] are given by Akk ¼
ð25Þ
ð24cÞ
2.2.3. Currents distribution along the grid At every calculation step, once the transient voltages are computed, the currents induced in interconnected conductors of grounding grid are obtained by numerically integrating the following telegrapher’s equation: @U @I þ RI þ L ¼ 0 @Z @t
Z ¼ x or y
ð28Þ
3. Numerical results In all computational examples the lightning current is expressed by the double exponential function with parameters: I0 ¼1.1043, a ¼0.07924 106 s 1 and b ¼0.07924 106 s 1. Fig. 5 shows the transient voltage at the feeding point calculated by AM and TLM approach, respectively, for all four scenarios for grid configuration and soil conductivity s1 ¼0.001 S/m, while Fig. 6 shows resulting transient impedance of the grounding systems. Comparing the results obtained by two different approaches a very good agreement can be observed for grid types 1 and 2, good
Fig. 4. (a) Conditions at the grid edges; (b) equivalent electrical network of grounding grid.
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
1105
Such a behavior is due to the fact that parts of the grid behave as single antennas and there are many reflections from discontinuities which MTL fails to take into account and consequently ensure accurate results. This effect is more evident in lower conductive soils then in higher ones, what can be seen in Figs. 7 and 8, that show transient voltages and impedances at injection point calculated for ground conductivity of s2 ¼0.01 S/m. The agreement between results obtained via two methods is found to be satisfactory especially for later time instants which correspond to lower frequency part of the spectrum. For very early times the results are very similar, although the values of
Fig. 5. Transient feeding-point voltage for dry soil: (a) linear scale and (b) log scale.
Fig. 7. Transient fed-point voltage for wet soil: (a) linear scale and (b) log scale.
Fig. 6. Transient impedance for dry soil: (a) linear scale and (b) log scale.
agreement for type 3, while major differences occur for type 4 grid. Some differences appear in the very early (10 8–10 7 s) time instants which corresponds to the high frequency content of the input signal spectrum, which cannot be accurately predicted by the MTL method. Consequently, there are some differences in rather early time instants, as it is visible from Fig. 5, i.e. MTL fails to accurately predict the early time behavior of grounding grids. After that very early stage, the same can be drawn for the transient impedance, as well. Greater the grid size, the worse agreement between results is achieved. One would expect that MTL method would work better if the wires are longer, but it is not a case in this particular grid configuration.
Fig. 8. Transient impedance for wet soil: (a) linear scale and (b) log scale.
1106
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
peak transient voltage, for all cases of grid configuration are somewhat higher. Comparing the results for both values of ground conductivities it can be noticed that peak value of voltage is advanced for the case of higher conductivity. Also, the values of the peak voltages are pretty much alike regardless of the grid size. It is well-known that for very early time instants the higher frequency part of the impedance spectrum is important [15]. As the grid density (mesh) remains the same for all grid configurations that part of the frequency spectrum is unchanged regardless of the grid size, as shown in Fig. 9. Depending on the conductivity that part of spectrum starts at different frequencies. In the case of s2 ¼0.01 S/m that frequency is about 3 MHz (Fig. 9), while for s1 ¼0.001 S/m is above 30 MHz (Fig. 10). This is because in high conductivity environment effective length of the grounding wires becomes very short at higher frequencies. Due to different frequency spectrum at lower conductivity up to 30 MHz peak values of voltage are therefore very different and decrease with the grid size. Very early time behavior is identical for all configurations since the content of the spectrum related to higher frequencies is pretty much the same. Note that Fig. 9 is in logarithmic scale while Fig. 10 in linear, just for sake of clarity. Also, it should be noted that this analysis is based only on the results arising from the wire antenna approach, since MTL results are obtained directly in the time domain.
4. Conclusions Transient behavior of various grounding grid configurations has been analyzed using both the wire antenna theory and the modified transmission line method. The formulation based on wire antenna theory is related to the set of corresponding Pocklington integro-differential equation for curved wires, while the modified transmission line (MTL) approach is based on the time domain telegrapher’s equation. The integro-differential relationships arising from the antenna theory are numerically treated via the Galerkin–Bubnov scheme of the indirect Boundary Element Method. The MTL equations are solved using the finite difference method. The analysis is undertaken for different ground conductivities. Generally, the MTL method fails to predict accurate results for the very early time instants of the transient impedance, especially for the lower conductivity of the soil. At later time instants, there is a good agreement between the methods, although differences are higher as the grid size increases for the case of low conductivity soil. For the higher conductivity scenario an agreement between results obtained by different methods is quite satisfactory.
Appendix A. Derivation of the Pocklington equations set The set of Pocklington equations for a configuration of interconnected buried wires can be obtained as an extension of the Pocklington integro-differential equation for a single buried wire of arbitrary shape [11]. An extension to the case of multiple wires is straightforward [19]. The Pocklington equation for a single grounding electrode can be derived by enforcing the continuity conditions for the tangential components of the electric field along the wire surface. For the sake of simplicity, wire is first placed in infinite lossy medium, and then formulation is extended to a half space problem. For the PEC wire the total field composed from the excitation !exc !sct field E and scattered field E vanishes !exc !sct þ E Þ¼0 s^ Uð E
Fig. 9. Frequency spectrum of impedance for s2 ¼ 0.01.
on the wire surface
ðA1Þ
where s^ is the unit vector tangent at the observation point. Starting from Maxwell’s equations and Lorentz gauge the scattered electric field can be expressed in terms of the vector ! potential A ! !sct E ¼ jo A þ
! 1 r ðr A Þ jomeefec
ðA2Þ
The vector potential is defined by the particular integral Z 0 Iðs0 Þg0 ðs,s0 Þs^ ds0 ðA3Þ
! m A ðsÞ ¼ 4p
C
0 where I(s0 ) is the induced current along the line, s^ is the unit 0 vector tangent at the source point and g0(s,s ) is the corresponding Green’s function of the form
g0 ðs,s0 Þ ¼
ejk1 R0 R0
ðA4Þ
where R0 is the distances from the source point to the observation point qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðA5Þ R0 ¼ ðxx0 Þ2 þ ðyy0 Þ2 þ ðzz0 Þ2 þ a2
Fig. 10. Frequency spectrum of impedance for s1 ¼0.001.
Combining Eqs. (A1)–(A3) follows the Pocklington integro-differential equation for the unknown current distribution along the
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
single arbitrary wire antenna insulated in unbounded medium Z 0 Eexc ðsÞ ¼ CU Iðs0 ÞUs^ Us^ U½k21 þ rrg0 ðs,s0 Þds0 ðA6Þ C0
The integral Eq. (A6) can be extended for a case of a wire located near the interface between half-spaces by modifying the kernel to account for the field reflecting from the interface. Although various reflection coefficients approximations exist, the rigorous Sommerfeld integral approach is used in this work for the sake of accuracy. The excitation field component can be written as the sum of !inc the incident (direct) field E and field reflected from the inter!ref face E !exc !inc !ref E ¼ E þE
ðA7Þ
Now Eq. (A1) can be rewritten !inc !ref !sct s^ U E ¼ s^ U E s^ U E
ðA8Þ
!inc where first term s^ U E is given by (A6) while the field !ref is component due to the interface E 2 R 3 0 ^n 2 RU C 0 Iðs ÞUs U½k1 þ rrgi ðs,sn Þds0 !ref 5 ðA9Þ E ðsÞ ¼ CU4 R ! þ C 0 Iðs0 ÞU G s ðs,s0 Þds0 where gi ðs,sn Þ is the Green function arising from the image theory gi ðs,sn Þ ¼
where D1 ðlÞ ¼
D2 ðlÞ ¼
2
g0 þ g1
ðA10Þ
2k21
ðA14aÞ
g1 ðk21 þk20 Þ
ðA14bÞ
and
g0 ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffi l2 k20 ;
g1 ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffi l2 k21
ðA15Þ
Finally, combining Eqs. (A6)–(A9) yields the Pocklington integro-differential equation for the unknown current distribution along the single wire antenna of arbitrary shape buried in a lossy ground 2R 3 0 ^ ^0 Us U½k21 þ rrg0 ðs,s0 Þds0 0 Iðs ÞUs 6 C R 7 6 þRU 0 Iðs0 ÞUs^ Us^ n U½k2 þ rrgi ðs,sn Þds0 7 ðA16Þ Eexc 1 7 C s ðsÞ ¼ CU6 4 R 5 ! þ C 0 Iðs0 ÞUs^ U G s ðs,s0 Þds0 To derive the set of coupled Pocklington integro-differential equations for NW wires of arbitrary shape (1) the influence of each antenna has to be summarized, i.e. one obtains " NW R X 0 ^ 2 0 0 exc ^0 Esm ðsÞ ¼ CU C 0 In ðs ÞUsm Usn U½k1 þ rr g0n ðsm ,sn Þds n
n¼1
þ RU
Z
Cn0
þ Cn0
with R1 is the distance from the image point to the observation point, respectively qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi R1 ¼ ðxx0 Þ2 þ ðyy0 Þ2 þðz þ z0 Þ2 ðA11Þ
2 2 k21 g0 þ k20 g1 g1 ðk21 þk20 Þ
Z
ejk1 R1 R1
1107
n Iðs0n ÞUb s m Ub s n U½k21 þ rrgin ðsm ,snn Þds0
! In ðs’ÞUs^ m U G s ðsm ,s0n Þds0
m ¼ 1,2,. . .NW :
#
ðA17Þ
References
n
and s^ is the unit vector tangent at the source point of the image wire. ! The kernel G s ðs,s0 Þ is a correction term (5) containing the Sommerfeld integrals and is constructed from the following components for horizontal and vertical dipoles [11]: GVr ¼
@2 2 R k V @r@z 0
ðA12aÞ
! @2 2 2 R þk 1 k0 V @z2
GVz ¼
GH r ¼ cos j
ðA12bÞ
@2 2 R k V þ k21 U R @r2 1
!
1 @ 2 R k1 V þ k21 U R GH f ¼ sin f r @r V GH z ¼ j4poeefec cos f Gr
ðA12cÞ
ðA12dÞ ðA12eÞ
The superscript on the G denotes a vertical (V) or horizontal (H) current element and the subscript indicates the cylindrical component of the field vector. The horizontal current element is oriented along x axis. The Sommerfeld integral terms are: Z 1 0 UR ¼ D1 ðlÞeg1 9z þ z 9 J0 ðlrÞldl ðA13aÞ 0
VR ¼
Z
1 0
D2 ðlÞeg1 9z þ z 9 J0 ðlrÞldl 0
ðA13bÞ
[1] Meliopoulos AP, Moharam MG. Transient analysis of grounding systems. IEEE Transactions on Power Apparatus and Systems 1983;PAS-102(2):389–99. [2] Ramamoorty M, Narayanan MMB, Parameswaran S, Mukhedkar D. Transient performance of grounding grids. IEEE Transactions on Power Delivery 1989;4(4): 2053–9. [3] Liu Y, Zitnik M, Thottappillil R. Animproved transmission line model of grounding system. IEEE Transactions on EMC 2001;43(3):348–55. [4] Lorentzou MI, Hatziargyriou ND, Papadias C. Time domain analyisis of grounding electrodes impulse response. IEEE Transactions on Power Delivery 2003;2:517–24. [5] Liu Yaqing, Theethayi Nelson, Thottappillil Rajeev. An engineering model for transient analysis of grounding system under lightning strikes: nonuniform transmission-line approach. IEEE Transactions on Power Delivery 2005;20(2): 722–30. [6] Poljak D, Doric´ V. Wire antenna model for transient analysis of simple grounding systems. Part I: the Vertical grounding electrode. Progress in Electromagnetics Research 2006;64:149–66. [7] Poljak D, Doric´ V. Wire antenna model for transient analysis of simple grounding systems. Part II: the horizontal grounding electrode. Progress in Electromagnetics Research 2006;64:167–89. [8] Grcev L, Dawalibi F. An electromagnetic model for transients in grounding systems. IEEE Transactions on Power Delivery 1990;5(4):1773–81. [9] Grcev L, Heimbach M. Frequency dependent and transient characteristics of substation grounding systems. IEEE Transactions on Power Delivery 1997;12(1): 172–8. [10] Bridges GE. Transient plane wave coupling to bare and insulated cables buried in a lossy half-space. IEEE Transactons on EMC 1995;37(1):62–70. [11] Burke GJ, Miller EK. Modeling antennas near to and penetrating a lossy interface. IEEE Transactions on Antenna and Propagation 1984;AP32(10):1040–9. [12] Poljak D, Gizdic I, Roje V. Plane wave coupling to finite length cables buried in a lossy ground. Engineering Analysis with Boundary Elements 2002;26(1): 803–6. [13] Poljak D, Lucic R, Roje V. Transient analysis of electromagnetic field coupling to buried cables. In: Proceedings of ICEAA 01 international conference on electromagnetics in advanced applications. Torino, Italy; 2001. p. 335–8. [14] Poljak D, Doric V, Lucic R, Antonijevic S, Murko V. Comparison of wire antenna and transmission line model for transient analysis of grounding electrodes. Kanazawa: ICLP; 2006.
1108
D. Cavka et al. / Engineering Analysis with Boundary Elements 35 (2011) 1101–1108
[15] Poljak D, Doric V, Drissi KE, Kerroum K, Medic I. Comparison of wire antenna and modified transmission line approach to the assessment of frequency response of horizontal grounding electrodes. Engineering Analysis with Boundary Elements 2008;32:676–81. [16] Poljak D, Brebbia CA. Boundary element methods for electrical engineers. Southampton-Boston: WIT Press; 2005.
[17] King RWP, Wu TT. Analysis of crossed wires in a plane-wave field. IEEE Transactions on Electromagnetic Compatibility 1975;17(4):255–65. [18] Sunde ED. Earth conducting effects in transmission systems. New York, N.Y: Dover publications, Inc; 1968. [19] Poljak D. Advanced modeling in computational electromagnetic compatibility. New York: John Wiley and Sons; 2007.