## Abstract

### Purpose

The simulation of eddy currents in laminated iron cores by the finite element method (FEM) is of great interest in the design of electrical devices. Modeling each laminate by finite elements leads to extremely large nonlinear systems of equations impossible to solve with present computer resources reasonably. The purpose of this study is to show that the multiscale finite element method (MSFEM) overcomes this difficulty.

### Design/methodology/approach

A new MSFEM approach for eddy currents of laminated nonlinear iron cores in three dimensions based on the magnetic vector potential is presented. How to construct the MSFEM approach in principal is shown. The MSFEM with the Biot–Savart field in the frequency domain, a higher-order approach, the time stepping method and with the harmonic balance method are introduced and studied.

### Findings

Various simulations demonstrate the feasibility, efficiency and versatility of the new MSFEM.

### Originality/value

The novel MSFEM solves true three-dimensional eddy current problems in laminated iron cores taking into account of the edge effect.

## Keywords

## Citation

Hollaus, K. (2019), "A MSFEM to simulate the eddy current problem in laminated iron cores in 3D", *COMPEL - The international journal for computation and mathematics in electrical and electronic engineering*, Vol. 38 No. 5, pp. 1667-1682. https://doi.org/10.1108/COMPEL-12-2018-0538

## Publisher

:Emerald Publishing Limited

Copyright © 2019, Karl Hollaus.

## License

Published by Emerald Publishing Limited. This article is published under the Creative Commons Attribution (CC BY 4.0) licence. Anyone may reproduce, distribute, translate and create derivative works of this article (for both commercial and non-commercial purposes), subject to full attribution to the original publication and authors. The full terms of this licence may be seen at http://creativecommons.org/licences/by/4.0/legalcode

## 1. Introduction

A laminated core represents a periodic micro-structure which is well suited for the multiscale finite element method (MSFEM). The aim of MSFEMs is to reduce the computational costs of eddy currents in very large laminated iron cores drastically without losing accuracy (Dular, 2008; Hollaus and Schöberl, 2018).

It can be stated that MSFEMs for eddy currents in laminated iron in two dimensions (2D) are well established. Problems in 2D have been solved very satisfactorily using a magnetic vector potential (MVP) ** A**, a current vector potential

**or a the mixed formulation with**

*T***and the current density**

*A***, see (Hollaus and Schöberl, 2018).**

*J*However, MSFEMs for three-dimensional (3D) problems are still far away from being a satisfactory solution. Analyzing the numerical examples in the literature it is very striking to see that there have been no real 3D MSFEM simulations presented up to now. Most of the examples are rotationally symmetric, for instance a toroidal transformer (Gyselinck *et al.*, 2006), or the magnetic flux is parallel to the laminates (Dular, 2008). Both kinds of problems exhibit no magnetic stray fields.

MSFEMs for 3D problems can be divided into methods solving real 3D problems and those considering 2D/1D-problems. The 2D/1D-MSFEMs are based on the assumption that the end effects of electrical machines, i.e. magnetic stray fields can be neglected; thus, each laminate is exposed to the same electromagnetic field distribution and therefore a simulation of a single laminate suffices, (Bottauscio and Chiampi, 2002; Rasilo *et al.*, 2011). A 2D problem is solved essentially reducing the computational costs compared to brute force 3D finite element method (FEM) models (Handgruber *et al.*, 2013; Schöbinger *et al.*, 2018).

The present paper deals with problems where this assumption is not applicable and 3D problems have to be solved. The aim of this work is to present a novel MSFEM for laminated iron stacks in 3D and its universal applicability and efficiency compared to the standard finite element method (SFEM). This MSFEM is based on ** A**. After recalling the analytic solution of eddy currents in an infinite slab the construction of the basic MSFEM approach with

**will be discussed in detail. Different aspects like averaging of coefficients and the edge effect are addressed. The method is capable to consider air gaps and the edge effect too.**

*A*Simulation results obtained by all specific MSFEMs will be shown and compared to reference solutions computed by the standard finite element method (SFEM) demonstrating the versatility of the new MSFEM for problems in 3D. The savings in computational costs using the new MSFEM instead of SFEM are presented at the end.

