Meshless Methods for Computation of Overhead Transmission Lines Electromagnetic Fields

Abstract

The meshless element-free Galerkin method (EFGM) is used to solve partial differential equations responsible for obtaining the electromagnetic fields generated by a transmission line. For this purpose, a 2D model based in a real transmission line is constructed and simulated. The results obtained by EFGM have good accuracy and they are compared against analytical solution and Finite Element Method in order to verify the effectiveness, advantages and disadvantages of the method.

Index Terms
meshfree/meshless; element-free Galerkin method (EFGM); finite element method; transmission lines

I. INTRODUCTION

Meshless methods are numerical procedures and are often used to solve partial differential equations (PDE). They are applied in different areas such as mechanical, civil and electrical as for example, in electromagnetic field problems. There are a wide variety of numerical techniques suitable for solving PDE, most of them depend on a structured connection between the nodes known as mesh-based techniques. In meshless methods, there is no need to explicitly establish connectivity relationships among the discretization nodes. These methods do not use a mesh, but rather, a cloud of points distributed throughout the problem domain [ ¹ ]–[ ⁴ ].

The element-free Galerkin method (EFGM) will be used in this work. It was originally proposed by [ ⁵ ] and it belongs to a class of global weak-form meshless methods, which are based on the global Galerkin weak-form for the problem equations. The Moving Least Squares approximation is often used to construct shape functions, and a set of background cells are required to evaluate the integrals resulted from the weak-form [ ⁶ ]–[ ⁹ ].

This work consists in the computational modeling of three-phase transmission line (TL) with uniform parameters. The EFGM is used to obtain the values of the electric and magnetic fields. A 2D domain surrounds a section of this line and the fields are evaluated 1 meter above the ground. To verify the accuracy of the Meshless Method in this kind of application, comparisons of EFGM results with FEM and analytical calculations were performed. It was observed that the EFGM has presented satisfactory results without the need of a great concentration of nodes near the conductors, as required by the finite element method (FEM) [ ¹⁰ ].

Calculations of electromagnetic fields are important, as well as calculating the maximum limits allowed on national and international regulations, since numerous researches have been taking place over the last decades, with the purpose of analyzing the relationship between the electric and magnetic fields and possible effects that may affect the health of people exposed to such fields. More specifically, there is a great concern on decreasing the intensity of low frequency magnetic fields, as these fields have a relationship with some cancers, such as childhood leukemia [ ¹¹ ]–[ ¹⁴ ].

We wish to present a new contribution through the modeling of transmission lines, verifying the advantages and limitations in the use of EFGM in the evaluation of electromagnetic fields in TLs, since there are still few studies using EFGM in this line of research.

II. Element-Free Galerkin Method

For the construction of the mathematical model, initially, a fictitious geometry representing the original model is determined. The geometry of the problem represents the domain (Ω) and the boundary problem (Γ). Then, the mathematical equations capable of representing the geometry problem valid in the domain and in the boundaries are determined [ ¹⁵ ].

The methodology does not require an explicit mesh. Only a node distribution across the problem domain, contours, and boundaries without the need of a pre-established connection between the nodes and an unknown function governed by a Boundary Value Problem (BVP) [ ¹⁶ ]. Fig. 1 shows a 2D node distribution at one point x_I = (x, y) over the entire domain.

Fig. 1
Meshfree discretization in a 2D domain.

Each node is associated with a closed subdomain called the node influence domain. The influence domain is responsible for the approximation function around node I, which uses the shape function (Φ_I). It represents the solution of the problem by a linear combination of the shape function and the solution of the nodal values [ ¹⁰ ]. Fig. 2 illustrates the domains of influence that may overlap each other and may have circular or rectangular shapes. The domains of influence need to surround the entire problem domain, once there is no connection between the nodes.

Fig. 2
Representation of the circular domains of influence of each node.

The element-free Galerkin method presented in 1994 by [ ⁵ ] is considered a variational method. It uses moving least squares (MLS) or interpolation moving least squares (IMLS) for the construction of the shape function. Galerkin’s weak-form is applied to develop the system of equations using the cells for numerical integration [ ⁸ ].

The EFGM applies MLS to approximate the u(x) function with u^h (x) It was initially developed by mathematicians for surface adjustments and data regression [ ¹⁷ ]. It is, currently, the most widely used method for constructing shape functions of meshless methods.

The MLS approximation begins with ( 1 ) [ ¹⁶ ]:

(1) $${\text{u}}^{\text{h}}\left(x\right)=\sum _{\text{i}=1}^{\text{m}}{p}_{\text{i}}\left(x_{\text{I}}\right)a_{\text{i}}\left(x\right)=p^{\text{T}}\left(x_{\text{I}}\right)a\left(x\right),$$

where p^T(x_I) is a complete polynomial base, in case of two dimensions [1 x y], m is the number of base functions, x^T =(x,y) is the position vector in two dimensions, and a^T ( x )={a₁(x) a₂(x) … a_m(x)} are the coefficient vectors.

