Introduction

The impulse earthing resistance of earthing electrodes is a very important factor which has to be taken into consideration in the lightning protection design of civil buildings and electric substations. An appropriate choice of protection measures against lightning overvoltage depends to a large degree on the knowledge of the values of the impulse earthing resistances of earthing electrodes. As the simplest form of earthing electrodes, vertical erathing electrodes are widely used in practical earthing systems. The impulse resistances of vertical earthing electrodes have been investigated for many years. Although a few algorithms were presented by different authors [1−3], their modeling for soil ionization effect under high impulse current is still a problem. These previous algorithms usually utilized a linear characteristic to represent the relationship between the electric field and current density in the ionization zone around the earthing electrode. In fact, this relationship has pronounced nonlinearity for typical kinds of soils in terms of the experimental investigation [4]. The nonlinear characteristic should therefore betaken into account in order to perform a more accurate calculation of the impulse earthing resistance.

The aim of this paper is to propose an efficient algorithm for calcu-lating the impulse resistances of vertical earthing electrodes. For the sake of engineering application, the inhomogeneity of soil is neglected. In the algorithm, the vertical earthing electrode is divided into a series of segments and the average potential method is employed to calculate the resistance matrix. The formula for the low current earthing resistance is derived from the resistance matrix. The soil ionization effect under high impulse current is further considered by using a nonlinear characteristic to represent the relationship between the electric field and current density in the ionization zone. With a simplified treatment made for the ionization zone, the effective radius of the equivalent earthing electrode is evaluated from the nonlinear characteristic. The impulse earthing resistance is then obtained by substituting the effective radius into the formula of the low current earthing resistance. Validity of the proposed algorithm has been verified by comparing calculated and measured results.

Low Current Earthing resistance

Consider a vertical earthing electrode in a homogeneous soil with the resistivity ρ, as shown in Figure 1. When a low current I diffuses from it into the soil, the soil ionization cannot be caused. According to the electromagnetic field theory, the presence of the interface between

air and soil can be taken into account by installing the image electrode[5], as illustrated in Figure 2. The image electrode depicted by the dotted line is installed above the earth surface; meanwhile the actual and image electrodes are symmetrical about the earth surface. In consideration of

the diffused current distribution along the electrode length, the actual electrode and its image are divided into N segments, respectively (see Figure 2）. At an arbitrary point on the surface of jth actual segment, as shown in Figure 3, the potential generated by the currents of jth actual segment and its image segment (j′ th segment) is calculated as [5]

${\varphi }_{\text{j}}(z_s)=\frac{\text{ρ}}{4\text{π}}\left[{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}\frac{{\tau }_{\text{j }}dz}{\sqrt{{(z-z_{\text{s}})}^2+r_0^2}}+{\displaystyle {\int }_{{\text{ z}}_{\text{1}}}^{z_2}\frac{{{\tau }^{\prime }}_{\text{j}}d{z}^{\prime }}{\sqrt{{(z^{\prime }-z_{\text{s}})}^2+r_0^2}}}}\right]\text{ (1)}$

where τ_j and τ_j^′ are the linear current densities of jth segment and its image, respectively. As the segment is rather short, the diffused current distribution on the segment is approximately considered to be uniform. Letting I_j and I_j^′ denote the diffused currents of jth segment and its image, respectively,

τ_j and τ_j^′ become τ_j= I_j/(z₂-z₁) and τ_j^′= I_j^′/(z₂-z₁). Owing to the image symmetry, we have I_j= I_j^′. Thus, (2.1) is rewritten as

${\varphi }_{\text{j}}(z_{\text{s}})=\frac{\text{ρ}{I}_{\text{j}}}{4{\text{π(z}}_{\text{2}}-z_1)}\left[{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}\frac{dz}{\sqrt{{(z-z_{\text{s}})}^2+r_0^2}}+{\displaystyle {\int }_{{\text{ z}}_{\text{1}}}^{z_2}\frac{d{z}^{\prime }}{\sqrt{{(z^{\prime }-z_{\text{s}})}^2+r_0^2}}}}\right]\text{ (2)}$

By taking the integral average of φ_j(z_s) over the length of jth segment, the average potential is obtained as

