G Meunier
email: gerard.meunier@g2elab.grenoble-inp.fr

Q.-A Phan

O Chadebec

J M Guichon

B Bannwarth

R Torchio

Unstructured PEEC Method with the use of Surface Impedance Boundary Condition

Keywords: Volume Integral Formulation, Unstructured PEEC, Surface Impedance Condition, Losses. 1

An unstructured-PEEC method for modelling electromagnetic regions with surface impedance condition (SIBC) is proposed. Two coupled circuits representations are used for solving both electric and/or magnetic effects in thin regions discretized by a finite element surface mesh. The formulation is applied in the context of low frequency problems with volumic magnetic media and coils. Non simply connected regions are treated with fundamental branch independent loop matrices coming from the circuit representation.

Purpose -An unstructured PEEC (Partial Element Equivalent Circuit) method for modelling electromagnetic regions with Surface Impedance Boundary Condition (SIBC) is proposed Design/Methodology/Approach -Thanks to the use of Withney face elements, two coupled circuits representations are used for solving both electric and/or magnetic effects in thin regions discretized by a finite element surface mesh. The air is not meshed.

Finding -The new surface impedance formulation enables the modeling of volume conductive regions in order to efficiently simulate various devices with only a surface mesh.

Research limitations implications -The propagation effects are not taken into account in the proposed formulation Originality/Value -The formulation is original and is efficient for modelling non simply connected conductive regions with the use of SIBC.

The unstructured PEEC SIBC formulation has been validated in presence of volume magnetic nonconductive region and compared with a SIBC FEM approach. The computational effort is considerably reduced in comparison with volume approaches.

Introduction

A major interest in using Volume Integral Methods for 3D magnetic and electric fields analysis is that the air region does not need to be meshed. Moreover, the development of new efficient matrix compression algorithms (e.g. FMM, ACA, which greatly improve memory storage and resolution time of fully dense matrix systems) renewed the interest of solving Maxwell's equations with integral methods based on the Green's function. In this work, we propose a new formulation dedicated to electromagnetic regions treated by the surface impedance condition [START_REF] Yuferev | Surface Impedance Boundary Conditions in Terms of Various Formalisms[END_REF]. This work is an extension of a previous volume integral approach based on 3D facet interpolations of the current density and of the magnetic flux density [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF], [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF] and SIBC approaches in the case of conductive but nonmagnetic regions [START_REF] De Grève | A Mixed Surface Volume Integral Formulation for the Modeling of High Frequency Coreless Inductors[END_REF]. The main interest is that volume conductive regions only require a surface mesh and that the air is not discretized (unlike finite element approaches combined with SIBC). In others words, only the interface between active and air region need to be meshed. This approach can be compared to Boundary Equation Method (BEM) which also requires the surface mesh of the interface between air and active regions. Some formulations have already been developed in order to take account conductive and magnetic regions [START_REF] Zheng | Three-dimensional eddy current analysis by the boundary element method[END_REF] [START_REF] Hiptmair | Boundary element methods for eddy current computation, in Boundary Element Analysis[END_REF]. These approaches are powerful and are more general than the SIBC formulation presented in this work because they are not based on any approximation. However, BEM implementations are not straightforward and their performances are penalized by the need to integrate complex Green's kernels which is always a difficulty and a can be source of numerical inaccuracies. So, when the skin depth allows the use of SIBC, our integral formulation becomes very attractive in terms of computing performances.

Moreover, thanks to the circuit representation, these formulations can be easily and naturally coupled with external circuits, PEEC-1D cables or coil which is not straightforward in the BEM context.

In this paper, we first present fundamental equations for building the circuit approach formulation for volume con-ducting regions treated with SIBC. The proposed formulation is applied to low frequency problems (capacitive effects and propagation effects will be neglected) in presence of coils and magnetic non conducting regions.

Formulation

Volume Integral Formulation