The vector of coefficients, a(x) in ( 1 ) is a function of x . The coefficients of a(x) are minimized using the norm L2 standard [ ¹⁶ ],

(2) $$\text{J}=\sum _{\text{i}=1}^{\text{n}}{\text{W}}_{\text{i}}{\left[{\text{u}}^{\text{h}}\left({\text{x}}_{\text{i}}\right)-{\text{u}}_{\text{i}}\right]}^2.$$

In ( 2 ), n is the number of nodes within the domain, associated with the weight function W_i=W( x - x _i) (i=1, 2 … n), and u_i is defined as the nodal parameter of u at x=x_i.

For related J stationary condition:

(3) $$\frac{\partial \text{J}}{\partial {\text{a}}_{\text{j}}}=0,\phi =1,2,\dots \text{n}$$

which results in the following linear system:

(4) $$A\left(x\right)a\left(x\right)=B\left(x\right)U_s,$$

where Us = (u₁ u₂ u₃ … u_n)^T is the vector responsible for storing the field function nodal parameters on all nodes of the support domain and A ( x ) is known as the momentum matrix defined by [ ¹⁶ ]:

(5) $$A\left(x\right)=\sum _{\text{i}=1}^{\text{n}}{W}_{i}p\left(x_{\text{i}}\right)p^T\left(x_{\text{i}}\right).$$

The matrix B represented in equation (4) is defined as

(6) $$B_{\text{I}}\left(x\right)=p^T\left(x_{\text{i}}\right)W_i$$

and W_i is a diagonal matrix formed from the weight functions.

Solving equation (4) in a(x) gives [ ¹⁶ ]:

(7) $$a\left(x\right)=A^{-1}\left(x\right)B\left(x\right)U_s$$

Substituting equation (7) into equation (1) we have:

(8) $${\text{u}}^{\text{h}}\left(x\right)=\sum _{\text{i}=1}^{\text{n}}{\text{Φ}}_{\text{i}}{\text{u}}_{\text{i}}={\Phi }^T\left(x\right)U_s,$$

where Φ(x) is the vector of the MLS shape function corresponding to the n nodes of the support domain of point x, which can be rewritten as:

(9) $${\Phi }^T\left(x\right)=p^T{A}^{-1}\left(x\right)B\left(x\right).$$

To obtain the derivatives of the u^h(x) function, one must first find the derivatives of the shape function, as described below [ ¹⁶ ]

(10) $${\Phi }^T\left(x\right)={\text{γ}}^T\left(x\right)B\left(x\right),$$

where

(11) $${\text{γ}}^T=p^T{A}^{-1},$$

γ^T can be found according to ( 11 ), provided that matrix A is symmetrical.

Partial derivatives of γ can be obtained by solving ( 12 ) and ( 13 ) [ ¹⁶ ]:

(12) $$A{\text{γ}}_{,\text{i}}=p_{,i}-A_{,i}\text{γ},$$

(13) $${\Phi }^T{}_{,i}\left(x\right)={\text{γ}}^T{}_{,i}B+{\text{γ}}^T{B}_{,i},$$

where i indicates the x or y coordinates, and the comma corresponds to the partial derivative in relation to the spatial derivatives. Partial derivatives of γ can be obtained by solving ( 12 ) and ( 13 ) [ ¹⁶ ].

MLS do not satisfy the Kronecker delta property, Φ_I (x_j) ≠ δ_IJ, so that u^h (x_i) ≠ U_I. That is, the nodal parameters u^h (x) are not equal to the values of u^h (x_I). The MLS approximation depends on all nodal parameters for node I, whose domain of influence involves node x_I. Therefore, to impose boundary conditions on the EFG method, it is necessary to use special techniques such as Lagrange Multipliers or IMLS to modify the weak-form [ ¹⁰ ], [ ¹⁸ ], [ ¹⁹ ].

In this work, a variation of MLS called Interpolating Moving Least Squares (IMLS) was used. IMLS main features consist on making W_i infinite at the points x_I, if necessary to interpolate in these points [ ⁸ ], [ ¹⁷ ].

To meet the properties of the Kronecker delta, the following weight function is used:

(14) $$\text{W}\left({\text{x}}_{\text{I}}\right)=\frac{1}{{\text{r}}^{2\text{n}}+\text{ε}},$$

where W(x_I) is the weight or window function, n is a positive integer and ε must be a positive real number. This weight function becomes infinite at point x_I, thus ensuring that the points near this node have some value and the others approach to zero [ ¹⁹ ].

III. Mathematical Modeling of the Problem

The following subtopics are used as electromagnetic and variable formulations of the mathematical model for the transmission line employing EFGM.

A. Equations for Magnetic Flux Density

