Transmission Line Modeling of Grounding

Electrodes and Calculation of their Effective Length

under Impulse Excitation
M. I. Lorentzou, N. D.Hatziargyriou, Senior Member, IEEE

Abstract-- Grounding electrodes, being the basic component literature [1],[18]. Grounding electrode is treated as in series
of any grounding system, need to be accurately modeled in connected π-circuits tending to the open-ended transmission
transient analysis studies. An easy-to-apply methodology for
calculation of voltage and current distribution along the line when the number of circuits increases. Analytical
electrode has been presented and validated in [1]. It is based on formulae are presented for current and voltage distributions
a distributed parameters transmission line model of the along the electrode. Impulse injection current has been
electrode. In the middle stages, closed–form mathematical
described by a double exponential waveform function.
expressions are used for the solution of telegraphy equations
that describe the propagation of voltage and current waves Grounding electrode model is constructed under the
along the electrode. The final expressions include few terms, so assumption that the soil is homogeneous. In this paper per
they are easy to apply. In this paper per unit length parameters unit length inductance conductance capacitance and
of the electrode are comparatively calculated using various
approaches. Calculation of the impulse response of grounding
resistance of the electrode are comparatively calculated using
electrodes shows that there is a limit in their length that various approaches and mathematical formulas.
seriously contributes in lowering raised potentials. This limit is Long grounding electrodes have lower grounding
the “effective length” of grounding electrodes and is calculated resistance and better behaviour under low and high frequency
in the paper.
excitation. In case a sinusoidal current source is applied at a
Keywords: Grounding Electrode, Transient Response,
point of a grounding conductor, maximum GPR (Ground
Effective Length, Ground Potential Rise (GPR). Potential Rise) decreases as its dimensions increase. However,
when high frequency components are injected as in the case
of lightning, there is an upper limit in the length of the
I. INTRODUCTION electrode that substantially affects the maximum GPR value.

M ODELING of grounding systems and grounding

electrodes is of primary importance for the analysis of
the transient behavior of any electrical installation. Many
This is called here ‘effective length”. In this paper, the effect
of effective length of grounding electrodes is demonstrated,
for various soil and impulse excitation characteristics.
attempts to accurately determine the transient response of a
grounding electrode are presented in literature [2]- II. THEORETICAL BACKGROUND
[5],[12],[15],[16]. These attempts are either based on A. Grounding Electrode Model
computer models, solved numerically or on analytical Correct calculation of per unit length parameters of the
expressions for current and voltage distributions based on electrode is essential for the accuracy of the results. Special
simplifications or special initial conditions [2],[3],[17]. attention is given to the assumptions done in the middle
In this paper a novel methodology applicable on the stages of the development of calculation formulas.
calculation of the behavior of a grounding electrode under
surge conditions is used. Attempts made in the past involve A.1. Circuit model of the electrode
techniques for Laplace inversion or prediction of the results As shown in literature at low frequencies the impedance of a
based on experimental data. Advantages of the new grounding electrode can be represented by a single resistance,
methodology are simplicity and accuracy, because it is based while at high frequencies by a lumped R-L-C circuit. Three
on a distributed model of the electrode and closed form sets of formulas for the parameters of the circuit model of the
mathematical expressions. This methodology is validated in electrode, resistance Re, inductance Le, conductance Ge and
capacitance Ce are often used in the existing literature.
M. I. Lorentzou is with the Hellenic Transmission System Operator, 17122 One is from the reference work by Rudenberg [6] and is
Athens, GREECE (e-mail: [email protected]). used for the vertical rod:
N. D. Hatziargyriou is professor at the Department of Electrical Engineering, ρ 2l 2l
Power Systems Division, of National Technical University of Athens, 15773 G −1 = log , C = 2πεl log and
Athens, GREECE (e-mail: [email protected]). 2πl α α
µ 0l 2l
Presented at the International Conference on Power Systems Transients L= log (1)
(IPST’05) in Montreal, Canada on June 19-23, 2005 2π α
Paper No. IPST05 - 239
 1