Let us consider a linear magnetoharmonic problem with electromagnetic regions Ω (with current density J and magnetization M), and source coils Ω0 (with imposed current density J0). Based on the solution of Maxwell equations, electric field E and magnetic field H can be written in terms of the Green's functions, the current density J, and magnetization M by means integral expression over Ω. In the frequency domain, by neglecting propagation effects, and thanks to Lorentz gauge, we have:

In the above expressions, r is the distance between observation point P and the integration point. Volume Integral Equations [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF] are then obtained by matching ( 1) and ( 2) with constitutive relationships E(J) and H(M) inside Ω and Ω , which, in the case of linear properties, are :

𝐄 = 1 𝜎 𝐉 𝐇 = 1 𝜇 -1 𝐌
where σ is the conductivity and µ = µ0µr is the permeability. Different formulations can be obtained by discretizing the regions and by choosing adequate unknowns and interpolation shape functions. The use of 2-form Whitney face interpolation for both current density J and magnetic flux density B leads to the unstructured PEEC method proposed in [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF].

Surface Impedance Boundary Condition (SIBC)

We consider that the frequency is sufficiently high to use the surface impedance first order approximation [START_REF] Yuferev | Surface Impedance Boundary Conditions in Terms of Various Formalisms[END_REF]. Thus, we consider that Then, we have:

𝐊 = 𝐉 𝒔 𝛿 1 + 𝑗 𝚽 = 𝐁 𝐬 𝛿 1 + 𝑗 (3)
From ( 1) and ( 2), the evaluation of electric and magnetic fields can be obtained from K and Φ on Γ. Integral equations are then obtained by matching the values of electric and magnetic fields on the boundaries of the pellicular region. This leads to two surface integral equations with K and Φ as variables:

1 + 𝑗 2𝜎𝛿 𝐊 ≈ -𝑔𝑟𝑎𝑑 𝑉 - 𝑗𝜔𝜇 4𝜋 𝐉 𝟎 𝑟 𝑑Ω - 𝑗𝜔 4𝜋 𝜇 𝐊 𝑟 𝑑Γ + 𝜇 -1 𝜇 𝚽 ∧ 𝒈𝒓𝒂𝒅 1 𝑟 𝑑Γ 1 + 𝑗 2𝜎𝛿 𝚽 ≈ -𝒈𝒓𝒂𝒅𝜑 + 1 4𝜋 𝐉 𝟎 ∧ 𝒈𝒓𝒂𝒅 1 𝑟 𝑑Ω + 1 4𝜋 𝐊 ∧ 𝒈𝒓𝒂𝒅 1 𝑟 𝑑Γ (4) 2.

Use of 2D face element on 3D surface

Region Γ is discretized by surface finite element meshes (composed of triangles or quadrilaterals for instance) on which K and Φ are interpolated with 2D first order surface face elements:

𝐊 = 𝒘 𝒔𝒋 𝐼 𝚽 = 𝒘 𝒔𝒋 𝜙
𝒘 𝒔𝒋 is the face function related to j-th face, whereas 𝐼 and 𝜙 are the current (in A) and magnetic flux (in Wb) flowing through the j-th face, respectively. It is worth notice that, with a surface mesh, face elements are associated to the edges of the mesh).

We have:

𝑑𝑖𝑣 𝒘 𝒔𝒊 = ± 1 𝑠 𝒘 𝒔𝒊 . 𝒏 = ± 1 𝑙 on face i (0 on other faces) ( 5
)
where 𝑠 is the area of the e-th element, and 𝑙 is the length of i-th face (geometric edge). The signs in [START_REF] Nguyen | An independent loops search algorithm for solving inductive PEEClarge problems[END_REF] depend on the choice of the global orientation of the face. These relations ensure the continuity of K.n and 𝚽.n between adjacent elements. Using [START_REF] Nguyen | An independent loops search algorithm for solving inductive PEEClarge problems[END_REF], by considering any scalar potential u and thanks to divergence theorem, we obtain: where 𝛥𝑢 is the difference between both mean values of the potential u on elements sharing face i, as shown in Fig. 2. These properties will be useful to construct an equivalent circuit representation based on the dual mesh, as developed in the next section.

