A Novel MSFEM Approach Based on the A-Formulation for Eddy Currents in Iron Sheets

The multiscale finite-element method (MSFEM) has been shown to be an efficient tool in the simulation of eddy currents in laminated cores. It allows for the accurate simulation of the electric and magnetic fields without needing to resolve the individual steel sheets in the finite-element mesh. This article presents a novel multiscale formulation for the low-frequency case which corrects some shortcomings of the approach used up to now. It enforces more physically accurate continuity and divergence conditions on the solution while also increasing its accuracy. The advantages of the new approach are studied both analytically and by means of a numerical example.


I. INTRODUCTION
T HE simulation of eddy currents in electrical devices is essential to predict the total losses, as well as to identify regions that are at risk of overheating. Since most industrial iron cores consist of a large number of thin iron sheets, a full simulation using the standard finite-element method (FEM) requires a prohibitive amount of finite elements and therefore too many degrees of freedom (DoFs) in the resulting equation system. This makes its solution a challenging task that requires a lot of time and computational resources [1].
A wide range of homogenization techniques have been developed, which address the problem by replacing the core with a homogenized bulk material that does not have to resolve the individual iron sheets. This approach was used, for example, in [2] and [3] where the parameters of the homogenized material were computed by an auxiliary 1-D problem using one iron sheet as a unit cell.
This work belongs to the multiscale FEM (MSFEM) that separates the solution into its coarse components on the scale of the entire core and its fine components on the scale of the thickness of a single sheet [4].
In [5], an MSFEM approach for the 3-D A-formulation of the eddy current problem is discussed, which is shown to yield satisfactory results at only a fraction of the computational effort required by a full FEM. A quick overview of this approach is given in Section III, followed by an analysis of its remaining shortcomings in Section IV. A novel MSFEM approach is presented in Section V, with Section VI discussing how it is able to overcome the shortcomings of the old approach.
The performance of both methods is compared in a numerical example in Section VII. It is shown that the novel approach is able to approximate the losses at a much higher accuracy Manuscript  while still being far more efficient than a full finite-element computation.

II. A-FORMULATION
Consider the eddy current problem in the frequency domain and the material relations given by the equations with the electric field strength E, the current density J, the magnetic field strength H, the magnetic flux density B, the magnetic permeability μ, the electric conductivity σ , and the angular frequency ω = 2π f , where f is the frequency. The MSFEM has been shown to work well for both nonlinear materials and materials with hysteresis [6]. The FEM simulation utilizes a magnetic vector potential A ∈ H (curl) that fulfills curlA = B. Using (1) leads to the second-order partial differential equation Multiplication of (2) with a test function v ∈ H (curl) and integration by parts yields the weak formulation: for all v ∈ H (curl). A more detailed discussion of formulations for the eddy current problem can be found in [7].

III. OLD MSFEM FORMULATION
The main idea of MSFEM is to separate the solution A of (2) into its fine-scale components (on the scale of the thickness of a single sheet) and its coarse-scale components (on the scale of the whole core). Using an analytic approximation Left: Fine mesh that resolves the iron sheets (red) and the nonconducting air and insulation layers (blue). Right: Coarse mesh that only distinguishes the air (blue) and the homogenized bulk medium of the core (red).
of the fine-scale components, it is possible to solve for the remaining coarse-scale components on a coarse mesh which does not have to resolve each sheet (see Fig. 1). This leads to a massive reduction in the required DoFs for the resulting system. Note that throughout this article the properties of the methods will be visualized for a 2-D example with an alternating prescribed magnetic field in the y-direction for ease of presentation. All properties also hold in the 3-D case, as shown in the numerical example at the end.
In the old first-order MSFEM approach discussed in [5], A is approximated by where A 0 , A 1 ∈ H (curl), and w 1 ∈ H 1 are defined on the coarse mesh. The microshape function φ 1 is piecewise linear across each sheet and each insulation layer between the sheets, as well as continuous over the material boundaries (see Fig. 4).
In the surrounding air domain, only A 0 is used. To obtain the weak formulation, (4) is inserted in (3) for both the trial function A and the test function v. All material parameters and shape functions are averaged over the thickness of one iron sheet and one insulation layer [5].

IV. SHORTCOMINGS
While the MSFEM formulation (4) has been shown to give good results for a wide range of applications, some of its properties can still be improved.

A. Continuity of J
For the current density J ∈ H (div) normal-continuity holds. This especially means that the z component of J needs to be continuous in the z-direction. More specifically, the z component of J vanishes at the boundaries between the sheets and the insulation layers (see J ref in Fig. 2). As can be seen, while the current density J MS obtained by the old MSFEM approach is able to include the edge effect in an averaged sense, it does not fulfill the required continuity. Even if the edge effect was disregarded, the main component A 0 already yields a current density that is not normal-continuous.