and the other [7][8] is used for horizontal electrodes: where γ c = [ jωµ c (σ c + jωε c )]2
ρ  4l   4l  and ω is the angular frequency. A simplified equation for
G −1 =  log − 1 , C = 2πε log − 1 , easier and accurate calculation of Zi can be found in [20].At
2πl  α   α 
very low and zero frequencies this impedance reduces to the
µ 0 l  2l  classical dc resistance of a cylindrical conductor:
L= log − 1 (2)
α  α  Zi =
1
(8)
where α the radius of the conductor is replaced by πα 2σ c
2 ⋅ a ⋅ z when the conductor is placed at depth z.
Alternatively the following equation proposed by Dwight et A.3.Comparative evaluation of calculation formulas
al. [9] can be used to calculate G: A drawback of most calculation formulas for estimation of
the lumped parameters Re, Le, Ce and Ge is that they are based
ρ   2l   2l  s s2 s4 
G−1 = ⋅ ln  + ln  − 2 + − + 4
on quasi-static approximation that limits their validity to
2πl   α   s  l 16⋅ (0.5l) 512⋅ (0.5l) 
2
some upper frequency that depends on the size of the
(3) grounding system and the electrical characteristics of the
where s=2z is two times the burial depth and L, C can be earth.
calculated from the equations found in [4]: More particularly, Sunde’s assumptions in derivation of
l the electrode parameters are that the potential due to current
L= , C = ε ⋅ ρ ⋅G (4) leaving a conductor element is the same for a point source at
c0 ⋅ ε ⋅ ρ ⋅ G
2
the axis of the conductor. Thus is acceptable when dealing
Other, less commonly used calculation formulas can be with long conductors. In order to solve the equations of the
found in literature [3][17] potential versus the leakage current along the entire
The electrode can then be considered as a transmission conductor it is assumed that the resistance of the conductor is
line open at the lower end. The input impedance at the negligible, so that the voltage drop may be disregarded and
energization end is [11]: the boundary condition to be satisfied at the surface of the
Z = Z 0 coth γl (5) conductor y=α is that dV ( x,α ) dx = 0
Considering the accuracy of (6), as noted in [7], when the
A.2. Frequency dependent transmission line model current varies rapidly, as for example in the case of lightning
For given conductor and soil data an electrodynamic strike, these are an improved approximation of the solution in
model can be evaluated that results in the frequency terms of an exponential mode of propagation. There is a
dependent characteristic impedance Zc(ω) and propagation logarithmic approximation of the variation vs. conductor
constant γ(ω)[13][14]. Sunde [7] developed the mathematical length of the main integral functions q, p which are added to
equations for the unit length admittance and self-impedance the internal impedance and admittance terms to give (6)
of the grounding electrode: Grounding resistance values calculated from formula (3)
−1 are closer to the results of finite elements method calculation
 −1 1 1.12 K 
Y (Γ ) =  Yi + log  and well known commercial software [21] than the results
 π (σ + jωε ) Γα  calculated from formula (2).
jωµ 1.18 K In Fig.1, the magnitude of the input impedance |Z| at the
Z (Γ ) = Z i + log (6)
2π (
α γ 2 + Γ2
12
) energization end obtained from the various calculation
formulas, vs. electrode length is plotted. Curves 1 and 2
where α is replaced by a + 4z when the conductor is
2 2
show the results of (2) and (3) when the parameters R, L, C,
placed at a depth z. The above equations (6) are solved G, are calculated for the whole length of the electrode and
Γ Z( Γ ) Y( Γ ) then they are divided to find the p.u. values. Curve 3 shows
iteratively for at every single frequency in the results from calculation formulas (6). Curves 4, 6 are the
the range of interest. In practical examined cases the internal input |Z| of n-pi-circuit model of the electrode where the
impedance Zi that represents energy within the conductor is parameters R, L, C, G, are calculated for the whole length of
nonzero. The internal impedance of cylindrical conductor is the electrode and then they are divided to find the p.u. values
derived in a classical paper by Schelkunoff [19] where for a for each pi-circuit. Curves 5, 7 are obtained as 4, 6 however
solid conductor with radius α and internal characteristics σs , the parameters put in every elementary pi-circuit are
µc and εc it is for p.u. length: calculated for a 1m long electrode.
jω c I 0 (γ cα ) In order to avoid drawbacks caused by assumptions in the
Zi = (7)
2παγ c I1 (γ cα ) middle stages of the development of the mathematical
equations (2) it is proposed here to use (3) and (4) for
calculation of R, L, C, G, for the whole length of the electrode B. Voltage and Current distributions along the grounding
and then divide to find the per unit length values. electrode
In this paragraph the transient response of grounding
electrode in terms of voltage and current distributions along
its length is calculated solving the equations of propagation
directly in time domain.

