Numerical treatment of floating conductors based on the traditional finite element formulation

A method to model a conductor with undefined potential (commonly known as floating conductor), is sometimes required in the electric field analysis. This paper presents and compares the main methods to deal with such issue, based on the traditional finite element formulation. The purpose is to guide the reader in the selection of a method under the following criteria: Accuracy, implementation and simplicity. The accuracy of each method was evaluated against the analytic solution of the capacitance matrix for a system of parallel cylindrical conductors. Based on the results of the simulations, the characteristics of the methods and this paper’s criteria, a qualification of each method performance is done. Additionally, particular cases in which a specific method could be the most suitable option to deal with floating conductors, are analyzed. In general, the Virtual Permittivity method exhibits the best performance. However, it is shown that the method’s accuracy is influenced by the rounding errors. This paper proposes an additional consideration on the method formulation in order to verify the accuracy of the results.


Introduction
The capacitance matrix [c] and the inductance matrix [l] are often required in the characterization of transmission lines and the analysis of printed circuit boards PCB, to study the electromagnetic compatibility and the electromagnetic interference of such systems [1].In many cases, it is not possible to obtain an analytical solution for [c] and [l], due to complex geometries and materials.Therefore, numerical methods have an important role in the analysis of electric and magnetic fields.The Finite Element Method FEM is a powerful tool to solve the electric field distribution of multiconductor wire lines, PBC, multilayer structures, etc. [2,3,4].An important feature of FEM, is that it enables to describe with great flexibility the geometry or the media of a given problem [2].Hence, it is ideal to deal with complex materials and shapes.
In the electric field characterization, a method to model an isolated perfect conductor is sometimes required [5].This topic has been object of study in many scientific papers [6,7,8].An isolated conductor is usually known as a floating conductor, and it refers to the case in which the potential of the conductor is not defined.Although it constitutes an equipotential surface, the conductor cannot be defined by a Dirichlet boundary condition.Therefore, traditional FE solvers are not suitable to deal with floating conductors.
This paper presents different methods based on the traditional finite element FE formulation to deal with floating conductors in the electric field analysis [9].The methods are compared based on their accuracy, implementation and simplicity.A short and detailed description of each method is presented.In order to evaluate the accuracy of the methods, the capacitance matrix of a system of parallel cylindrical conductors was estimated by each method and the results were compared with the analytical solution; the procedure to calculate the capacitance matrix is presented for each method.Based on the results of the simulations, the characteristics of the methods and this paper's criteria, a qualification of each method performance is done.Typically, the virtual permittivity method VPM is the preferred option to deal with floating potentials.However, in particular cases, another method could be the most suitable option due to its characteristics.Furthermore, it is shown that the VPM accuracy is influenced by rounding errors.For this reason, an additional consideration on the VPM formulation is proposed in order to avoid inaccurate results.
This paper is presented as follows: Section 2 presents the problem formulation for floating conductors.In section 3, the different methods to deal with floating conductors are described.Section 4 shows the analytical solution and the simulation results of [c] for a configuration of cylindrical conductors.In section 5, the performance of the methods is analyzed.Also, it is shown that the VPM is influenced by rounding errors.Therefore, additional considerations must be taken into account in the method formulation, to verify the accuracy of the results.Finally, the main conclusions are presented.

General Formulation
Fig. 1 shows the model of N perfect electric conductors (PEC) in an electrostatic field.Ω is the solution domain, Ω C the domain of each conductor, Γ is the boundary of each domain with subscript i for i = 1, 2, .., f, ..., N, Ω and f ∈ i is the floating conductors subscript for f = 1, 2, ..., M .
Figure 1: Model of a N conductors system [7].
The model is described by the electromagnetic set of equations presented in Eq. (1).
Where ϕ is the potential function; ρ is the volume charge density; ε is the material permittivity.According to the electromagnetic theory, a PEC in the presence of an electrostatic field becomes an equipotential surface.If the equipotential value is known, a Dirichlet boundary condition ϕ i = u 0i is obtained and Ω C may be excluded from the problem domain because the potential is constant throughout the conductor's volume.A floating potential is referred to the case when the u 0 value is unknown and it makes part of the solution [6,10].In addition, the charge of a PEC is redistributed in the presence of nearby charges due to the electrostatic induction.Nevertheless, law of charge conservation establishes that the initial charge in an isolated conductor must be the same before and after it is exposed to an electric field.Therefore, all the methods to deal with floating potentials are based on the following: A PEC represents an equipotential surface with unknown value and it also satisfies the law of charge conservation with Q i = q 0i .