B. Divergence Condition
The current density J also fulfills the condition divJ = 0. The respective equation with the approach (4) reads as Note that, because J MS is not an element of H (div) due to the incorrect continuity conditions, (5) only strictly holds in the interior of the finite elements.
Because the constant function, the shape function φ 1 and its derivative φ 1 are linearly independent, each term in (5) needs to be equal to zero individually. As can be seen, this requires a strong coupling between A 1 and w 1 , which is not directly enforced by the MSFEM approach.

C. Separation of Effects
The function A 1 in (4) can be readily interpreted as the coefficient of the laminar currents. However, A 0 and w 1 cannot be separated from each other. The z component of J MS is given as Numerical results confirm that the term containing φ 1 is insignificant due to the weak dependence of the solution on z. This mainly leaves a weighted sum of A 0,z and w 1 . As can be seen in Fig. 3, both terms individually are much greater than the resulting currents and they almost cancel each other out, with only the coefficient of the edge effect remaining.
Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.
The new microshape functions φ * 0 and φ * 2 are of polynomial order 0 and 2, respectively. The asterisk is introduced to separate them from the microshape functions φ 0 and φ 2 used in the conventional MSFEM setting [8]. The function φ * 0 is chosen to be piecewise constant, fulfilling the conditions that its average is equal to 1 and that σ φ * 0 is continuous (see Fig. 4). The function φ * 2 is the integral of φ 1 with respect to z. To achieve continuity of the tangential trace in an averaged sense, A 0 is split into A 0,air and A 0,core . These components are coupled by the condition on the boundary of the core where n is the outward-pointing normal vector. Here and in the following, overlined expressions are averaged over the thickness of one iron sheet and one insulation layer. Details for the nonlinear and hysteretic case can be found in [9]. Note that the middle term of (7) does not contribute to (8) The complete weak MSFEM formulation is given as: Find A 0 ∈ H (curl) and w 1 , w 2 ∈ H 1 fulfilling given boundary conditions and the coupling condition (8) for all v 0 ∈ H (curl) and q 1 , q 2 ∈ H 1 . Note that the appearance of second other derivatives necessitates the use of either a suitable finite element order or a mixed formulation.

VI. ADVANTAGES
The proposed method is able to overcome the previously listed shortcomings.

A. Continuity of J
Because φ * 2 is zero on the material boundaries, the continuity of the third term in (7) follows immediately. As an additional benefit, due to φ * 2 being of second order, the quality of the edge effect approximation is significantly better (see Fig. 5). Furthermore, φ * 0 was constructed in a way that ensured J 0 would become continuous.

B. Divergence Condition
Due to the explicit coupling of the second and third term in (7), several terms in the divergence of J MS cancel out. The remaining equation is The divergence, therefore, depends only on the components of single functions of the approach, implying that there is no more implicit coupling. While divA 0 and [(∂divw)/∂z] are already numerically close to zero, they may be controlled more rigorously by including respective penalty conditions in the weak formulation.

C. Separation of Effects
In (7), every term has a straightforward physical interpretation. As can be seen in Fig. 6, A 0 is the magnetic vector potential of the mean magnetic field. The function w is the coefficient of the laminar currents and its divergence is the coefficient of the edge effect.

VII. NUMERICAL EXAMPLE
Consider the 3-D laminated core consisting of four conducting square-shaped sheets (see Fig. 7). Each sheet is 0.5 mm thick and the ratio of the iron sheets to the insulation layers is 97.5%. The material is assumed to be linear with a magnetic permeability of μ = 1000μ 0 and an electric conductivity of   σ = 2 MS/m. To ensure the system is regular and avoid using an explicit gauge condition, an artificial conductivity of σ 0 = 1 S/m is prescribed in the insulation layers and the surrounding air domain. The frequency is f = 50 Hz. For the sake of efficiency, the three planes of symmetry were used to model only one-eighth of the total geometry. The problem is excited via the Biot-Savart field of two filament currents with an amplitude of I = 1 A. All simulations were done in the open-source software Netgen/NGSolve [10] using a direct solver for the linear systems. Fig. 8 compares the current densities calculated by a reference finite-element solution, by the old MSFEM approach, and the new one in detail of the cross section of the core. It can be seen that while the laminar currents are similar in all three methods, the old MSFEM approach is unable to resolve the edge effect in a satisfactory way. The new approach would still require higher-order terms to resolve the edge effect perfectly, but is already far better.
This also affects the total losses, as can be seen in Table I. The new MSFEM approach is able to approximate the losses up to an error of less than 1%, while the old approach only achieved 3%-4%. Table I also shows one tradeoff concerning the new approach: Because second derivatives of w are required, it is not possible to use the lowest-order finiteelement spaces for all unknowns, which increases the number of required DoFs. However, as can be seen, the old approach is unable to compete with the new one even at a much higher number of DoFs.

VIII. CONCLUSION
The MSFEM has been shown to be an efficient tool to achieve an approximate eddy current solution in laminated cores at high accuracy. A novel MSFEM was presented which corrects some of the shortcomings of the old approach while further increasing the accuracy. This new approach will be the basis for reliable error estimation for the 3-D eddy current problem in future work.