B.1.Impulse excitation model

Transient response of grounding electrodes is calculated
considering impulse current excitation of a double
exponential waveform. It is:
(
I source (t ) = I o ⋅ e at − e βt ) (9)
Where α and β are real negative constants.

B.2.Derivation of Mathematical Expressions for Calculation

of the Impulse Response of the electrode

Derivation of mathematical expressions for the proposed

methodology is presented in detail in [1]. In this paper the
mathematical background will be briefly described.
Time domain solution of the propagation equations equals
to the sum of a partial solution and the solution of the
homogeneous equation for the total network.

Fig. 2. Voltages and Currents at lumped elements of the equivalent network that
represents grounding electrode.

B.3. Partial solution of the differential equations

By expressing voltages and currents in the middle nodes of

the electrode in Fig. 1, as a function of voltages and currents
at their right side, a differential equation of kth order is
obtained:

i0 = a1 D kVn + a2 D k −1Vn + a3 D k −2Vn + ...ak DVn (10)

where A=R+LD, B=G+CD, D=d/dt k>n, and ai real
constants. Differential equations of the form (10) can be
written for every unknown current i1,… in or voltage V0, V1,
V2,… Vn~1 shown in the circuit of Fig. 1.
For known current i0 (t ) = I 0 e − e ( at
) , equation (10)
βt

accepts the partial solution Vn = I (C e − C e )

0 a ,n
at
β ,n
βt

where Ca ,n , C β ,n are real constants. Assuming the current

βt
distribution along the electrode ik = Ca ,k e − C β ,k e
at
we
Fig. 1. Impedance magnitude of the grounding electrode as “seen” from one end, can uniquely determine constants Ca ,k , C β ,k in order to
versus electrode length, when various calculation formulas are used.
satisfy Kirchoff’s laws. The following expressions are
obtained from calculation of the limit of expressions for
 current and voltage distributions as the number of segments n
tends to infinity:
( )
+C 2 (k ) ⋅ coshγ x 2 (k ) (l − x ) ⋅ e x 2 (k )⋅t }
(14b)
C1 (k ) C2 (k )
sinh (γ a (l − x )) at sinh (γ β (l − x )) βt and
are suitably calculated to satisfy the
I p ( x, t ) =
sinh (γ a l )
e −
sinh (γ β l )
e (11)
boundary condition
V x , x Le C e = 0 (
at any point x of the
)
− cosh(γ a (l − x )) at cosh (γ β (l − x )) βt
Gt
− e electrode. That is to assure velocity of propagation
V p ( x, t ) =
LeCe . Furthermore C1 (k ) = C2 (k ) , so the expressions
e − e + Const ⋅ e C

sinh (γ β l )
e

sinh (γ a l ) =1
(14a) and (14b) are further simplified.

where γα = (Re + aLe )(Ge + aCe )

B.5.Complete time-domain solution of the differential
equations
B.4.General solution of the homogeneous differential
equations
Complete time domain solution of the equations of
Full expressions for voltage and current distributions along
propagation of voltage and current waves along the electrode,
the electrode are the (11) plus the general solutions of the
is:
I (x ,t ) = I h (x ,t ) + I p (x ,t )
homogeneous differential equations. The differential equation
that relates source current to the voltage at the end of the (15)
electrode is of the form: V (x ,t ) = V h (x ,t ) + V p (x ,t )
n +1 n+1
 a   β 