Capacitance Matrix Approach
The Capacitance Matrix Approach (CMA) is based on these principles, but it deals with floating conductors in an indirect way.This approach proposes to change the charge conditions of the system, instead of modeling the behavior of the floating conductor, in order to solve a problem with traditional Dirichlet boundary conditions.System of equations (2) shows the charges vector [Q] (for N conductors), in terms of the capacitance matrix [c] and the conductors voltage vector (2) For each conductor, Q i or V i is set as a system condition.The other variable is unknown and it might be required as part of the solution.A floating conductor represents the case when Q i is known and V i is unknown.The CMA replaces the charge conditions Q i by its respective u 0i , which only can be done if [c]'s values are known [11,12].It must be emphasized that [c] can be found by any of the proposed methods.However, only a few methods use the CMA to deal with floating potentials.The advantage of the CMA is that [c] contains all the information of the system, so if any voltage or charge condition is changed, [Q] and [V ] can be obtained solving Eq.( 2).For other parameters, like the potential or electric field distributions, a traditional FE with u 0if for i = 1, 2, ..., N must be performed.The drawback of this approach is that at least (N + 1)/2 simulations are required to find [c]'s values.The reason is that using [c]'s symmetry, the number of unknown coefficients is N (N +1)/2; divided by N , which is the number of equations that can be obtained in each simulation for most of the methods, give the minimum number of simulations that must be performed to solve [c].In this paper, a simulation will be referred to the process of finding [Q] for a given [V ] or vice versa, using FEM.

Methods for Traditional FE Formulation
In this section, the traditional FE methods to deal with floating potentials are presented.A concise explanation of each method's formulation is given.
Also, the main advantages and disadvantages of the methods are mentioned.For a detailed overview of each method, the corresponding references can be consulted.

Charge Representation
The Charge Representation Method (CRM) depicts the surface charge of a conductor as an equivalent charge distribution in the conductor's domain [7].Equation (3) shows the charge of the conductor i.
Where ρ ik is the charge density and v ik is the volume for k = 1, 2, ..., nt, with nt being the number of the finite elements in the conductor.The CRM can deal with charge conditions Q i = 0 but it has several disadvantages.Firstly, it does not consider the conductor as an equipotential domain.Secondly, the method depends on the user's expertise to select the number and the location of the simulated charges in order to model the electrostatic induction effect.Finally, there is an increase in the number of nodes due to the inclusion of Ω C .

Homogeneous Boundary Condition
The Homogeneous Boundary Condition Method (HBCM) deals with floating conductors with Q f = 0.This method consists of reducing the FE matrix based on the fact that the equipotential surface of a floating conductor constitutes a single variable with unknown value [10].System of equations ( 4) shows the equivalent FE matrix equation to the Poisson's equation [9].
Where [ϕ] is the vector of unknown potentials, [S] is the overall stiffness matrix and [B] is a vector related to the initial conditions of voltage u 0 and to the charge is related to the number of nodes en that conform the finite element, to the charge distribution ρ k and to the volume v k of the adjacent elements k to the node ι. [u] is related to the known potentials u 0j of the adjacent nodes j to the node ι.
The procedure for each floating conductor f , consist in selecting a reference node r that belongs to the boundary Γ f to represent the unknown potential ϕ f = ϕ r .Once r is selected, the system is reduced in the following way: All the [S] rows that belong to Γ f , are added to the S r row and then these are removed from the matrix.In the same way, the [B] rows that belong to Γ f are added to the B r row and then are removed from the vector.Next, the same procedure is applied to the columns of [S] that belong to Γ f .Once the reduced system is solved, [ϕ] is reconstructed replacing ϕ Γ f = ϕ r .A detailed explanation of this procedure is presented in [7].The advantages of this method are that it guarantees accurate results, and the fact that the conductor domain Ω C is excluded of the problem.On the other hand, the matrix reduction requires more computational resources.In addition, a traditional FE solver for the Poisson's equation can require modifications to perform the matrix reduction.Usually, the user does not want to deal with the FE package code because a major comprehension of the code and the capacity to modify it, are needed.