$\begin{array}{l}{\varphi }_{\text{ j}}^{\text{a}}=\frac{1}{{\text{(z}}_{\text{2}}-z_1)}{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}{\varphi }_{\text{j}}(z_{\text{s}})dz}\\ =\frac{\text{ρ}{I}_{\text{j}}}{4{\text{π(z}}_{\text{2}}-z_1{)}^2}\left[{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}\frac{dzd{z}_{\text{s}}}{\sqrt{{(z-z_{\text{s}})}^2+r_0^2}}+{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}{\displaystyle {\int }_{{\text{ z}}_{\text{1}}}^{z_2}\frac{d{z}^{\prime }d{z}_{\text{s}}}{\sqrt{{(z^{\prime }-z_{\text{s}})}^2+r_0^2}}}}}}\right]\\ =\frac{\text{ρ}{r}_{\text{0}}{I}_{\text{j}}}{4{\text{π(z}}_{\text{2}}-z_1{)}^2}\left[2+2\Psi \left(\frac{z_2-z_1}{r_0}\right)-2\Psi \left(\frac{z_1+z_2}{r_0}\right)+\Psi \left(\frac{2z_1}{r_0}\right)+\Psi \left(\frac{2z_2}{r_0}\right)\right]\text{ (3)}\end{array}$

where

$\Psi \left(\xi \right)=\xi {\text{ sinh}}^{-\text{1}}\xi -\sqrt{1+{\xi }^2}$

The self resistance of jth segment is given as

$R_{\text{jj}}=\frac{{\varphi }_{\text{ j}}^{\text{a}}}{I_{\text{j}}}\text{ (4)}$

In a similar manner, the average potential generated by the currents of kth segment and its image on the surface of jth segment, as shown in Figure 4, is found as

$\begin{array}{l}{\varphi }_{\text{jk}}^{\text{a}}=\frac{\text{ρ}{I}_{\text{k}}}{4{\text{π(z}}_{\text{2}}-z_1)\left({\text{z}}_{\text{4}}-z_3\right)}\left[{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}{\displaystyle {\int }_{-{\text{z}}_{\text{4}}}^{-z_3}\frac{dzd{z}_{\text{s}}}{\sqrt{{(z-z_{\text{s}})}^2+r_0^2}}+{\displaystyle {\int }_{-{\text{z}}_{\text{2}}}^{-z_1}{\displaystyle {\int }_{{\text{ z}}_{\text{3}}}^{z_4}\frac{d{z}^{\prime }d{z}_{\text{s}}}{\sqrt{{(z^{\prime }-z_{\text{s}})}^2+r_0^2}}}}}}\right]\\ =\frac{\text{ρ}{r}_{\text{0}}{I}_{\text{k}}}{4{\text{π(z}}_{\text{2}}-z_1)\left({\text{z}}_{\text{4}}-z_3\right)}\left[\Psi \left(\frac{z_3-z_2}{r_0}\right)-\Psi \left(\frac{z_3-z_1}{r_0}\right)+\Psi \left(\frac{z_4-z_1}{r_0}\right)-\Psi \left(\frac{z_3-z_1}{r_0}\right)\right.\\ \left.\Psi \left(\frac{z_2+z_3}{r_0}\right)-\Psi \left(\frac{z_1+z_3}{r_0}\right)+\Psi \left(\frac{z_1+z_4}{r_0}\right)-\Psi \left(\frac{z_4+z_2}{r_0}\right)\right]\text{ (5)}\end{array}$

where I_k is the current of kth segment. The mutual resistance between jth and kth segments is given as

$R_{\text{jk}}=\frac{{\varphi }_{\text{ jk}}^{\text{a}}}{I_{\text{k}}}\text{ (6)}$

By taking the self and mutual resistances as the matrix elements, the resistance matrix of the vertical earthing electrode is formed as

$R=\left[\begin{array}{cccc}{R}_{11}& R_{12}& \cdot \cdot \cdot & R_{1N}\\ R_{21}& R{}_{22}& \cdot \cdot \cdot & R_{2N}\\ ⋮& ⋮& ⋮& ⋮\\ R_{N1}& R_{N2}& \cdot \cdot \cdot & R_{NN}\end{array}\right]\text{ (7)}$

Using R gives the relationship between voltage and current

