Impact of Grounding Modeling on Lightning-Induced Voltages Evaluation in Distribution Lines

: Lightning-induced voltages are one of the main causes of shutdown in distribution lines. In this work, attention is focused on the effects of wideband modeling of electric grounding in the overvoltage calculation along insulator strings due to indirect lightning strikes. This study is done directly in the time-domain with the grounding being represented with an equivalent circuit accounting for its dynamics. Results show that the adoption of commonly adopted simpliﬁed grounding models, such as low-frequency resistance, may lead to an underestimation of the overvoltage. According to the results, differences in the order of 25% can be found in some studied cases.


Introduction
Transmission and Distribution Systems are highly affected and damaged by direct and indirect lightning events. Direct events occur when lightning directly strikes the line; such events are hazardous but rare and are typically studied and analyzed in Transmission System (TS). On the other hand, indirect events occur when lightning strikes the ground in the proximity of a power system; these events are much more frequent with respect to direct ones, but the overall voltage induced in the power system is usually much lower. For this reason, indirect events are not of interest for TS since the induced voltages are generally lower than the line Critical FlashOver voltage (CFO), but they are vital when dealing with Distribution Systems (DS), which are characterized by a low CFO.
Given the above, this work presents an evaluation of the impact of grounding modeling on lightning-induced voltage. Thus, the main original contribution of this paper is to include, in the time domain type simulations, an equivalent electric circuit that reproduces the complete frequency response of grounding, with full inclusion of the aforementioned effects. The Hybrid Electromagnetic Model (HEM) is used to determine the wideband grounding frequency response Z(ω) [18,19]. To implement the Z(ω) in silico, the Vector Fitting (VF) technique is applied to generate an equivalent electric circuit that is easily inserted in EMT-type software [20,21]. In the following, the grounding circuit will be implemented in the software developed in [22]. In this paper, as commonly proposed in the IEEE Standard [5], the coupling between the tower and the lightning channel and the coupling between the lightning channel and the grounding electrodes are neglected.
The results illustrate that the induced voltages, considering the grounding modeled via R LF are quite different from those results using Z(ω), with perceptual differences reaching values of around 25%. It is noticeable that the differences increase with the soil resistivity and with the point of occurrence of the lightning (lightning striking closer to the DS increase the perceptual differences) for both first and subsequent return strokes. The paper is organized as follows: Sections 2-4 show the lightning field-to-line coupling problem equations, the tower and the grounding modeling, respectively; while Sections 5 and 6 present the test cases and the results. Section 7 is dedicated to the conclusions.

Induced-Lightning Modeling
The lightning-induced voltages occurring in a DS are here evaluated, recalling the procedure presented in [22,23]. This procedure is usually divided into two steps: (i) the ElectroMagnetic (EM) fields computation and (ii) the field-to-line coupling.

EM Fields Computation
The EM fields are computed analytically considering the approach proposed in [24] and validated in [25]. The method requires as input the knowledge of the channel-base current, the return stroke height and the return stroke velocity. It can be applied both to perfect electric conductor ground and soil characterized by a finite conductivity. The only assumption required is the Transmission Line model for the attenuation of the current along the channel. The main advantage of this approach consists of the possibility of dealing with analytical formulas, which guarantee a fast solution and a low computational effort.

Field-to-Line Coupling
The field-to-line coupling computation is obtained considering the well-known Agrawal model [26], which is here presented in its extended version taking into account the presence of a finite-conducting ground and a multi-conductor line.
where V s i (x, t), I i (x, t) and E inc,x,i (x, t) are the scattered voltage, the current and the tangential component of the exciting electric field (computed in the previous subsection) on the ith conductor at distance x from the beginning of the line. As expressed in Equation (1), the knowledge of the inductance and capacitance matrices (L and C) is required. Please note that ξ g i is the time-domain expression for the ground impedance [27]. The total voltage occurring on the i-th conductor at the point x can be then expressed as the sum of the scattered voltage and the incident voltage, whose value depends on the vertical electric field (computed in the previous subsection).
The proposed methodology is adapted to an EMT-type software (in this framework Simulink-Simscape is used), through the finite-difference time-domain (FDTD) technique. In this case, a second-order scheme is adopted with dt = 10 ns and dx = 9 m, which satisfies the well-known Courant stability condition. Further details can be found in [22].