a2  − 2  − β 2  − 2 
2 2

i0 =  β2   α2  Vn
Only a few of the infinite terms in (14) are needed to form
(12) the solution (15) with high accuracy.
β −α
d
where D = , A=R+LD, B=G+CD, β2=-Β-α2,
dt
C. Determination and estimation of the Effective Length
− AB + AB ( AB + 4 )
a2 =
2A In case a sinusoidal current source is applied at the
The general solution of (4) where the right part tends energization end of a grounding electrode, when its length
- (Ge + xC e Re + xLe ) ⋅ sinh( γ x l ) increases, maximum GPR decreases almost continuously.
to , is a linear However, there is an upper limit in the length of the electrode
combination of exponential terms. Each of them corresponds that seriously effects the max GPR value observed in case of a
to a root of the characteristic polynomial. Roots of (12) are: lightning strike. The length of the conductor value, beyond
x = − Ge C e γ r l = kπj , k = ±1,2,3....∞ ⇒ which no serious decrement is observed, is called here
and

r1,2 ( k ) =
− Re C e l 2 − Le Ge l 2 ± (R C l e e
2
− Le Ge l 2 ) 2
− 4 Le C e l 2π 2 k 2
“effective length”. The max.GPR value decrement is constant
versus length increment, and it is rapid for lengths below the
2 Le C e effective length and slow when conductors are longer than the
The general solution to the homogeneous equation for the effective length.
voltage at the end is: Another definition of the effective length determines it as