$\left[\begin{array}{c}{U}_1\\ U_2\\ ⋮\\ U_N\end{array}\right]=\left[\begin{array}{cccc}{R}_{11}& R_{12}& \cdot \cdot \cdot & R_{1N}\\ R_{21}& R{}_{22}& \cdot \cdot \cdot & R_{2N}\\ ⋮& ⋮& ⋮& ⋮\\ R_{N1}& R_{N2}& \cdot \cdot \cdot & R_{NN}\end{array}\right]\left[\begin{array}{c}{I}_1\\ I_2\\ ⋮\\ I_N\end{array}\right]\text{ (8)}$

where U_j (j=1, 2, ..., N) is the voltage of jth segment. As the conducing surfaces are equipotential, we have U₁=U₂= ··· =U_N=U [5]. Calculation of the inverse matrix R⁻¹ gives the following expression

$\left[\begin{array}{c}{I}_1\\ I_2\\ ⋮\\ I_N\end{array}\right]=\left[\begin{array}{cccc}{G}_{11}& G_{12}& \cdot \cdot \cdot & G_{1N}\\ G_{21}& G{}_{22}& \cdot \cdot \cdot & G_{2N}\\ ⋮& ⋮& ⋮& ⋮\\ G_{N1}& G_{N2}& \cdot \cdot \cdot & G_{NN}\end{array}\right]\left[\begin{array}{c}U\\ U\\ ⋮\\ U\end{array}\right]=\left[\begin{array}{c}{\displaystyle \sum _{k=1}^N{G}_{1k}}\\ {\displaystyle \sum _{k=1}^N{G}_{2k}}\\ ⋮\\ {\displaystyle \sum _{k=1}^N{G}_{Nk}}\end{array}\right]U\text{ (9)}$

where G_jk (j, k=1, 2, ..., N) is the element of inverse matrix R⁻¹. In terms of (9), the total current is found as

$I={\displaystyle \sum _{j=1}^N{I}_{\text{j}}}=U{\displaystyle \sum _{j=1}^N{\displaystyle \sum _{k=1}^N{G}_{\text{jk}}}}\text{ (10)}$

As a result, the earthing resistance can be obtained as

$R_{\text{e}}=\frac{U}{I}=\frac{1}{{\displaystyle \sum _{j=1}^N{\displaystyle \sum _{k=1}^N{G}_{\text{jk}}}}}\text{ (11)}$

As shown in Figure 1, a vertical earthing electrode of r₀=0.015m is considered here. The earthing electrode is made of plain carbon steel (con-ductivity σ=0.5×10⁷ S•m^-1). The soil is kind of clay. Its resistivity is ρ=42Ω·m and water content is about 20%. The earthing resistances calculated from (11) are given in Figure 5. The corresponding measured values are given simultaneously in Figure 5, which was obtained by the fall of potential method [6]. As indicated in Figure 5, the calculated values are close to the measured ones.

Impulse Earthing Resistance

For a high impulse current with crest value I, representative of lightning, the impulse resistance R_ie of the earthing electrode is defined as a ratio of the crest value U_m of the impulse potential on the earthing electrode to I, i.e. R_ie=U_m/I [7]. The impulse current can produces great current density and high electric field intensity near the earthing electrode. When the electric field intensity on the surface of the earthing electrode exceeds the critical value E_c of soil ionization gradient, the breakdown will occur. This process can be illustrated by Figure 6 [8].

As the current increases, streamers are developed and in turn arcs are generated. Within the streamer and arc zones, the resistivity decreases from its original value to a limit approaching conductor [8]. In addition, there is a semiconductive zone between the streamer zone and the constant resistivity zone. For the purpose of simplifying calculation, this process can be described by an equivalent model shown in Figure 7. In the equivalent model, the semiconductive zone is neglected since it is small.

The streamer and arc zones are equivalently modeled as an ionization zone [8,9]. As the vertical earthing electrodes used in the actual cases usually have a shorter length (less than 5m), the non-uniformity of the ionization zone along the electrode length, as shown in Figure 7(a), is not pronounced and may be neglected. Therefore, the ionization zone is approximately considered to be uniform along the electrode length, as shown in Figure 7(b) [8−10]. The border of theionizationzoneis delimitedbythecriticalelectric field value E_c. A nonlinear characteris-tic is introduced to represent the relationship between the electric field andcurrent density in theionizationzone. The nonlinear characteristic is given as [4]