Tower Modeling
The modeling of the tower is usually neglected in lightning-induced voltages studies. However, in this framework, the tower is included in the model according to [28,29] and is modeled as a lossless transmission line, whose characteristic impedance is: h being the tower height and r the tower radius.

Electrical Grounding Modeling
In this paper, the grounding transient behavior is modeled by HEM [18,19]. This model is an electromagnetic computational method developed for the numerical solution of lightning problems and, according to CIGRÉ [30], it is classified as a hybrid electromagneticcircuit approach. The main motivations for using HEM are as follows: (i) it is accurate and flexible, i.e., it can be used in different types of grounding configuration; (ii) its results have been extensively validated experimentally, such as measurements in TS [18,19], horizontal electrodes [18,31], vertical rods [31,32], and typical substation grounding grids [33] and (iii) it is faster than traditional full-wave methods (without losing accuracy). It is worth mentioning that the usage of HEM has increased significantly recently [30,34,35].
Basically, HEM consists of subdividing the actual system (in this case, electrical grounding) into N small conductive cylindrical segments and, for each segment, the electromagnetic theory is applied. After that, by using the circuit theory, it is possible to obtain a matrix system that computes the wideband response of the electrical grounding. For the sake of clarity, we present a brief overview of HEM below. More details about HEM are described in [18,19].
It is worth commenting on the fact that HEM corresponds to an electromagnetic model developed in the frequency domain. Thus, it is necessary first to determine the frequency spectrum (depending on the phenomenon of interest). After that, the electrical grounding is divided into N segments, where the length of each segment is equal to 10 times its radius (thin wire approximation). A discussion about segmentation length is presented in [36].
Each segment is considered a source of two currents, one longitudinal that flows along the electrode (I L ) and another transversal that flows from the electrode to the surrounding soil (I T ). It is worth noting that I L generates a non-conservative electric field and I T a conservative one. With the aid of the magnetic vector and electric scalar potentials, both voltage drops (∆V) and electric potentials (V) in each pair of segments (transmitter and receiver) are determined. Additionally, double integral equations are established for ∆V and V. These integrals depend on the frequency, geometry, soil parameters and I T and I L distributions. However, the distributions of I T and I L are not known and are integrands of the integrals. From this point on, the Method of Moments (MoM) is applied to solve these integral equations [37]. The effect of the air-soil interface is included using the method of images, similar to [38,39].
The I T and I L distributions considered in this paper are of the piecewise-constant function type [18,19]. MoM makes it possible to transform integral equations into algebraic ones, the solution of which allows determining all the quantities of interest (in the frequency domain), such as I T , I L ,∆V and V distributions; transverse (capacitive and conductive couplings) and longitudinal (resistive and inductive couplings) impedances (self and mutual); electromagnetic field; harmonic grounding impedance (Z(ω)); low-frequency grounding resistance (R LF ), etc.
Also, it has been documented in the literature, over almost one hundred years [40], that the soil is a dispersive medium, i.e., the response is not instantaneous. Several researchers have presented a numerical solution to consider this dispersivity in the frequency domain, such as [31,[40][41][42][43][44][45][46][47][48][49]. According to [50], for the values of conductivity considered in this paper, the frequency-dependent soil electric parameters can be neglected in the computation of the EM fields that illuminate the line. On the other hand, according to [17], the impact of considering the frequency-dependence of the soil in the grounding modeling generates sensible differences. Thus, in this paper, the frequency-dependence of the soil parameters is considered in the grounding modeling. Use is made by the formulations proposed in [31], since a CIGRE Brochure has suggested them [51]. Equations (4) and (5) illustrate the formulation.
where σ( f ) is the frequency-dependent soil conductivity (in mS/m), σ 0 is the low-frequency soil conductivity (in mS/m), r ( f ) is the frequency-dependent soil permittivity and 0 is the vacuum permittivity. To obtain mean results (more details about it in [31]), one should use h(σ 0 ) = 1.26σ −0.73 0 , γ = 0.54 and r∞ = 12. As specified before, the calculation of induced voltages is performed directly in the time domain; however, Z(ω) is a frequency domain quantity. Thus, the well-known Vector Fitting (VF) approach is used for fitting the calculated frequency domain grounding response with rational function approximations [20]. The passivity is enforced by perturbation [21].
Finally, based on the obtained rational function, it is possible to synthesize an electric network that can be promptly included in the time-domain simulation. It is important to note that this electric circuit generates the same frequency response as the harmonic grounding impedance provided by HEM. Thus, it includes reactive and electromagnetic wave propagation effects.

Test Cases
This section presents the test cases related to the comparison between two different grounding modeling, i.e., the low-frequency grounding resistance (R LF ) and the harmonic grounding impedance (Z(ω)).
Let us consider a 1.2 km matched three-phase DS (Figure 1). The three-phase conductors' heights are 10, 11 and 12 m, respectively, while the shield wire height is 14 m. The horizontal distance between each conductor and the shield wire is 2.4 m. The conductors' diameter is 1.83 cm, while the shield wire diameter is 0.72 cm.
The span between each tower is 50 m. To consider the influence of the adjacent towers, a total of five towers are modeled in detail, because the towers in longer distance have little impact on the overvoltage. According to [52], for lightning-related phenomena adjacent towers place a moderate influence. Each tower is 14 m high and with a base diameter of 0.5 m. According to Equation (3), a value of Z c = 244.17 Ω is considered. The propagation velocity along the tower is considered to be 0.8c [53,54]. Moreover, the insulators are modeled taking into account their parasitic capacitance. The value of parasitic capacitance of each insulator is 7.68 nF, and they were calculated considering the information in [55][56][57].
Each tower is grounded with a grounding system as shown in Figure 2. This is a typical configuration for grounding distribution networks in the State of Minas Gerais, Brazil. It consists of three vertical rods 2.5 m long interconnected by a horizontal galvanized steel cable 6 m long. The vertical rods are copper-plated steel, with a diameter of 15 mm.
The equivalent circuit of the system composed of a three-phase distribution line, tower and grounding system is shown in Figure 3.  Two cases for the soil parameters will be considered. The soil conductivity and permittivity will be frequency-dependent according to [31], where σ 0 = 10 mS/m and 1 mS/m, respectively. These cases correspond to two different grounding harmonic responses according to the grounding modeling proposed in Section 4. Figures 4 and 5 show Z(ω) and R LF of the two considered cases. Based on the behaviors described in these figures, it is possible to verify that: (i) grounding can only be represented by R LF in the low-frequency range, where Z(ω) tends to R LF ; (ii) the limit frequency of the low-frequency range increases with a reduction in conductivity; (iii) in the intermediate-frequency range there is a predominance of capacitive behavior of the grounding, verified by the decrease of Z(ω) in relation to the R LF ; (iv) the limit frequency of the intermediate-frequency range also increases with the reduction in conductivity and (v) only in the high-frequency range inductive effect is predominant, mainly for higher conductivity values. Thus, the response of the system under study (DS and grounding) will be a direct function of the frequency spectrum of the electromagnetic signal that requests it. As a consequence, it is expected that the overvoltages in the insulator string are sensitive to grounding modeling.  When we consider a grounding model described by R LF , the implementation in the EMT-type software is trivial, while when we consider the harmonic grounding impedance, it is possible to obtain the synthesis of the electric circuit to be implemented in the EMT-type software by using the approach presented in Section 4.
The general layout of the circuit obtained from the Vector fitting approach is described in Figure 6, while the values of the passive circuit are proposed in Tables 1 and 2 for σ = 10 mS/m and in Tables 3 and 4 for σ = 1 mS/m. It is worth mentioning that these equivalent circuits are mathematical models that have a frequency response very close to Z(ω), but their electrical parameters do not have physical consistency, it is also important mentioning that although some elements may have negative values, the circuit is passive in overall. Hence the existence of negative values for resistance, inductance and capacitance in Tables 1-4 do not mean that the circuit is not passive. The i-index appearing in Tables 1-4 refers to the electrical branch.     To compare the grounding modeling, 12 different tests have been implemented (Table 5), each one differing for the soil conductivity, stroke location and stroke type (first or subsequent). The stroke location is always placed in front of the middle of the line, Figure 7 illustrates the distance between the lightning-channel and the tower under study. It is important to highlight that if the closest point of the DS to the stroke location is in the mid-span, it would, naturally, change the maximum overvoltage but not the perceptual differences between the approaches. The lightning return stroke channel is characterized by a height of 8 km and a speed equal to one-half the speed of light in a vacuum. The channel-base current is modeled as a sum of two Heidler's functions as in Equation (6), with parameters reported in Table 6. The representation of the channel-base current is proposed in Figures 8 and 9.
being  Figure 7. Sketch that illustrates the distance between the tower under study and the lightning stroke Location.

Results
In this section, the results for the test cases of Table 5 are presented, showing the voltage across the phase B insulator string (V insulator in Figure 3) and the voltage difference occurring on the grounding system (V grounding in Figure 3). Figures 10-15 show the results for tests T1-T6, corresponding to a typical first stroke. The main differences in terms of voltage across the insulator can be observed considering a low soil conductivity (Figures 13-15) and near stroke locations (60 m). This is extremely important because the closer the stroke location, the higher (and the more dangerous) the induced voltage. For example, let us consider Test T4 (Figure 13). If we use the lowfrequency grounding resistance (R LF ) as grounding model, the maximum induced voltage across the insulator string is 115.12 kV, while if we consider the harmonic grounding impedance (Z(ω)), which represents in a better way the reality, the voltage is 119.40 kV. This shows how the difference in the modeling could lead to either a fault or not across the insulator strings.
On the other hand, when the harmonic grounding impedance model presents a voltage across the insulator higher with respect to the R LF case, the voltage on the grounding system is lower. This can be explained as follows: let us consider Figure 3; the voltage difference occurring on the insulator string is It is reasonable to assume that the voltage on the conductor does not change in a meaningful way. Considering the two different approaches (based on grounding system modeling), the only difference is the current flowing in the shield wire conductor causing a different coupling with the phase conductor. Even if not negligible, the coupling between conductors does not represent the dominant aspect in the lightning-induced voltages (which is the electric field illuminating the conductor). Consequently, V insulator + V sw is almost constant. The shield wire voltage is: with the same current, V tower is constant in the two cases but V grounding varies because the impedance varies according to Figures 4 and 5 for σ = 10 mS/m and 1 mS/m, respectively. Let us consider the most critical case, i.e., σ = 1 mS/m: from Figure 5 it is clear that for each considered frequency Z(ω) < R LF , thus with the same current the voltage on the grounding system is lower if we consider the harmonic impedance Z(ω) and consequently also V sw is lower. Since V insulator + V sw = constant, if V sw decreases , V insulator increases. This aspect is confirmed in Tests T4-T5-T6, T10-T11-T12.
The results for subsequent strokes can be observed in Figures 16-21. The results are in agreement with the previous ones, confirming a significant increase of the maximum voltage if the equivalent circuit (Z(ω)) is taken into account, especially if the soil conductivity is low. Moreover, the percentage increase considering the harmonic grounding impedance with respect to the low-frequency resistance is much more significant with respect to the first stroke case considering σ = 1 mS/m. It is expected since the subsequent strokes are faster. Thus, it has a higher frequency spectrum (the region where there are the highest differences between R LF and Z(ω)). Additionally, it is important to highlight that after a while, both overvoltages, considering R LF and Z(ω), tend to the same value (first and subsequent strokes). For instance, if we consider the T10 and use the R LF as grounding model, the maximum induced voltage across the insulator string is 66 kV, while if we consider the Z(ω), the voltage is 81 kV. This shows how the modeling difference could lead to either a fault or not across the insulator strings, especially for subsequent strokes.            Finally, Table 7 shows the percentage increase in the maximum voltage across the phase B insulator considering the harmonic grounding impedance (Z(ω)) with respect to the low-frequency grounding resistance (R LF ). According to the previous considerations, the differences are almost negligible if the soil conductivity is high (tests T1-T3 and T7-T9), but they become consistent when the soil conductivity decreases (tests T4-T6 and T10-T12). This behavior is more evident for close stroke location (test T4 and T10).

Conclusions
Lightning induced-voltages are usually computed considering only the low-frequency grounding resistance when one considers the grounding system of the distribution tower. This work presented the impact of two different models for the grounding system of distribution line towers on the lightning-induced voltage on the phase insulators computation. The comparison between the low-frequency grounding resistance (R LF ) and the equivalent circuit corresponding to the wideband grounding frequency response (Z(ω)) shows that considering only R LF may lead to non-negligible underestimation of the maximum induced voltage. This aspect is more evident for subsequent strokes in the case of close stroke locations and low soil conductivities, which represents, by the way, one of the configurations when the lightning-induced voltages on a distribution line are high and potentially dangerous. On the other hand, for high soil conductivity, the differences between the two models are negligible. Future work will extend this analysis to the evaluation of a distribution line lightning performance to check whether this trend is also confirmed when dealing with statistical calculations. Additionally, in future works it is expected to include the coupling between the lightning channel and the grounding electrodes, similarly, as presented in [58].