( )
∞ the length of the horizontal conductor, beyond which the
V h (l ,t ) = ∑ C 1 (k ) ⋅ e x 1 (k )⋅t + C 2 (k ) ⋅ e x 2 (k )⋅t maximum transient voltage at the injection point is length
k =1
(13) independent.
where C1(k), and C2(k) are real constants and k is integer. An example of effective length is given in Fig. 3. The
Considering the propagation of each of these exponential maximum GPR produced by a 9kA 1.4/17µs impulse current
terms we obtain: strike and also by a sinusoidal current with the same
∞
I h (x ,t ) = ∑{C1(k )⋅
(
sinh γ x1 (k ) (l − x ) )⋅e x1 ( k )⋅t magnitude at the injection end of horizontal conductors of
k =1 Z x1 (k ) various lengths are comparatively plotted. The 1.4/17µs
(14a) impulse current wave shape has short rising time to its peak.
+C2 (k ) ⋅
(
sinh γ x2 (k ) (l − x ) )⋅ e x 2 (k )⋅t
} For this reason there is larger deviation of the results from
Z x2 (k ) those obtained when the sinusoidal current is injected.
and
( )
∞
Vh (x ,t ) = ∑{C1 (k ) ⋅ coshγ x1 (k ) (l − x ) ⋅ e x1 (k )⋅t

k =1
 8.0E+5
injection point is decreasing and finally goes to an asymptotic
7.0E+5 value when the length of the conductor is increasing.
6.0E+5 According to the definitions of the effective length of
5.0E+5
horizontal grounding conductor, the maximum transient
voltage at the injection point will not decrease further, when
Voltage/V

4.0E+5 Lightning
Response the length of the conductor exceeds a certain value. It is
3.0E+5
50 Hz Response observed from the simulations that the effective length of the
2.0E+5
horizontal grounding conductor for the given current impulse
1.0E+5 and soil conditions is approximately equal to 40m.
0.0E+0
0 50 100 150 200 250 300

Conductor Length/m

Fig. 3. Max Ground Potential rise values vs. conductor length

Effective length is dependent on frequency or the rise time

of the impulse energization. It can be estimated using a
diagram of |Z| vs. length or |Z| vs frequency. Effective length
of an electrode in 500 Ohm-m soil in 100 kHz is 20 m as it
can be seen in Fig. 4a. In Fig. 4b, |Z| of a 20 m long electrode
vs. frequency has been plotted. It remains to a low value until
100 kHz where it begins to increase. Plot of |Z| of a longer
electrode would show the expected increase to a lower
frequency value, as it would be more reactive. Consequently it
appears that 20m corresponds to the longer electrode with low
impedance in 100 kHz, so that is the effective length value in
100 kHz.
Fig. 5. Voltage distribution vs. time at various points of a 140 m long electrode
buried in low relative permittivity soil (εr= 1).

(a) (b)
Fig. 4. Impedance magnitude of the grounding electrode as “seen” from one end,
versus electrode length, when various calculation equations are used.

III. APPLICATION EXAMPLES Fig. 6. Voltage distribution vs. time at injection (start) point of electrodes of
The proposed model is applied to the transient analysis of various lengths (εr= 4)..
a 140m long electrode with radius 1.5mm, buried in 0.9m in
300 Ohm-m soil. Injection current has a 7/28 µs waveform. The ratio of peak voltage to peak current has been plotted
Voltage profiles at various points of the electrode are shown when a fast 1/4µs (Fig. 7) and a slow 8/20µs (Fig. 8) impulse
in Fig. 5. The soil relative permittivity is considered low so current has been injected to show the influence of effective
the effect of the capacitive component weakens. In this case, length of electrodes placed in soil with various characteristics.
the electrode shows a reactive behaviour. This results in faster In Figs 7 and 8 it can be observed that the faster the
appearance of the voltage peak at the injection point. injection current rises to its maximum value, the shorter the
In Fig. 6 the transient responses of horizontal grounding effective length is.
conductors with different lengths buried in soil with εr=4,
ρ=1000Ωm, and µr=1 are simulated. Impulse current is given
from the formula I(t)=I(e-500000t –e-1000000t). The lengths of the
conductors are 20m, 40m, 50m, 80m, 100m, and 280m. The
radius of the conductors is 4mm and the burial depth is
0.75m. In Fig. 6 the magnitude of the transient voltage at the
 10000 [5] S.Cattaneo, A.Geri, F.Mocci, G.Veca: “Transient Behavior of Grounding
Systems Simulation: Remarks on the EMTP’s and special code’s use”
Proceedings of the 21st EMTP Users Group Meeting, 5-7 June 1992
1000
Kolympari Crete,Greece.
Vmax/Imax

[6] R.Rudenberg: “Electrical shock waves in power systems”, Cambridge

3000 m Harvard Univ Press, 1968, pp.273-277.
100
1000 m [7] E. D. Sunde, "Earth Conduction Effects in Transmission Systems" D. Van
300 m Nostrand Company Inc., New York 194
100 m
10 [8] G.F.Tagg, "Earth Resistances" Pittman Publishing Corporation, London
50 m
1964.
[9] W. Dwight “Calculation of Resistances to Ground” Electrical Engineering,
1
0 20 40 60 80 100
pp.1319-1328, 1936.
Length (m)
[10] N.D. Hatziargyriou, M.Lorentzou: 'Grounding Systems Design Using
Fig. 7. Vmax / Imax of horizontal conductors of various lengths under fast impulse. EMTP', 23rdEuropean EMTP Users Group Meeting, Barcelona, 9-11
November 1997.
10000
[11] M.Lorentzou, N.Hatziargyriou: “Modeling of Long Grounding Conductors
Using EMTP”, IPST '99 - Budapest, 20-24 June 1999.
[12] R.G.Olsen, M.C.Willis : “A comparison of exact and quasi-static methods
1000 for evaluating grounding systems at high frequencies” IEEE Trans., on
Power Delivery, vol. 11, 1996, pp.1071-1081.
Vmax/Imax

[13] F.Menter,' Accurate Modeling of conductors imbedded in Earth with

100 3000 m Frequency Dependent Distributed Parameter Lines' (Proceeding of the 1st
1000 m
EMTP User Group Meeting, 1992)
300 m [14] F.Menter: 'EMTP Based Model for Grounding System Analysis' IEEE
10
100 m Trans on PWRD, Vol.9, no.4, October 1994.
50 m [15] L.Grcev, V.Arnautovski: "Frequency Dependent and Transient Impedance
1
of Grounding Systems: Comparison between Simulation and
0 20 40 60 80 100 Measurement" Lightning & Mountains, Chamonix, France 1997
Length (m)
[16] L.Grcev F.Dawalibi: “An electromagnetic model for transients in
Fig. 8. Vmax / Imax of horizontal conductors of various lengths under slow grounding systems” IEEE Trans. on Power Delivery, vol.5, 1990,
impulse. pp.1773-1781
[17] Devgan S., Whitehead E.R: “Analytical models for distributed grounding
Systems” IEEE Trans. On Power Apparatus and Systems, Vol. PAS-92,
Part III, pp.1763-1770, Sep./Oct. 1973
IV. CONCLUSIONS [18] M.I. Lorentzou, N.D. Hatziargyriou: “Accurate Determination of Current
In this paper transmission line modeling of grounding and Voltage Distributions at Grounding Electrodes Under Surge
th
Conditions”, Proceedings of the 26 International Conference on Lightning
electrodes is presented and used for calculation of the
Protection ICLP, Krakow, Poland, 2-6 September 2002
transient response of grounding electrodes. The method used [19] S.A. Schelkunoff: “The Electromagnetic Theory of Coaxial Transmission
to solve for voltage and current distribution along the Lines and Cylindrical Shields”, Bell System Technical Journal, 1934, 13,
electrode is based on closed form solution of the telegraphy pp.532-579.
[20] L.M. Wedepohl, D.J. Wilcox: “Transient analysis of underground power
equations. Solution is achieved directly in time domain, so transmission systems” Proceedings of the IEE, Col. 120, No.2, February
any transformation to and from the frequency domain is not 1973.
needed. The main advantage is that results are obtained using [21] CYMGRD User’s guide and Reference Manual, March 1995.
only a few terms in the final expressions, making the method
simple and easy to apply. VI. BIOGRAPHIES
In the paper special attention is given to the accuracy of
calculation of per unit length parameters of the electrode.
M.I.Lorentzou received a PhD from the Electrical and Computer Engineering
Finally the effect of effective length of grounding Department, NTUA in 2001. Her research activity concerns analysis of the
electrodes that limits their performance in higher frequency transient behavior of grounding systems, Lightning Protection of Wind Turbines,
energization is demonstrated, for various soil and impulse and calculation of power system transients. She works at the Power System
Plnaning Department of the Hellenic Tranmsission System Operator. She is
excitation characteristics. member of Technical Chamber of Greece.