The mathematical model of the magnetic field starts from Maxwell’s Ampere equation, given by ( 15 ) in the differential form [ ²⁰ ]:

(15) $$\nabla \times H=J_t+\text{j}\omega D.$$

Since magnetic field analysis is performed from a quasi-static point of view, the displacement current density $\frac{\partial D}{\partial \text{t}}=\text{j}\omega D$ can be neglected because the displacement current j ω D is much smaller than the conduction current density (jω D ≪ J_t ) [ ²¹ ], [ ²² ].

In a homogeneous, linear and isotropic medium, characterized by μ (permeability), ε (permittivity) and σ (conductivity), the constitutive relations are valid. For the magnetic field it is necessary to use the constitutive relation ( 16 ) [ ²¹ ], [ ²³ ]

(16) $$H=\nu B,$$

where B is the magnetic flux density (Wb/m²), H is the magnetic field intensity (A/m) and υ is the material reluctivity, $\text{ν}=\frac{1}{\text{μ}}$ (A.m/Wb) [ ¹⁹ ] [ ²¹ ]–[ ²³ ].

Considering J_t as the present total current

(17) $$J_t=J_s+J_e,$$

where J_s is the current density vector (A/m²) of TL and J_e is the induced current density vector (A/m²). In this case, there is no J_e portion, that is J_t = J_s , the current density in each TL conductor.

(18) $$\nabla \times vB=\nabla \times v\nabla \times A,$$

doing B =∇ × A where A is the magnetic vector potential, and with J_t and A perpendicular to the plane of interest, we have A = A k and J = J k .

Substituting ( 19 ) into ( 18 ) results the Poisson equation (20) [ ²³ ]:

(19) $$B={\text{B}}_{\text{x}}i+{\text{B}}_{\text{y}}j=\frac{\partial \text{A}}{\partial \text{y}}i-\frac{\partial \text{A}}{\partial \text{x}}j,$$

(20) $$\frac{\partial }{\partial \text{x}}\left(v\frac{\partial \text{A}}{\partial \text{x}}\right)+\frac{\partial }{\partial \text{y}}\left(v\frac{\partial \text{A}}{\partial \text{y}}\right)=-{\text{J}}_{\text{t}},$$

Rewriting ( 20 ), results in ( 21 )

(21) $$\nabla \cdot v\nabla \text{A}=-{\text{J}}_{\text{s}}.$$

It is necessary to apply the respective boundary conditions in the problem: Dirichlet (Essential Boundary Condition - EBC) ( 22 ) and Neumann (Natural Boundary Condition) ( 23 ):

(22) $$\text{A}={\text{A}}_{\text{u}}\text{in}{\text{Γ}}_{\text{u}},$$

(23) $$-v\frac{\partial \text{A}}{\partial \text{n}}=\tilde{\text{t}}\text{in}{\text{Γ}}_{\text{t}},$$

where n is the outward normal direction on the Neumann boundary, $\tilde{\text{t}}$ and A_u are the values at the boundaries Γ_t, and Γ_u, respectively and Γ= Γ_t U Γ_u.

In the calculation of magnetic flux density, Fig. 3 (a) illustrates the Dirichlet boundary indicated in red with zero value and J_sA, J_sB, J_Sc are phase conductors in a three phase balanced system. Fig. 3 (b) shows nodes in blue that are located at the surface of one of the conductors. For the magnetic field, A is only imposed in the boundary domain and it is considered null.

Fig. 3
(a) Dirichlet condition null at domain boundary. (b) Distribution of nodes in part of the domain and around the phase-A conductor.

Equations (21) , (22) and (23) represent the mathematical modeling of electromagnetic fields of transmission lines. Since this formulation presents a difficult solution, it is essential to obtain the weak-form and use the Galerkin method.

The variational form of the problem and the use of Galerkin method is presented below:

(24) $${\displaystyle {\int }_{\Omega }^{}\left(-\nabla \cdot v\nabla \text{A}-J_s\left(\text{t}\right)\right)}\text{ω}\partial \Omega =0\forall \text{ω}\in {\text{H}}^1.$$

Where ω is the weight function and H¹ is the Hilbert space of degree 1.

Manipulating the equation, making use of vector identities and applying the divergence theorem:

(25) $$\sum _{i=1}^n{\text{ωi}}_{}\left\{\sum _{j=1}^n\left[{\displaystyle {\int }_{\Omega }^{}\left(\nabla {\text{Φ}}_{\text{j}}\nabla {\text{Φ}}_{\text{i}}{\text{A}}_{\text{j}}\right)}\partial \Omega +{\displaystyle {\int }_{{\Gamma }_{\text{t}}}^{}\left({\text{Φ}}_{\text{i}}\check{\text{t}}\right)}\partial \text{Γ}-{{\displaystyle {\int }_{\Omega }^{}{J}_s\left(\text{t}\right)\text{Φ}}}_{\text{i}}\partial \Omega \right]\right\}=0.$$