## 2. The multiscale finite element method – MSFEM

A laminated iron core exhibits two very different scales. The large scale is determined by the overall dimensions of a laminated core, e.g. the length *L* and the height *H* of a transformer core as shown in Figure 1 on the left, and the small (micro-) scale determined by the tiny dimensions of the thickness of the laminates *d* and the width of an air gap *d*_{0} in between, see Figure 1 on the right. The ratio of these scales is about 10^{6}, and thus very large. Modeling each laminate and air gap of large electrical devices would yield a very large finite element model and consequently an extremely large equation system impossible to solve with reasonable computational effort. However, a laminated iron core represents a quasi-periodic structure, this means not strictly periodic because of the finite overall dimensions, which can be exploited by the MSFEM advantageously. One period *p* is composed of *d* and *d*_{0}.

To substantiate the MSFEM approach, the exact solution of eddy currents in an infinite slab is highly relevant and for convenience the main results are summarized in the following (Stoll, 1974). A Cartesian coordinate system is used.

A single component MVP * A = Ae_{x}* is assumed to be selected at the surfaces

*z = ±d/2*of the slab prescribing a magnetic flux per unit length, such that the magnetic field

**points into the**

*H**y*-direction inducing eddy currents pointing into the

*x*-direction, i.e.

*and*

**H**= H**e**_{y}*, respectively, see Figure 2. On the other hand*

**J**= J**e**_{x}*holds.*

**J**= −jωσ**A**Provided the linear problem is given in the frequency domain, the quasi-static magnetic field with the phasor convention *e ^{jωt}* reads as:

Holds. The solution (2) is described by a hyperbolic sine which is an odd function. Therefore, odd Gauss–Lobatto polynomials

*to construct MSFEM approaches with*

_{i}(s)**. The transformation**

*A**s = 2z/d*holds with

*s∈[−1,1]*and

*z∈[−d/2,d/2]*. The micro-shape functions

*ϕ*are extended by zero in [−(

_{i}*d*] and

_{0}+d)/2,−d/2*[d/2,(d*including the air gap, except

_{0}+d)/2]*ϕ*

_{1}which is extended linearly and becomes zero in {−(

*d*)/2}. Figure 3 shows how the micro-shape functions ϕ