Virtual Permittivity
The Virtual Permittivity Method (VPM) consist in representing a PEC with Q f = 0, as an extreme form of dielectric material.The representation is accomplished by defining an infinite electric permittivity ε ∞ for the conductor [10].In the numerical approach, infinite means "very large".The principal disadvantage of the VPM, consist in defining the value of ε r in order to obtain accurate results [8].Suggested ε r values are among 10 5 − 10 7 [10,13].However, the selection of ε r is not that simple as it is shown in section 5. Another drawback of this method is that Ω C must be considered.Therefore, the number of nodes increases according to the mesh size and the geometry of the conductor.

Minimum Energy
The Minimum Energy Method (MEM) deals with the floating potentials using the CMA [13].[c] is calculated based on the energy equation of the system [13].Equation (5) shows the energy of the system W .
The equation is given first in terms of the electric field E, and later in terms of the capacitance matrix [c] and the vector of the conductors voltages [V ].The MEM consist in performing N (N + 1)/2 FE simulations, in which W is calculated using the E solution.Then based on the first and third part of Equ.( 5), an equation system must be propose and solve.The simplest approach to find the self-capacitance coefficients c ii , consists in perform N simulations in which all but one conductor are defined with a boundary condition V j = 0 for j = i.Equation ( 6) defines c ii coefficients for this approach.
Where W i is the system energy found in the simulation with potential V i = 0.In the same sense, the boundary conditions V i = 0, V j = 0 and V k = 0, for k = i ∧ k = j, can defined for the conductors in order to calculate the mutual-capacitance coefficients c ij = c ji by Eq. (7).
As it can be observed besides W ij , the self-coefficients c ii and c jj are required to find mutual-coefficient c ij .Although the number of simulations required for this method increases at a rate of N (N + 1)/2, the calculation of W is simple once E is solved.Therefore, this method is useful when N has a low value.

Gauss Law
The Gauss Law Method (GLM) also calculates the capacitance matrix of the system, with a different procedure, in order to deal with floating potentials.In this case, as for most of the methods, an equation system based on Equ.( 2) must be propose and solve.A simple procedure to solve [c] by the GLM, consists on performing N simulations in which the conductors are defined with a boundary condition V k = 0 for k = j.Equation ( 8) is solved by numerical integration to find the charge Q i for i = 1, 2, ..., N in each conductor.
Equation ( 9) calculates the capacitance coefficients c ij for the proposed procedure.
Where Q ij is the charge of each conductor, due to the voltage excitation V j = 0. Using the symmetry of [c] and choosing the right combination of voltage excitation, the number of simulations can be reduced to a minimum of (N +1)/2.In such case, the resulting equation system must be solved (instead of Eq.( 9) ), to calculate [c].Although this method required fewer simulations than the MEM, the definition and the numerical integration of s could be an exhausting task, specially for multiple conductors and 3D geometries.

Simulations and results
To evaluate the accuracy of the proposed methods, the numerical solution and the analytical solution of [c] were compared with a configuration of parallel cylindrical conductors above ground.The test case is shown in Fig. 2 (not scaled).The distance between the center of each conductor and ground is 150, the distance between the center of two adjacent conductors is 100 and the radius of each conductor is r = 6.Since traditional FE simulations were performed, the system of conductors was enclosed by a rectangle of 4200 × 2100 with V = 0 to avoid significant errors.In addition, the number of nodes used in each simulation is 5 × 10 6 in order to obtain reliable results.The simulations were performed on the Partial Differential Equation toolbox of Matlab.In the case of the HBCM, the matrices [p], [e] and [t], were imported from the toolbox to calculate matrices [B] and [S] [14].Regarding the matrices, [p] contains the (x, y) coordinates of the nodes, [e] contains the information of the geometry's edges (used to set the boundary conditions), and [t] contains the indices of the nodes that conform each element.The matrix calculation and the solution of [V ], are done with a code developed on Matlab which performs the matrix reduction presented in section 3.2.The results are expressed in percent error between the analytical solution and the solution obtained by the corresponding method.

Analytical Solution
Fig. 3 shows a general model for N parallel conductors above ground with their respective image conductors.Where H is the distance between the conductor's center and ground, x ij and y ij are the horizontal and the vertical distances between the conductors i and j , respectively, for i = 1, 2, ..., N and j = 1, 2, ..., N .

Maxwell's Potential Coefficient Method is used to calculate the analytical value of [c],
given by Eq. (10). Where: The method treats the system of conductors as a multi-line charge system under the assumption that the surface charges are uniformly distributed around each conductor.This assumption can be made if H r. Nevertheless, r must be considered in [c]'s calculation to obtain an accurate expression [15,16,17].The analytical solution of [c] for the test case is shown in Eq.( 11).