$E=a{J}^b\text{ (12)}$

where a and b are constants and J is the current density

$J=\frac{I}{2\text{π}{r}_{}l+2\text{π}{r}_{}^2}=\frac{I}{2\text{π}r\left(r+l\right)}\text{ (13)}$

Substituting (12) into (13), the electric field intensity is rewritten as

$E=a\frac{I^b}{{\left(2\pi \right)}^b{r}^b{\left(r+l\right)}^b}\text{ (14)}$

On the border of the ionization zone, E should satisfy the boundary condition.

${\left.E\right|}_{r=r_i}=E_{\text{c}}\text{ (15)}$

The values of a, b and E_c can be found in [4] for typical kinds of soils. Thus, the boundary radius of the ionization zone can be derived from (14) and (15)

$r_i=\frac{-l+\sqrt{l^2+4Q}}{2}\text{ (16)}$

where

$Q=\frac{I}{2\text{π}}{\left(\frac{a}{E_{\text{c}}}\right)}^{\frac{1}{b}}\text{ (17)}$

The soil ionization is basically equivalent to an increase in the dimension of the earthing electrode, which can be taken into account by the effective radius r_e of the equivalent earthing electrode in constant resistivity zone, asshowninFigure 8. The voltage between the earthing electrode and the border of the ionization zone, as shown in Figure 8(a), is expressed as

$U_i={\displaystyle {\int }_{r_0}^{r_i}\frac{a{I}^b}{{\left(2\text{π}\right)}^b}\frac{1}{r^b{\left(r+1\right)}^b}dr=\frac{a{I}^b}{{\left(2\text{π}\right)}^b}}W\text{ (18)}$

where

$W={\displaystyle {\int }_{r_0}^{r_i}\frac{1}{r^b{\left(r+1\right)}^b}dr}\text{ (19)}$

On the other hand, the voltage between the equivalent earthing electrode and the border of the ionization zone, as shown in Figure 8(b), is expressed as

$U_{\text{e}}={\displaystyle {\int }_{r\text{e}}^{r_i}\frac{\rho I}{2\text{π}r\left(r+l\right)}}dr=\frac{\rho I}{2\text{π}l}\ln \frac{r_i\left(r_{\text{e}}+l\right)}{r_{\text{e}}\left(r_i+l\right)}\text{ (20)}$

Because U_i and U_e must be equal, the effective radius of the equivalent eathing electrode can be determined from (18) and (20)

$r_{\text{e}}=\frac{r_{i}l}{r_i\left(\lambda -1\right)+\lambda l}\text{ (21)}$

where λ=exp(χ) and

$\chi =\frac{al{I}^{b-1}}{\rho {(2\text{π})}^{b-1}}W\text{ (22)}$

As the integration for (19) is difficult to be evaluated analytically, a numerical solution is given to obtain W

$W\approx \frac{\Delta r}{2}{\displaystyle \sum _{n=1}^N\left[\frac{1}{r_{\text{n}-1}^b{\left(r_{\text{n}-\text{1}}+l\right)}^b}+\frac{1}{r_n^b{\left(r_n+l\right)}^b}\right]}\text{ (23)}$

where Δr is the radial step and N=(r_i-r₀)/Δr.

By replacing r₀ with r_e in (11), the impulse earthing resistance can be obtained. To check the validity of the proposed algorithm, two numerical examples are given here. The data for the first example are: r₀=0.0127m, l=3.05m, ρ=87.2Ω·m, E_c=127kV/m, a=3094.6 and b=0.51 [4]. At different current crest values,the calculated impulse earthing resistances are shown in Figure 9.

The data for the second example are: r₀=0.025m, l=1m, ρ=43.5Ω·m, E_c=350kV/m, a=219.05 and b=0.82 [4]. The calculated crest potentials on the earthing electrode (U_m=R_ieI) are shown in Figure 10. Furthermore, the corresponding measured results [11,12] are given in Figures 9 and 10, respectively, for comparison. It can be seen from Figures 9 and 10 that a better agreement appears between calculated and measured results.