_{0}+d)/2,(d_{0}+d*fit into the periodic structure with*

_{i}*d*and

*d*

_{0}.

Thus, the polynomials facilitate the required continuity of the unknown solution, *ϕ _{i}* (–1) = 0 and

*ϕ*(1) = 0 with

_{i}*i*=

*3, 5.*

### 2.1 Construction of a multiscale finite element method approach with A

To construct a MSFEM approach with ** A**, the eddy currents of a reference solution, detailed in Figure 4, are studied.

Based on the eddy current distribution, the approach:

*in (4) describes the large-scale behavior of the solution, whereas the others that of the fine scale. The large-scale behavior takes account of the large eddy current loops induced by the magnetic stray flux perpendicular to the lamination (Figure 5). At the fine scale, the main magnetic flux parallel to the lamination induces eddy currents confined to flow in very narrow loops shown in Figure 6. These currents are assumed to be split into two parts. The laminar part which is parallel to the laminates and represented by the second term in equation (4), the third term includes the edge effect, i.e. the part where the currents turn around to form closed loops (Figure 6).*

**A**_{0}The boundary value problem to be solved is the ECP:

*represents the conducting domain (iron) and Ω*

_{c}_{0}the non-conducting domain (air). The weak form is:

Find

*, where*

**v**_{h}∈V_{0}*V*.

_{h}⊂H(curl,Ω)For a unique solution the regularization with 0<σ_{0}⪡σ is applied (Ledger and Zaglmayr, 2010). The solution of equation (6) with the SFEM serves as reference solution for the MSFEM. To end up with a weak form for the MSFEM, equation (4) becomes the trial function and

*consists of the iron laminates and the air gaps.*

_{m}

**A**_{1},

*w*

_{1}and

*ϕ*

_{1}are restricted to Ω

*, whereas*

_{m}

**A**_{0}is valid in the entire domain Ω=Ω

_{m}∪Ω

_{0}. Dirichlet and thus essential boundary conditions are prescribed by means of

**A**_{0}exclusively, and only natural boundary conditions are provided for

**A**_{1}and

*w*

_{1}. This is especially true for planes of symmetry. To obtain the weak form for the MSFEM, simply speaking,

*and*

**A**_{h}*in the weak form of the SFEM [equation (6)] are replaced by*

**v**_{h}*and*

**A**˜_{h}*, respectively, resulting in:Find*

**v**˜_{h}*(*, where

**v**_{0h},**v**_{1h}, q_{1h})∈V_{0}*U*and

_{h}⊂H(curl,Ω), V_{h}⊂H(curl,Ω_{m})*W*have been selected. The micro-shape function

_{h}⊂H^{1}(Ω_{m})*ϕ*

_{1}is a periodic, piecewise linear and continuous function, i.e.

*ϕ*.

_{1}∈H_{per}(Ω_{m})### 2.2 Averaging of the highly oscillating coefficients

The arising highly oscillating coefficients in equation (8) make the finite element assembling very expensive. To overcome this problem these coefficients are averaged. Here, the stiffness term is treated representing the mass term too. Writing the stiffness term in detail yields:

Analogue operations are carried out also for the mass term of equation (8).

Coefficients
*σ*, σϕ* _{1,z}, σϕ_{1}*, etc. are averaged over the period

*p*=

*d*+

*d*

_{0}:

*λ*means either

*σ*in iron or air. The bar marks averaging.

Highly oscillating coefficients are replaced by averaged ones, whereby the bilinear form and the linear form in

The averaged coefficients are constant and therefore a rather coarse FE-mesh suffices to get an accurate approximation of the solution. In fact, equations (8) and (10) are solved.

### 2.3 Biot–Savart field and multiscale finite element method

Fields due to currents in coils can be considered by the Biot–Savart field. Rearranging of the linear form yields:

*is the Biot-Savart-field*

**h**_{S}*remains as test function for the linear form.*

**v**_{0h}## 3. Numerical example

The single phase transformer shown in Figure 7 is used to study various simulations of different MSFEMs. The core consists of 183 laminates yielding a fill-factor of *k _{f}* = 0.9734. An electric conductivity of

*σ*= 2.0 · 10

^{6}S/m and a relative permeability of

*μ*= 1,000 in the linear case have been selected. The cross-section of a cylindrical coil is shown in Figure 8. It consists of two layers (dark rings), 60 turns per layer. The length of the coil equals 192 mm. The arrangement of the core with the coils exhibit three planes of symmetry.

_{r}A handmade mesh was created by means of hexahedral FEs to simplify the modeling of each laminate for the reference solution. The Biot–Savart field was used to avoid the modeling of the cylindrical coils. Due to the symmetry one eighth of the problem has been considered in the simulations.

## 4. Different simulations

### 4.1 Frequency domain and higher-order multiscale finite element method

We start with the linear case in the frequency domain and show how to cope with small penetration depths making use of higher-order MSFEM (HMSFEM). To this end, the basic approach [equation (4)] is extended by adding higher order micro-scale terms (Hollaus and Schöberl, 2015), leading to:

Figure 9 shows the eddy current losses computed by SFEM and MSFEM. The relative error of MSFEM presented in Figure 10 is obtained by comparing to SFEM results which have been obtained by a brute force finite element model discretizing each steel sheet. The lowest-order MSFEM approach [equation (4)] is valid as long as the variation of the MVP across the laminate thickness *d _{fe}* can be approximated by a linear function well. For decreasing penetration depths

*δ*approach [equation (4)] starts to fail.

By adding higher-order terms, the accuracy is clearly improved. The reason why the fifth-order approach does not show a better accuracy than the third-order one is that the reference solution is not reliable for high frequencies. Reference solutions for an order higher than 3 could not be solved on the available server with 4 times 16 cores (Intel(R) Xeon(R) CPU E7-8867 v3) and 2 TByte RAM.

### 4.2 Time stepping method and time stepping method multiscale finite element method

To deal with nonlinear materials simulations with the time stepping method (TSM) and MSFEM (TSMSFEM) have been carried out using implicit Euler as time stepping scheme and the fixed point method Bíró and Preis (1995) has been exploited to solve the nonlinear system. Iron is highly nonlinear, but assumed to be isotropic. The magnetization curve used in the simulations is determined by measurement points and linear interpolation (Figure 11). The curve is convex-concave. Input currents are selected with 1.0, 2.0 and 3.0 A (peak value) to deal with different states of saturation. Simulations with these currents have been carried out at 50 and 500 Hz. For the reasons of comparability, the eddy current losses in the laminated core presented in Figures 12 to 17 are scaled to the current in the wire of the coils *I* squared.

#### 4.2.1 Results, 50 Hz.

The Figures 12 to 14 show the losses at *f *=* *50 Hz. The agreement of the losses with respect to time obtained by SFEM and MSFEM is excellent. The influence of the saturation due to different input currents is clearly visible.

#### 4.2.2 Results, 500 Hz.

Eddy current losses obtained at *f *=* *500 Hz (Figures 15 to 17) show a clear transient initial phase behavior before the steady state is reached. A very satisfactory agreement between SFEM and MSFEM is obtained.

### 4.3 Harmonic balance method and MSFEM

Most of the sources of ECPs alternate harmonically in time and only the solution of the steady state has to be calculated. However, in case of nonlinear materials, the solution is not harmonic any more, but still periodic. Thus, the solution can be represented as a Fourier series. This can be exploited advantageously by the harmonic balance method (Yamada and Bessho, 1988) or as also called the multi-harmonic ansatz (Bachinger *et al.*, 2005), i.e. a truncated Fourier series expansion at a finite number. Only a few harmonics are required for a sufficiently accurate approximation. That’s why the harmonic balance method is superior to the time stepping method particularly in case of a transient that takes a long time. The harmonic balance finite element method (HBFEM) saves mainly computation time in simulations of large devices with harmonic excitation and nonlinear material properties. A rigorous estimate for the total error due to the use of truncated Fourier series is presented in Bachinger *et al.* (2005). The successful use of HBFEM in simulations of electromagnetic devices in the frequency domain can be found in De Gersem *et al.* (2001) or in Gyselinck *et al.* (2002). A 2D FEM considering the main magnetic flux with a 1D diffusion equation across the lamination and using a multi-harmonic ansatz of the MVP including hysteresis is shown in Bottauscio *et al.* (2000).

The harmonic balance method is combined with the MSFEM (HBMSFEM) to exploit the advantages of both methods.

For nonlinear problems with time harmonic excitation and steady state, the harmonic balance method is preferably used (Bíró and Preis, 2006). The steady state solution * u(x,t)* is periodic in time with period

*T*:

An approximated solution can be written as a truncated Fourier series:

*c*and

*s*for cosine and sine, respectively, and with the upper bound

*N*.Based on the basic MSFEM approach [equation (4)], the HBMSFEM approach can be written as:

with the coefficient functions:

*α*=

*c*,

*s*and

*k∈ℕ, k≤N*holds. The time average

*in equation (15) is not used in this work. The hat indicates truncated Fourier expansion of the HBMSFEM approach.*

**A**_{0}(**x**)To compare the results of TS and HBFEM, the losses obtained by TS are averaged over the first and second period. Simulation results for the losses are summarized for *f *=* *50 and *f *=* *500 Hz in Appendix 1 in Tables I and II. There is a very satisfactory agreement.

## 5. Computational costs

The number of required degrees of freedom (DOFs) is valid for one-eighth of the single phase transformer and is given in Appendix 2 in Tables III. In general, using MSFEM the number of DOFs can be reduced by a factor of about ten for the studied example of the single phase transformer in the present work.

## 6. Discussion

The presented MSFEM fits very well to ECPs in transformers. Electrical machines can be treated in two ways. The assumption that each iron sheet is exposed to the same electromagnetic field pattern is often permitted. In this case, a 2D/1D MSFEM can be exploited advantageously (Schöbinger *et al.*, 2019; Hollaus *et al.*, 2018). However, when the stray field cannot be neglected or its influence is of interest, a method which copes with 3D problems is absolutely necessary. The presented MSFEM should also work for 3D problems of electrical machines.

The reduction of computational costs grows with the number of iron sheets in the laminated core. Although, applying MSFEM reduces large problems essentially, the remaining complexity is still too large to be solved conveniently. Further methods based on MSFEM needs to be developed.

Large equation systems resulting from problems with many very thin iron sheets compared to the overall dimensions of the core are extremely ill-conditioned and except small problems impossible to solve iteratively due to the lack of an appropriate preconditioner. Similarly, an appropriate preconditioner is missing for equation systems from the MSFEM too. Therefore, both the reference problems and the MSFEM problems have been solved using a direct solver.

## 7. Conclusions

Based on the results in this work, it can be concluded that the MSFEM presented here is very powerful because it reduces the complexity of the ECP in laminated iron cores essentially compared to SFEM and discretizing each sheet, copes with any penetration depth, considers the edge effect, allows to include nonlinear material properties in a straightforward way and is capable to deal with real 3D problems.

## Figures

### Figure 3.

Gauss–Lobatto polynomials are the micro-shape functions, scaled s ∈ (Dular, 2008), the laminate is grey

Eddy current losses in W at f = 50 Hz

I in A |
1 | 2 | 3 | |||
---|---|---|---|---|---|---|

Period \TS | SFEM | MSFEM | SFEM | MSFEM | SFEM | MSFEM |

TSM 1 Period^{st} |
3.51 | 3.52 | 8.92 | 8.97 | 13.1 | 13.3 |

TSM 2 Period^{nd} |
3.58 | 3.526 | 9.062 | 9.091 | 13.73 | 13.79 |

HBMSFEM^{a} |
3.979 | 9.112 | 14.05 |

Up to the 5th harmonic.

Eddy current losses in W at f = 500 Hz

I in A |
1 | 2 | 3 | |||
---|---|---|---|---|---|---|

Period \TS | SFEM | MSFEM | SFEM | MSFEM | SFEM | MSFEM |

TSM 1 Period^{st} |
104 | 108 | 401 | 412 | 632 | 646 |

TSM 2 Period^{nd} |
135 | 141 | 537 | 550 | 853 | 869 |

HBMSFEM^{a} |
137 | 539 | 801 |

Up to the 5th harmonic.

No. degrees of freedom DOF

Method | SFEM | MSFEM | |||
---|---|---|---|---|---|

FE order | DOF | FE order | MSFEM order | DOF | |

Time harmonic | 3 | 8,739,144 | 2 | 3 | 310,082 |

Time stepping | 1 | 1,116,860 | 1 | 1 | 103,879 |

Harm. balance^{a} |
LO^{b} |
874,836 | LO^{b} |
1 | 95,256 |

Harm. balance^{a} |
1 | 6,701,160 | 1 | 1 | 623,274 |

Up to the 5th harmonic.

Lowest order Nédélec elements

Source: Schöberl and Zaglmayr (2005).

## Appendix 2. Computational costs

## References

Bachinger, F., Langer, U. and Schöberl, J. (2005), “Numerical analysis of nonlinear multiharmonic eddy current problems”, Numerische Mathematik, Vol. 100 No. 4, pp. 593-616.

Bíró, O., and (1995). and Preis, K. “Finite element calculation of time-periodic 3d eddy currents in nonlinear media”, in Homna, T. (Ed), Advanced Computational Electromagnetics, Elsevier, Budapest Hungary, pp. 62-74.

Bíró, O. and Preis, K. (2006), “An efficient time domain method for nonlinear periodic eddy current problems”, IEEE Transactions on Magnetics, Vol. 42 No. 4, pp. 695-698.

Bottauscio, O. and Chiampi, M. (2002), “Analysis of laminated cores through a directly coupled 2-D/1-D electromagnetic field formulation”, IEEE Transactions on Magnetics, Vol. 38 No. 5, pp. 2358-2360.

Bottauscio, O., Chiampi, M. and Chiarabaglio, D. (2000), “Advanced model of laminated magnetic cores for two-dimensional field analysis”, IEEE Transactions on Magnetics, Vol. 36 No. 3, pp. 561-573.

De Gersem, H., Sande, H.V. and Hameyer, K. (2001), “Strong coupled multiharmonic finite element simulation package”, Compel - the International Journal for Computation and Mathematics in Electrical and Electronic Engineering, Vol. 20 No. 2, pp. 535-546.

Dular, P. (2008), “A time-domain homogenization technique for lamination stacks in dual finite element formulations”, Journal of Computational and Applied Mathematics, Vol. 215 No. 2, pp. 390-399.

Gyselinck, J., Sabariego, R. and Dular, P. (2006), “A nonlinear time-domain homogenization technique for laminated iron cores in three-dimensional finite-element models”, IEEE Transactions on Magnetics, Vol. 42 No. 4, pp. 763-766.

Gyselinck, J., Dular, P., Geuzaine, C. and Legros, W. (2002), “Harmonic-balance finite-element modeling of electromagnetic devices: a novel approach”, IEEE Transactions on Magnetics, Vol. 38 No. 2, pp. 521-524.

Handgruber, P., Stermecki, A., Bíró, O., Belahcen, A. and Dlala, E. (2013), “Three-Dimensional Eddy-Current analysis in steel laminations of electrical machines as a contribution for improved iron loss modeling”, IEEE Transactions on Industry Applications, Vol. 49 No. 5, pp. 2044-2052.

Hollaus, K. and Schöberl, J. (2015), “A higher order Multi-Scale FEM with a for 2D eddy current problems in laminated iron”, IEEE Trans. Magn, Vol. 51 No. 3.

Hollaus, K. and Schöberl, J. (2018), “Some 2-D multiscale Finite-Element formulations for the eddy current problem in iron laminates”, IEEE Transactions on Magnetics, Vol. 54 No. 4, pp. 1-16.

Hollaus, K., Schöberl, J. and Schöbinger, M. (2018), “Air gap and edge effect in the 2D/1D method with the magnetic vector potential a using MSFEM”, Proc. 18th IEEE CEFC, Hangzhou, 28-31 October.

Ledger, P. and Zaglmayr, S. (2010), “hp-Finite element simulation of three-dimensional eddy current problems on multiply connected domains”, Computer Methods in Applied Mechanics and Engineering, Vol. 199 Nos 49/52, pp. 3386-3401.

Rasilo, P., Dlala, E., Fonteyn, K., Pippuri, J., Belahcen, A. and Arkkio, A. (2011), “Model of laminated ferromagnetic cores for loss prediction in electrical machines”, IET Electric Power Applications, Vol. 5 No. 7, pp. 580-588.

Schöberl, J. and Zaglmayr, S. (2005), “High order Nédélec elements with local complete sequence properties”, Compel - the International Journal for Computation and Mathematics in Electrical and Electronic Engineering, Vol. 24 No. 2, pp. 374-384.

Schöbinger, M., Schöberl, J. and Hollaus, K. (2018), “An error estimator for multiscale FEM for the Eddy-Current problem in laminated materials”, IEEE Transactions on Magnetics, Vol. 54 No. 3, pp. 1-4.

Schöbinger, M., Schöberl, J. and Hollaus, K. (2019), “Multiscale FEM for the linear 2-D/1-D problem of eddy currents in thin iron sheets”, IEEE Transactions on Magnetics, Vol. 55 No. 1, pp. 1-12.

Stoll, R. (1974), The Analysis of Eddy Currents, Ser. Monographs in Electrical and Electronic Engineering, Oxford University Press, Oxford.

Yamada, S. and Bessho, K. (1988), “Harmonic field calculation by the combination of finite element analysis and harmonic balance method”, IEEE Transactions on Magnetics, Vol. 24 No. 6, pp. 2588-2590.

## Acknowledgements

This work was supported by the Austrian Science Fund (FWF) under Project P 27028.