Charge Representation
Based on the CRM, [c] can be found with N simulations as follows: All the conductors will have an initial charge Q ik , where i = 1, 2, ..., N is the subscript of the conductor and k = 1, 2, ..., N is the subscript of the simulation.The voltages V jk are found in each simulation k for j = 1, 2, ..., N .Finally, Eq.( 12) shows the form of the N sets of equations for i = 1, 2, ..., N that must be solved, in order to find the capacitance coefficients c ij .
In the test case, a constant charge density ρ through the domain of the conductors was considered.Therefore, Q ik = πr 2 ρ ik .The [c] value obtained by the CRM in percent error against the analytical solution, is shown in Eq.( 13).

Homogeneous Boundary Condition
The proposed procedure to find [c] by the HBCM, consist in defining the conductors with where k is the subscript of the simulation.In this way, for each simulation k a conductor ik is defined with an initial condition V ik .Thus, the variables that must be found in each simulation k are the charge Q ik and the voltages V jk for j = 1, 2, ..., N .To find Q ik , Equ.( 8) can be used.However, if it is convenient, Q ik 's equation can be replace by Equ.( 5) in each simulation.This is useful for the cases in which calculating the system energy W is easier than calculating the charge Q ik .For k = N simulations, Eq.( 14) shows the form of the N sets of equations for i = 1, 2, ..., N , that must be solved to find [c].
The [c] value obtained by the HBCM percent error rate, is shown in Eq.( 15).

Virtual Permittivity
The same procedure for the HBCM can be used to calculate [c] with the VPM.The relative permittivity used in the simulations was ε r = 1 × 10 7 .The [c] value obtained by the VPM in percent error rate, is shown in Eq.( 16).

Minimum Energy
The procedure to find [c] by the MEM is explained in section 3.4.The [c] value obtained by the MEM in percent error rate, is shown in Eq.( 17).

Gauss Law
The procedure to find [c] by the GLM is explained in section 3.5.The [c] value obtained by the GLM in percent error rate, is shown in Eq.( 18).

Accuracy Analysis
Table .1 shows the RMS error between the analytical solution and each method for the test case.Also, it present the maximum error obtained in each case.Except for the CRM, the results indicate that all methods are a suitable option to approximate the capacitance matrix [c].The CRM presents the highest error order.This fact highlights the importance of a good representation of the charge distributions, even for simple geometries.For this reason, the basic formulation of this method is not recommended to deal with floating potentials.