𝒘 𝒊𝒔 𝒈𝒓𝒂𝒅 𝑢 = Δ𝑢 (6)
It is worth notice that the surface divergences of K and Φ are generally different from zero. Indeed, these divergences are related to the current / flux flowing out from Γ. Thus, assuming that the faces of the elements are oriented towards the outside, we have:

𝑑𝑖𝑣 𝐊 = 1 𝑠 𝐼 𝑑𝑖𝑣 𝚽 = 1 𝑠 𝜙

Unstructured PEEC model for thin regions

Applying a Galerkin projection with facet functions wi to equations (4), and by considering (5), we get a circuit representation on the dual mesh [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF]: The equivalent circuit (dual mesh) is made by the branches connecting two adjacent elements (i.e. face of primal mesh), where elements are the nodes of the dual mesh. {ΔVs} and {Δφs} represent the differences between electric and magnetic potentials on the branches of the dual mesh. R and Y are sparse finite element matrices while L, C and C' are fully populated integral matrices.

Δ𝑉 Δ𝜑 = 𝑅 + 𝑗𝜔𝐿 𝑗𝜔 𝐶 𝐶 𝑌 𝐼 𝜙 + 𝑆 𝑇 (7)
In order to take into account magnetic fluxes flowing out from magnetic regions and to consider the presence of capacitive effects in the electric regions (if any), equations ( 5) must be completed. The value of electric and magnetic potentials on each face element of Γ can be obtained from equations [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF]. Considering that div J = 0 in ΩJ, we have :

𝐉. 𝒈𝒓𝒂𝒅 1 𝑟 𝑑𝛺 = 𝐉. 𝒏 𝜞 𝑑Γ = 1 𝑠 𝐼 𝑑𝛤
where 𝐼 is the capacitive current flowing out from element l. Moreover, we have:

𝐌. 𝒈𝒓𝒂𝒅 1 𝑟 𝑑𝛺 = 1 𝜇 𝜇 -1 𝜇 1 𝑠 𝜙 𝑑𝛤
where 𝜙 is the magnetic flux flowing out from element l.

The average value of potentials on each element of Γcan then be written as:

{𝑉 } = 𝑃 {𝐼 } {𝜑 } = 𝑄 {𝜙 } with 𝑃 = 1 4𝜋𝜀 1 𝑠 1 𝑠 1 𝑟 𝑑𝛤 𝑑𝛤 𝑄 = 1 4𝜋𝜇 1 𝑠 𝜇 -1 𝜇 1 𝑠 1 𝑟 𝑑𝛤 𝑑𝛤
{𝜙 } and {Ie} are the magnetic fluxes and electric currents flowing out from the surface elements of Γ to the air. By considering that electric and magnetic potentials are null at infinite, the electric and magnetic circuits are completed by adding branches which connect nodes (i.e. surface elements) of Γ to the infinite (Fig. 3). The solution can then be obtained using a circuit solver based on an independent loops search technique for instance.

Fundamental circuit equations to be solved are:

[𝑀 ]{∆𝑉} = 0 𝑤𝑖𝑡ℎ {∆𝑉} = ∆𝑉 ∆𝑉 𝑀 {∆𝜑} = 0 𝑤𝑖𝑡ℎ {∆𝜑} = ∆𝜑 ∆𝜑
where [𝑀 ] and 𝑀 are the branch-fundamental independent loop matrices of the electric and magnetic equivalent circuit representations, respectively. The unknowns of the system to solve are the mesh currents IM and mesh flux density 𝜙 :

{𝐼} = [𝑀 ]{𝐼 } 𝑤𝑖𝑡ℎ {𝐼} = 𝐼 𝐼 {𝜙} = 𝑀 {𝜙 } 𝑤𝑖𝑡ℎ {𝜙} = 𝜙 𝜙
The final system be solving is:

[𝑀 ][𝑍 ][𝑀 ] [𝑀 ] 𝑍 [𝑀 ] 𝑀 𝑍 [𝑀 ] 𝑀 𝑍 [𝑀 ] 𝐼 𝜙 = 𝑈 𝑈
With:

[𝑍 ] = 𝑅 + 𝑗𝜔𝐿 0 0 𝑃 𝑍 = 𝑗𝜔𝐶 0 0 0 𝑍 = 𝐶 0 0 0 𝑍 = 𝑌 0 0 𝑄 {𝑈 } = [𝑀 ] 𝑆 0 𝑈 = [𝑀 ] 𝑇 0
In the case of only conducting regions, the magnetic unknowns are not needed and the system to solve becomes:

[𝑀 ][𝑍 ][𝑀 ]{𝐼 } = {𝑈 }

Coupling with volume magnetic regions

In this section, capacitive effects are neglected. We present how to adapt the proposed formulation in order to consider volume magnetic regions ΩM (see second problem presented in the section "Results"). In the context of a circuit representation, we have to consider a magnetic equivalent circuit for volume magnetic region which contains internal and external branches [START_REF] Meunier | A magnetic flux -electric current volume integral formulation based on facet elements for solving electromagnetic problems[END_REF]. The magnetic fluxes flowing through internal magnetic branches are related to the flux density by the use of volume face elements, i.e.:

𝑩 = 𝒘 𝒗𝒋 𝛷

where 𝒘 𝒔𝒋 is the volume face function associate to the face j of the volume discretization and 𝛷 the magnetic flux flowing through the j-th face. External branches connect nodes on the external surface of the volume magnetic region with the infinite, as in the case of regions treated by surface impedance condition. The previous equations can then be adapted by considering the fluxes {𝜙 } and {𝜙 } flowing in the magnetic regions, the external fluxes {𝜙 } (treated by surface impedance or not) and surface currents {𝐼 }. We have:

{𝜙} = 𝜙 𝜙 𝜙 {𝐼} = {𝐼 } [𝑍 ] = [𝑅 + 𝑗𝜔𝐿] 𝑍 = [𝑗𝜔𝐶 𝑗𝜔𝐶 0] 𝑍 = 𝐶 𝐶 0 𝑍 = 𝑌 0 0 0 𝑌 0 0 0 𝑄 with 𝑌 = 𝒘 𝒔𝒊 𝒘 𝒔𝒋 𝜇 𝜴 𝑴 𝑑𝛺 𝐶 = 1 4𝜋 𝒘 𝒔𝒊 𝜇 -1 𝜇 𝒘 𝒗𝒋 ∧ 𝒈𝒓𝒂𝒅 𝟏 𝒓 𝜴 𝑴 𝑑Ω 𝑑𝛤 𝐶 = 1 4𝜋 𝒘 𝒗𝒊 𝒘 𝒔𝒋 ∧ 𝒈𝒓𝒂𝒅 𝟏 𝒓 𝚪 𝑑𝛤 𝑑Ω
Note that matrix Q must consider all the external faces of magnetic regions, treated by surface impedance condition or not. If two volume regions are in contact, a small change has to be done to determine Q. Otherwise, we can notice that C and C' matrices need a volume integration on the magnetic region. However, previous works have shown that a relative coarse meshing of magnetic regions leads to good results when volume integral formulation based on circuit representation is used [START_REF] Le-Van | A Volume Integral Formulation Based on Facet Elements for Nonlinear Magnetostatic Problems[END_REF].

In order to get a good accuracy, the integration of Green's kernels of L, P, Q, C matrices are computed by using the analytical integration technique proposed in [START_REF] Fabbri | Magnetic flux density and vector potential of uniform polyhedral sources[END_REF] for the self-interactions.

Non simply connected domains