After obtaining the variational formulation, the EFGM transfers the problem from the continuous domain to the discrete domain and yields:

(26) $${\text{K}}_{\text{ij}}={\displaystyle {\int }_{\Omega }^{}\nabla }{\text{Φ}}_{\text{j}}\text{ν}\nabla {\text{Φ}}_{\text{i}}\partial \Omega ,$$

(27) $${\text{F}}_{\text{i}}={\displaystyle {\int }_{\Omega }^{}{J}_s}\left(\text{t}\right){\text{Φ}}_{\text{i}}\partial \Omega -{\displaystyle {\int }_{{\Gamma }_{\text{t}}}^{}\left({\text{Φ}}_{\text{i}}\check{\text{t}}\right)}\partial \text{Γ},$$

The linear system is:

(28) $$KA=F,$$

where K and F are given by their equations (26) and (27):

B. Equations for Electric Field

The modeling for obtaining the electric field is similar to the magnetic field formulation, but without the current density term, which in this case is null. As it is a quasi-static model, the following approach is used: [ ¹⁰ ], [ ¹⁶ ], [ ²¹ ], [ ²² ]

(29) $$\nabla .\text{β}\nabla \text{U}=0.$$

Expanding ( 29 ), where U electric potential is the solution vector and considering β described below, results the Laplace equation [ ²³ ]:

(30) $$\frac{\partial }{\partial \text{x}}\left(\text{β}\frac{\partial \text{U}}{\partial \text{x}}\right)+\frac{\partial }{\partial \text{y}}\left(\text{β}\frac{\partial \text{U}}{\partial \text{y}}\right)=0.$$

Applying the divergence on both sides of the Ampère equation (15) , together with equation (30) and using the constitutive relations D =ε E and J =σ E , gives the following PDE:

(31) $$\nabla .(\text{σ}+\text{j}\overline{\text{ω}}\text{ε})\nabla \text{U}=0,$$

where ε is the permittivity of the medium, σ is electrical conductivity, $\overline{\text{ω}}=2\text{πf}$ , where f is the frequency, E is the electric field and D is the electric displacement vector.

Considering $\text{β}=\text{σ}+\text{j}\overline{\text{ω}}\text{ε}$ and applying the respective boundary conditions: Dirichlet ( 32 ) and Neumann ( 33 ):

(32) $$\text{U}={\text{U}}_{\text{u}}\text{in}{\text{Γ}}_{\text{u}},$$

(33) $$-\text{β}\frac{\partial \text{U}}{\partial \text{n}}=\tilde{\text{t}}\text{in}{\text{Γ}}_{\text{t}},$$

where n is the outward normal direction on the Neumann boundary, $\tilde{\text{t}}$ and U_u are the values at the boundaries Γ_t, and Γ_u, respectively and Γ= Γ_t U Γ_u.

The essential boundary condition (EBC) for the electric field computation is zero at the problem boundary and equal to the voltage potentials on the conductor surfaces of each phase ( V_nA, V_nB, V_nC ) with voltage values lagged by 120°, as shown in Fig. 4 (a) and Fig. 4 (b) .

Fig. 4
(a) Null Dirichlet condition at domain boundary. (b) Dirichlet condition on the conductor surface with potential value of phase-A.

The ( 30 ), ( 32 ) and ( 33 ) represent the strong form of the problem. Since this formulation presents a difficult solution, it is essential to obtain the weak-form by the Galerkin method.

Getting the variational solution of the problem and using the Galerkin method:

(34) $${\displaystyle {\int }_{\Omega }^{}\left(-\nabla \cdot \text{β}\nabla \text{U})\right)}\text{ω}\partial \Omega =0\forall \text{ω}\in {\text{H}}^1.$$

Manipulating the equation, making use of vector identities and applying the divergence theorem:

(35) $$\sum _{i=1}^n{\text{ω}}_{\text{i}}\left\{\sum _{j=1}^n\left[{\displaystyle {\int }_{\Omega }^{}\left(\nabla {\text{Φ}}_{\text{j}}\text{β}\nabla {\text{Φ}}_{\text{i}}{\text{U}}_{\text{j}}\right)}\partial \Omega +{\displaystyle {\int }_{{\text{Γ}}_{\text{t}}}^{}\left({\text{Φ}}_{\text{i}}\check{\text{t}}\right)}\partial \text{Γ}\right]\right\}=0\forall \text{ω}\in {\text{H}}^1$$

One arrives to the following linear system, where K is given by the equations:

(36) $${\text{K'}}_{\text{ij}}={\displaystyle {\int }_{\Omega }^{}\nabla }{\text{Φ}}_{\text{j}}\text{β}\nabla {\text{Φ}}_{\text{i}}\partial \Omega ,$$

(37) $${\text{F'}}_{\text{i}}={\displaystyle {\int }_{\Omega }^{}{\text{Φ}}_{\text{i}}}\check{\text{t}}\partial \text{Γ}=0.$$

The linear system is:

(38) $$K\text{'}U\text{'}=F\text{'}$$

where K’ and F’ are given by equations (36) and (37):

C. External Border Truncation

This approach assumes that within a given limit, the potential or normal potential derivative is close to zero. Thus ( 39 ) is valid for Dirichlet or Neumann conditions at the dummy boundary, having null values [ ²¹ ], [ ²² ]:

(39) $$\text{ϕ}=0\dots \text{or}\dots \frac{\partial \text{ϕ}}{\partial \text{n}}=0\text{em}{\text{Γ}}_{\text{u}}.$$

The distance from the outer boundary (Γ_u) must be at least five times the distance to the center of the problem [ ²⁴ ]. In the case of TL, there are two scenarios, in the first one, the height in the center of the tower is considered as the distance between the ground and the most distant conductor in the tower and in the second one, the height is considered between the lowest level of the catenary and the ground.

IV. Results

A real stretch of a TL that was interconnected to the Hydro-Québec Company System was chosen as a case study. The main characteristics of this transmission line are presented in Table I [ ²¹ ], [ ²² ].

TABLE I
Characteristics of Hydro-Québec TL

Data	TL Hydro-Québec
Average Current	4.68 kA
Operating Voltage	500 kV
Number of Phases	3
Number of Conductors per Phase	5
Geometric Mean Radius (GMR)	0.40 m
Conductors Radius	25.16 mm
Conductors Type	CAA Grosbeak
Phase Distance	5.84 m
Height of Measurement Point Above Ground Level	1 m
Vertical Spacing Between Center of Phase A and Ground	10.848 m
Vertical Spacing Between Center of Phase B and Ground	12.891 m
Vertical Spacing Between Center of Phase C and Ground	10.848 m

The TL presented in Table I is approximated by the geometry showed in Fig. 5 where A, B and C are equivalent conductors. To simplify, the shield wires are disregarded, since their effect in the electric evaluation at ground level is negligible [ ²⁵ ]. In the adopted model represented in Fig. 5 , the height considered is from the center of the conductors at the lowest level of the catenary to evaluate the maximum electromagnetic fields of the conductors in relation to the ground level.

Fig. 5
Structure of Hydro-Québec TL 500 kV.

Having the variational model presented in topic III, it is necessary to solve the problem through the EFGM and validate the results through the analytical calculation and the FEM, making the appropriate comparisons, analysis and verification of the error rates.

A. Magnetic Flux Density Result

The magnetic flux density and electrical field obtained with the EFGM, as well as the analytical solution and FEM results are presented in this section. The values were obtained at 1 m above ground level and the results can be seen in Fig. 6 , Fig. 7 and in Table II . The voltage, current and Geometric Mean Radius (GMR) values of the conductors are show in Table I .

Fig. 6
Magnetic flux density at 1 m from the ground (TL 500 kV Hydro-Québec).

Fig. 7
(a) Magnetic field distribution (TL 500 kV Hydro-Québec). (b) Logarithmic scale magnetic field distribution (TL 500 kV Hydro-Québec).

TABLE II
Magnetic field error analysis at 1 m from the ground (TL 500 kV Hydro-Québec).

Magnetic Flux Density		Error Calculation
Analytical	EFGM	Global Error	Relative Error
71.2 μT	68.38 μT	1.47 %	4.2 %

Magnetic fields are evaluated in the transmission line crossing range at 1 m from the ground. In this situation, the field is maximum at the center of the conductors, suffering a decay in the lateral conductors, this effect is observed in Fig. 6 . In Fig. 7a and Fig. 7b , we can see that the values of the magnetic field are high in the vicinity of the conductors and decrease as the distance from the conductors rises.

It is seen that the results found by EFGM are satisfactory. The maximum magnetic flux density value of 68.38 µT is found to be within the limit adopted by the International Commission for Non-Ionizing Radiation Protection (ICNIRP) [ ²⁶ ]. In this case, the proposed method achieved a better result than the FEM, based on the analytical solution.

B. Electric Field Result

The electric field obtained with the EFGM, as well as the analytical solution can be observed in Fig. 8 , Fig. 9 and Table III .

Fig. 8
Electric field at 1 m from the ground (TL 500 kV Hydro-Québec).

Fig. 9
EFGM Electric field distribution along the logarithmic servitude range (TLs 500 kV Hydro-Québec).

TABLE III
Electrical field error analysis at 1 m from the ground (TL 500 kV Hydro-Québec).

Electric Field		Error Calculation
Analytical	EFGM	Global Error	Relative Error
10.0 kV/m	10.5 kV/m	0.9 %	1.3 %

The electric field curve shown in Fig. 8 indicates a higher value in the lateral conductors, decreasing at central conductor, in the measurement range at 1 meter from the ground. It is possible to see that at the center of the TL the errors inherent to the numerical procedure are greater compared to the analytical model.

The distribution of the electric field next to the TL is shown in Fig. 9 . A logarithmic scale (log10) was used to obtain a better visualization, once the fields decrease rapidly.

Table III indicates the maximum electric field and error calculation.

C. Parameter Analysis

EFGM has several configuration parameters such as: number of Gauss points for numerical integration, size of influence domain (d_max) and number of nodes (shown in Fig. 10 and Fig. 11 ). Thus, some simulations are performed in order to certify these parameters and identify a better configuration for the application in TLs.

Fig. 10
Analysis of the variation of the size of the electric field influence domain (TL 500 kV Hydro-Québec).

Fig. 11
Analysis of the number of Gaussian integration points for the electric field (TL 500 kV Hydro-Québec).

According to [ ¹⁴ ] the size of the influence domain is represented again by ( 40 ):

(40) $${\text{d}}_{\text{I}}={\text{d}}_{\text{max}}{\text{C}}_{\text{I}}$$

where d_max is a proportionality constant with values between 2 ≤ d_max ≤ 4, and C_I is the average nodal distance between nodes near the node of interest. For example, if d_max = 2.0, this indicates an influence domain whose radius is 2.0 times the average value of distances between the nodes [ ¹⁹ ].

Fig. 11 shows the approximation of the solution by employing 2 and 4 Gauss points per cell for the numerical integration. It is found that when using two Gauss points the accuracy is satisfactory and the computational cost is approximately the half time compared to four Gauss points.

Table IV indicates the maximum magnetic field error calculation for 2 and 4 Gauss points per cell.

TABLE IV
Magnetic field error analys for sensitivity analysis of number of Gauss.

Magnetic Field			Error Calculation
Analytical	2 points EFGM	4 points EFGM	Relative Error 2 Gauss point	Relative Error 4 Gauss points
71.2 μT	68.38 μT	65.72 μT	7.7 %	4.2 %

D. Case study – Furnas TL

Table V describes the electrical characteristics of the Furnas TL analyzed by [ ²¹ ], [ ²² ].

TABLE V
Geometric characteristics (TL 525 kV Furnas).

Data	TL Furnas
Average Current	675.23 A
Operating Voltage	525 kV
Number of Phases	3
Number of Conductors per Phase	4
Geometric Mean Radius (GMR)	0.35 m
Conductors Radius	31.96 mm
Conductors Type	CAA Ruddy
Phase Distance	7.5 m
Height of Measurement Point Above Ground Level	1 m
Vertical Spacing Between Center of Phase A and Ground	17.975 m
Vertical Spacing Between Center of Phase B and Ground	25.475 m
Vertical Spacing Between Center of Phase C and Ground	17.975 m

Fig. 12 shows the curves calculated by [ ²¹ ], [ ²² ] and by the EFGM for the Furnas TL. Note that all solutions are symmetrical with respect to the x axis from -20 to 20 meters. Such behavior is indicated by the fact that the geometric configurations of the TL under study are balanced. In this case the height considered is the distance from the center of the conductors at the tower to the ground.

Fig. 12
Electric field distribution along the logarithmic servitude range (TL 525 kV Furnas).

Table VI indicates the maximum electric field and error calculation.

TABLE VI
Electrical field error analysis at 1 m from the ground (TL 525 kV Furnas).

Electric Field		Error Calculation
Analytical	EFGM	Global Error	Relative Error
4.42 kV/m	4.14 kV/m	4.79 %	6.3 %

It is seen that the results found by FEM in this case are better than EFGM. The maximum electric field value of 4.42 kV/m is found to be within the limit adopted by the International Commission for Non-Ionizing Radiation Protection (ICNIRP) for occupational exposure that is 8.3 kV/m [ ²⁶ ].

E. Comparison between FEM and EFGM algorithms

Fig. 13 represents the flowchart computational modeling for FEM and EFGM [ ¹⁶ ].

Fig. 13
Flowchart of the EFGM and FEM algorithms [ ¹⁶ ].

The main difference between the two algorithms is that EFMG does not need a mesh, being necessary only the distribution of nodes throughout the problem domain. The shape function can be built for a point of interest and can change through the local nodal distribution. In the FEM, as shown in Fig 13 , the shape function has pre-defined elements and is the same throughout the problem region.

The EFGM algorithm is shown more detailed in Fig.14 [ ¹⁰ ]:

Fig. 14
Flowchart for EFGM.

The step by step to EFGM algorithm is described below:

1. The geometry of the problem and its physical characteristics are defined;

2. Then it is necessary to configure the nodal coordinates;

3. Generate the integration cells;

4. Insert the points of integration, weight, Jacobian for each cell;

5. Determine the domain of influence of each node;

6. Start the loop over the integration points;

   1. Determine the nodes of the neighbor of the Gaussian points.

   2. Determine the shape function and its derivatives.

   3. Store matrix B .

   4. Add the contribution of matrix K .

7. Enforce boundary conditions;

8. System solution, solving for nodal parameters;

9. Draw the graphs.

Table VII compares FEM and EFGM characteristics for the purpose of this work [ ¹⁶ ].

TABLE VII
Differences between fem and efgm

Item	FEM	EFGM
Configuration	Mesh Generation	Node Generation
Shape function	Pre-defined elements	Local support domains
Discretized system	Banded, symmetric	Banded, symmetric
Imposition of essential Boundary condition	Standard and easy	Standard and easy
Adaptive analysis	Difficult for 3D cases	Easier
Commercial software packages availability	Many	Few
Accuracy	Accurate compared to Finite Difference Method (FDM)	More accurate than FEM in some cases
Computation speed	Fast	Slower (not optimized)
Stage of development	Well developed	With many challenging problems

V. CONCLUSIONS

The computational model EFGM has been successfully applied in the evaluation of electromagnetic fields generated by a real TL with four and five conductors per phase. As seen in Fig. 6 , the EFGM indicates a better approximation (to the analytical calculation curve) and a lower error rate when compared to the FEM. However, in Fig. 12 , the FEM modeling is better than EFGM in relation to the analytical curve. This indicates that in some cases the use of EFGM is better than the use of FEM. Therefore, due to the nature of this work, several aspects related to EFGM can be further studied and/or improved. For example: study a non-uniform node distribution throughout the domain region, using specific softwares for node creation, and implement a 3D model for the catenary region between the spans of the towers, with towers being located at different levels.

The use of FEM in calculations of electromagnetic fields is already consolidated. The EFGM method code is not yet optimized and the shape function calculation is more complex, hence, the processing time is longer than the FEM. However, the application of the Meshless technique in TLs is a promising tool, as there is no need for mesh generation as in Finite Elements, thus allowing great flexibility to deal with complex geometry problems. In contrast, Meshless technique has a complexity involved in the construction of shape functions.

Through EFGM, a weak global form is applied, thereby avoiding a number of implications found in meshless methods that make use of local domain. Another characteristic of EFGM is its great accuracy, in some cases greater than the Finite Element Method.

Calculations of electromagnetic fields in the vicinity of the TL are of great importance. With the method presented in this work it was possible to evaluate the electromagnetic field values considering international standards and guidelines.

REFERENCES

• ^[1]
  E. H. R. Coppoli, B. Ramdane, Y. Maréchal, “Meshless Local Radial Point Interpolation Method for Electromagnetic Device Modeling”, 19th International Symposium on Electromagnetic Fields in Mechatronics , Electrical and Electronic Engineering (ISEF), vol. 01, Nancy, France, 2019.
• ^[2]
  B. M. F. Gonçalves, M. M. Afonso, E. H. R. Coppoli, B. Ramdane and Y. Maréchal, “Application of Mixed Finite Element and Natural Element Method in anti-periodic electromagnetic devices”, International Journal o Applied Electromagnetics and Mechanics , vol. 1, pp. 1-10, 2018.
• ^[3]
  A. L. Paganotti, R. R. Saldanha, M. M. Afono, M. A. O. Schroeder, G. C. Ribeiro, “New Arrangements of Overhead Power Line Phase Conductors Achieved by Multiobjective Non Dominated Sorting Evolutionary Method”, 19th International Symposium on Electromagnetic Fields in Mechatronics , Electrical and Electronic Engineering (ISEF), Nancy. 2019.
• ^[4]
  A. F. L. Nogueira, “Analysis of magnetic force production in slider actuators combining analytical and finite element methods”, Journal of Microwaves and Electromagnetic Applications , vol. 10, no. 1, pp. 243-250, 2011.
• ^[5]
  T. Belytschko, Y. Y Lu, L. Gu, “Element-Free Galerkin Methods”, International Journal for Numerical Methods in Engineering , vol. 37, pp. 229-256, John Wiley & Sons, Ltd. Technological Institute, Northwestern University, Evanston, USA, 1994.
• ^[6]
  E. H. R. Coppoli, R. C. Mesquita, and R. S. Silva, “Periodic boundary conditions in element free Galerkin method”, COMPEL , vol. 28, pp. 922-934, Aug. 2009.
• ^[7]
  G. F. Parrreira, E. J. Silva, A. R. Fonseca, R. C. Mesquita, “The Element-Free Galerkin Method in 3-Dimensional Electromagnetic Problems”, IEEE Transactions on Magnetics, COMPUMAG , vol. 42, no. 4, April 2006.
• ^[8]
  U. C. Resende, E. H. R. Coppoli, M. M. Afonso, “A Meshless Approach Using EFG Interpolating Moving Least-Squares Method in 2-D Electromagnetic Scattering Analysis”, IEEE Transactions on Magnetics , vol. 51, pp. 1-4, 2015.
• ^[9]
  E. H. R. Coppoli, Ramdane, Y. Maréchal, M. Afonso, “Meshless Local Radial Point Interpolation Method for Electromagnetic Devices Modeling”, International Symposium on Electromagnetic Fields in Mechatronics , Electrical and Electronic Engineering (ISEF), IEEE Xplore , May 2020.
• ^[10]
  G. C. Ribeiro, “Use of Meshfree Methods for Evaluation of Electromagnetic Fields - Generated by Three-phase Air Transmission Lines”, Belo Horizonte, Brazil. Dissertation (Master). Graduate Program in Electrical Engineering. Federal Center of Technological Education of Minas Gerais Federal (CEFET-MG) and University of São João Del Rei (UFSJ), 2017.
• ^[11]
  A. C. Pinho, “The 2D Electric Field Effect Generated by Lines of Transmission and the Living Ones”, Florianópolis, Brazil. Dissertation (Master). Graduate Program in Electrical Engineering. Federal University of Santa Catarina (UFSC). Florianópoles, 1994.
• ^[12]
  E. L. Cartensen, “Magnetic Fields and Cancer”, IEEE Engineering in Medicine and Biology , 1995.
• ^[13]
  M. Havas, “Biological Effects of Low Frequency Electromagnetic Fields”, Electromagnetic Environments and Health in Buildings . London, England. Spon Press, pp. 535, 2004.
• ^[14]
  V.W. Filho, et. al. “Exposure to Magnetic Fields and Childhood Acute Lymphocytic Leukemia in São Paulo, Brazil”, The International Journal of Cancer Epidemiology, Detection, and Prevention. Cancer Epidemiology, vol. 35, pp. 534–539, Aug. 2011. São Paulo, Brazil.
• ^[15]
  S. A. Viana. “Study of Meshless Methods in Electromagnetic Problem Solving”, Belo Horizonte, Brazil. Dissertation (Master). Graduate Program in Electrical Engineering. Federal University of Minas Gerais (UFMG), 1998.
• ^[16]
  G. R. Liu, Y. T. Gu, “An Introduction to Meshfree Methods and Their Programming”, Sprinter, 2005.
• ^[17]
  P. Lancaster, P. Salkauskas, “Surfaces Generated by Moving Least Square Methods”, Mathematics of Computation , vol. 37, no. 155, pp. 141-58, 1981.
• ^[18]
  J. Dolbow, T. Belytschko, “An Introduction to Programming the Meshless Element Free Galerkin Method”, Evaston, United States of America. Archives of Computational Methods in Engineering , vol. 5, no. 3,pp. 207-241, 1998.
• ^[19]
  E. H. R. Coppoli, “Modeling of Electromagnetic Devices through Meshless Methods”, Petrópolis, Brazil. (Ph.D. Thesis). Graduate Program in Computational Modeling of the National Laboratory for Scientific Computing (LNCC), 2010.
• ^[20]
  J. Jim, “The Finite Element Method in Electromagnetics”, Jhon Wiley & Sons, Inc., 2002.
• ^[21]
  A. A. M. Farah, “Calculation of Electric and Magnetic Fields in Transmission Lines by the Finite Element Method”, Belo Horizonte, Brazil. Dissertation (Master). Graduate Program in Electrical Engineering. Federal Center of Technological Education of Minas (CEFET-MG), 2014.
• ^[22]
  A. A. M. Farah, M. M. Afonso, J. A. Vasconcelos, M. O. Schroeder, “A Finite-Element Approach for Electric Field Computation at the Surface of Overhead Transmission Line Conductors”, IEEE Transactions on Magnetics , December 2017. https://www.researchgate.net/publication/321738833
  » https://www.researchgate.net/publication/321738833
• ^[23]
  C. A. Balanis, “Advanced Engineering Electromagnetics”, New York, United States of America. John Wiley & Sons, Inc., 2ª Edition, 2012.
• ^[24]
  Q. Chen, A. Konrad, “A Review of Finite Element Open Boundary Techniques for Static and Quasi-Static Electromagnetic Field Problems”, IEEE Transactions on Magnetics , vol. 33, no. 1, January 1997.
• ^[25]
  M. P. Sarma, W. Janischewskyj, “Electrostatic Field of a System of Parallel Cylindrical Conductors”, IEEE Transactions on Power Apparatus and System , pp. 1069-1079, July 1969.
• ^[26]
  ICNIRP. International Commission on Non-Ionizing Radiation Protection: “Guidelines for Limiting Exposure to Time-Varying Electric, Magnetic and Electromagnetic fields”, pp. 494-522, vol. 74, April 1998.

Publication Dates

• Publication in this collection
  18 June 2021
• Date of issue
  June 2021

History

• Received
  20 Oct 2020
• Reviewed
  26 Oct 2020
• Accepted
  15 May 2021