Selection of a Method
The idea of this paper is to guide the reader in the selection of a method to deal with floating potentials.However, this is not a simple task.First, the performance of each method is influenced by several aspects such as its accuracy, computational resources, grid size, the user's expertise, etc.If all these aspects were considered in the comparison of methods for each particular problem (homogeneous or not homogeneous materials, complex or simple geometries, low or high number of conductors, etc.), the selection would be almost impossible.Additionally, different modifications could be made in the formulation of each method in order to improve its performance.For these reasons, this paper only considers the basic formulation of the methods, and the comparison is made based just on the accuracy, simplicity and implementation of each method.The method's qualification according to its characteristics (shown in section 2), the results of simulations and this paper's criteria, are shown in Table .2.The qualification is done from 1 to 5, being 5 the best rating.For instance, the VPM formulation is the easiest to understand and to implement.Therefore, it obtains 5 in simplicity and implementation.Regarding the accuracy parameter, most of the methods show similar results.Nevertheless, the HBCM obtain 5 because once it is implemented, it is the most reliable of the 5 methods.This means that the chances of inaccurate results due to a wrong selection of a parameter or a miscalculation, are the lowest.On the other hand, the traditional formulation of the VPM obtain 4 because its accuracy is highly influenced by the value of ε r that is chosen (as it is discussed later on this section).The next point to be treated, is the cases in which a specific method is the most suitable option.For instance, the only method presented in this paper that can deal directly with initial conditions Q i = 0, is the CRM.Nevertheless, the accuracy of this method depends on the representation of the charge distributions (main conclusion of section 4).As these are usually unknown, the CRM is not recommended to deal with floating potentials.Therefore, the recommended way to deal with initial conditions Q i = 0, is to calculate the matrix [c].All methods can perform this task.However, the MEM is recommended for problems with a low number of conductors (N ≤ 3).The reason is that the MEM presents the simplest procedure to find [c].Still, the required number of simulations for the MEM increases at a rate of N (N + 1)/2, which is a problem.In addition, if the electric field or the voltage distributions are desired, an extra simulation is required for all the methods based on the CMA (which is the case of the MEM).Therefore, if N value is too high, the remaining options are the HBCM, the VPM and the GLM.The procedure to find [c] is very similar for these methods, thus the selection depends on the user's ability to implement each method.Another important case, is when the matrix [c] is not required and the initial charge of the floating conductors is Q j = 0.For this case, the MEM and the GLM are viable options only if N ≤ 3. Otherwise, the selection must be made between the HBCM and the VPM.In general, these two methods are the best options to deal with floating potentials.The HBCM is a simple method that guarantees the accuracy of the results.However, modifications to the FE package may be required in order to implement it.If the FE package contains this feature or the user can modify the package code, the HBCM should be selected.Often, none of these conditions is met.Therefore, the recommended option to deal with floating is the VPM.The main issue with the VPM, is to define the value of ε r for the floating conductors.A study of how this parameter affects the accuracy of the method, was performed.The first step consisted of the 2D simulation of a parallel plates configuration, as shown in Fig. 4a, with a floating conductor in the middle of the plates.To consider the effect of the geometry, size and position of the floating conductor, six geometries were simulated (see Fig. 4).The location of the floating conductor is given in terms of the coordinates of the representative points of each geometry.For each case, simulations were performed with ε r = 10 m for m = 1, 2, ..., 25.As it guarantees accuracy, the results were compared with the HBCM solution.Fig. 5 shows the magnitude of the maximum potential ϕ error between HBCM and the VPM, for the different ε r 's values.It can be observed that the error order changes with the geometry and the ε r value.Nevertheless, the error has the same behavior for all geometries : The initial error decreases when the ε r value is increased.It reaches a minimum at certain ε r value that depends on the geometry.After this, if the ε r value is increased, the error also increases until a local maximum is reached.The problem is that in some cases, the maximum is equal or higher than the value of the initial error (up to a 10% of the voltage applied to the plates).
Relative permittivity ε r  Based on the error behavior, it is proposed that the method's accuracy is affected by the rounding error.To verify this hypothesis, simulations for the same geometries were performed with precision type single instead of double.The results are shown in Fig. 6.In Fig. 6a, the results with single precision are compared with the HBCM.As it can be observed, the error behavior is very similar to the double precision case (Fig. 5).However, a critical peak value is presented before the error stabilizes for the triangle geometry.This is an additional reason why the selection of ε r should not be taken lightly.The difference in the VPM results, for the precision types: double and single, are shown in Fig. 6b.It can be observed that the rounding error hypothesis is proved.For low values of ε r , both precision types present similar results (error ≤ 10 −1 ).Due to the fact that the single precision works with fewer significant figures, it is affected early by the rounding error.For this reason, a significant difference between the results can be obtained for certain values of ε r .However, as ε r 's value increases, both precision types will be affected by the rounding error.Therefore, at some point both precision types should have similar results; As it is shown in Fig. 6b.In this case, an error ≤ 10 −1 is obtained for an ε r ≥ 10 19 .the right choice.This will depend on the features of the problem, the characteristics of the method and the user's expertise.For instance, if the matrix reduction can be easily performed, the HBCM should be selected above the VPM.Moreover, expert users can try advanced methods like the one proposed in [8].Also, advanced formulations of the proposed methods could be tried.For example: improving the HBCM to deal with Q if = 0, could be an option.

Conclusions
Different methods to deal with floating conductors, were studied.Based on their characteristics, the results of the simulations and this paper's criteria, each method performance was qualified (see Table .2).In general, the Virtual Permittivity method shows the best performance.Nevertheless, it is shown that the method's accuracy is influenced by the rounding errors.Therefore, additional considerations in the method formulation are proposed in order to guarantee accurate results.

Figure 2 :
Figure 2: System of parallel cylindrical conductors.The radius of each conductor C i is r i = 6 for i = 1, 2, 3, 4

Figure 3 :
Figure 3: Model of parallel conductors and their image.

Figure 4 :
Figure 4: Parallel plate configuration and geometries used in the simulations.

Figure 5 :
Figure 5: Maximum ϕ error between the HBCM and the VPM (with precision type Double), for the different ε r 's values .

Table 1 :
RMS error and maximum error between the methods and the analytical solution.