V. REFERENCES N.D. Hatziargyriou received a diploma in Electrical and Mechanical

Engineering from NTUA, and an MSc and PhD degree from UMIST, U.K. He is
[1] M.I. Lorentzou, N.D. Hatziargyriou, B.C.Papadias: “Time Domain
currently Professor at Power Division of the Electrical Engineering Department
Analysis of Grounding Electrodes Impulse Response”, IEEE Trans. on
of NTUA. His research interests include power system transients, switching
Power Delivery, vol.5, 2003.
transients, dynamic security and renewable energy sources. He is Senior member
[2] R.Velasquez, D.Mukhedkar 'Analytical Modelling of Grounding
of IEEE and member of CIGRE SC 38 and of Technical Chamber of Greece.
Electrodes Transient Behavior' IEEE Transactions on Power Apparatus
and Systems, Vol. PAS-103, no.6 Jun. 1984
[3] C.Mazzetti,G.Veca: "Impulse Behavior of Grounding Electrodes" IEEE
Transactions on Power Apparatus and Systems, Vol. PAS-102, no.9 Sept.
1983
[4] A.P. Meliopoulos, M.G.Moharam: "Transient Analysis of Grounding
Systems" IEEE Transactions on Power Apparatus and Systems.Vol.PAS-
102, no.2, February 1983.