It is worth notice that non-simply connected domains are naturally treated thanks to the use of a circuit solver (see the second example presented in the section "Results"). In order to efficiently obtain the incidence matrices 𝑀 and 𝑀 , the circuit solver uses a specific algorithm, which combines the use of connectivity matrix between nodes and edges [START_REF] Bossavit | Yee-like schemes on staggered cellulagrids: A synthesis between FIT and FEM approaches[END_REF][7] of the primal mesh and independent loops search algorithm based on the circuit representation [START_REF] Nguyen | An independent loops search algorithm for solving inductive PEEClarge problems[END_REF]. A first connectivity matrix is easily obtained from the finite element mesh and provides a set of small loops (i.e. the loops around each node) which is then provided to the circuit solver to be completed. Then, a specific algorithm builds independent loop matrices by adding missing loops and eliminating superfluous ones. This technic, based on [START_REF] Douglas | A Polynomial-Time Algorithm to Find the Shortest Cycle Basis of a Graph[END_REF], leads to an efficient determination of 𝑀 and 𝑀 , with the advantages to prioritize small loops.

Post-treatment

The computation of the electric and magnetic fields inside the regions can be achieved tanks to the face finite element solution obtained for K and Φ. Fields in air and flux in coils can be obtained thanks to integral computed on the active regions [START_REF] Huang | A Highly Efficient Post-Processing Method for Computing Magnetic Flux Considering Magnetic and Conductive Regions[END_REF]. Finally, thanks to the solution obtained for K, losses are computed by integrating on Γ, i.e. A more complex problem is considered with multiply connected regions as proposed by IEEJ (Fig. 5) [START_REF] Nakata | Comparison of various methods of analysis and finite elements in 3-D magnetic field analysis[END_REF]. The conductivity and relative permeability of both plates are 5.5E7S/m and µr = 100, respectively. The ferrite region has a relative permeability of 3000.

The excitation coil is fed by an alternating current 1000A, 1000Hz. In order to determine the independent loop matrices, the connectivity matrix from the primal mesh is used. For each conductive part, two supplementary loops are found by the circuit solver (Fig. 7) for the electric circuit representation, allowing to take into account the eddy currents flowing around the hole and along the section (Fig. 7), if any. The density of the equivalent surface current K (Fig. 8) is compared with the result obtained by a FEM-SIBC formulation implemented in the Altair Flux TM software (Fig. 9). The difference in terms of Joule losses is about 0.04% for a fine mesh and of 1 % for a coarse mesh (Fig. 10). The new surface impedance formulation enables the modeling of volume conductive and/or magnetic regions in order to efficiently simulate various devices with only a surface mesh. The computational effort is considerably reduced in comparison with volume approaches. The unstructured PEEC SIBC formulation has been tested in presence of volume magnetic nonconductive region and compared with a SIBC FEM approach.

Fig. 1 .

 1 Fig. 1. Typical electromagnetic problem

  the current and flux density in the volume regions are mainly tangential to the boundary and located in a thin pellicular region placed of Ω. Tangential current density J and flux density B along the perpendicular direction (z) of the boundary are express by : Js and Bs being the tangential surface current and flux densitiy on Γ (the external boundary of Ω) respectively, whereas δ is the skin depth. Tangential surface current K (in A/m) and surface flux density Φ (in Wb/m) are obtained by integrating J and B along (z) i.e.:

Figure 2 .

 2 Figure 2. Primal and dual mesh (dotted lines)

Figure. 3 :

 3 Figure. 3: Example of magnetic equivalent circuit (dotted lines) for a cubic volume region treated by SIBC. External branches allow to take into account magnetic fluxes going outside of the region. Analog circuit is used for electric circuit with or without external branches depending of the presence of capacitive effects.Then, the complete set of equations representing the equivalent circuit interpretation of the electromagnetic problem is given by:

  We first present an academic problem composed of a conductive sphere (radius R = 1m, µr = 1, σ = 5.5E+7 S/m) placed in uniform external field B = 1T. Fig.1shows eddy current losses obtained by the analytical solution and the proposed PEEC SICBC formulation.

Fig. 4 .

 4 Fig. 4. Eddy current losses at different frequencies.

Fig. 5 .Fig. 6 .

 56 Fig. 5. The model proposed by IEEJ.

Fig. 8 .

 8 Fig. 8. Surface current K (real part) on conductive regions treated by SIBC

Fig. 9 .Fig. 10

 910 Fig.9. Comparison of the equivalent surface current K.