Nodal and Edge Finite Element Analysis of Eddy Current by a76m823ik


									                                                                                                                        Miklós KUCZMANN
                                                                                                                     Széchenyi István University

                                                 Nodal and Edge Finite Element Analysis of
                                                             Eddy Current Field Problems
Abstract. Eddy current field problems can be solved by different potential formulations based on the “quasi-static” Maxwell’s equations. The
potential formulations are obtained using a vector potential and a scalar potential, the most widely used techniques are the A,V-A, the T,Φ-Φ
formulations and their combinations. Vector potentials can be approximated by nodal or edge finite elements, scalar potentials are approximated by
scalar elements. The paper presents and compares these formulations through TEAM Problem No. 7 containing a multiply connected region.

Streszczenie. Problemy związane w prądami wirowymi mogą być rozwiązywane przy wykorzystaniu różnych sformułowań i różnych potencjałów,
bazujących na quasistatycznych równaniach Maxwella. Wykorzystywany jet w tych sformułowaniach potencjał wektorowy i potencjał skalarny,
najczęściej są to sformułowania A,V-A, oraz T,Φ-Φ, a także ich kombinacje. Potencjał wektorowy jest aproksymowany przez elementy węzłowe i
krawędziowe, a potencjał skalarny przez elementy skalarne. Artykuł przedstawia różne sformułowania i prowadzi ich analizę komparatystyczną na
bazie problemu No. 7 z katalogu TEAM Workshop. (Węzłowe i krawędziowe elementy w analizie problemów wiroprądowych)

Keywords: eddy current field, potential formulations, finite element method
Słowa kluczowe: prądy wirowe, sformułowania potencjałów, metoda elementów skończonych

Introduction                                                                    There are well-known potential formulations [1-5], here a
     A linear eddy current field problem contains a source of              short summary is presented.
the electromagnetic field, which can be a current driven or a                   The static magnetic field in Ω n can be described by the
voltage driven coil, and it generates a time varying magnetic
                                                                           reduced magnetic scalar potential Φ or by the magnetic
field in the vicinity of the coil. This field induces eddy
                                                                           vector potential A. The combination of these formulations is
currents inside a conducting material (e.g. inside an
                                                                           also valid.
aluminum plate), which currents generate a magnetic field,
                                                                                The source current density can be represented by the
and this field modifies the magnetic field supplied by the
                                                                           so-called impressed current vector potential T0 by
     In this case, not only the magnetic field, but the eddy                J 0 = ∇ × T0 , since ∇ ⋅ J 0 = 0 , and the magnetic field
current field must be simulated.                                           intensity can be represented by
     Eddy current field problems can be solved by using
vector potentials combined with scalar potentials in the                   (2) H = T0 − ∇Φ
frame of the Finite Element Method (FEM).
                                                                           from the first equation in (1). Using the magnetic Gauss law
Governing equations                                                        results in the following partial differential equation:
   The following linear static and “quasi-static” Maxwell’s
equations can be used to simulate linear eddy current field                (3) − ∇ ⋅ [μ∇Φ ] = − μ∇ ⋅ T0     in   Ωn .
                                       ⎧∇ × H = J,                         Here, the impressed current vector potential T0 can be
                                       ⎪                                   calculated by the Biot-Savart law.
    ⎧∇ × H = J 0 ,                     ⎪∇ × E = − B,                           The tangential component of magnetic field intensity can
    ⎪                                  ⎪
(1) ⎨∇ ⋅ B = 0 ,   in Ω n ,        and ⎨∇ ⋅ B = 0 ,  in Ωc .               be prescribed by the Dirichlet boundary condition Φ = Φ 0
    ⎪ B = μH,                          ⎪ B = μH,
    ⎩                                  ⎪                                   on ΓH n , where Φ 0 is constant [2,5], and a Neumann
                                       ⎪J = σE
                                       ⎩                                   boundary condition μ [T0 − ∇Φ ]⋅ n = 0 on ΓB prescribes
                                                                           B⋅n = 0 .
Here H, J0, J, E, B, σ and μ are the magnetic field
                                                                             Magnetic vector potential is introduced by the magnetic
intensity, the source current density, the eddy current                    Gauss law, i.e.
density, the electric field intensity, the magnetic flux density,
the conductivity, and the permeability, respectively.                      (4) B = ∇ × A ,
    The domain Ω n is the air region where source current
                                                                           and from the first equation in (1) the partial differential
is placed and μ = μ 0 , furthermore Ωc is the region where
eddy currents are present. The two regions are coupled by
interface conditions, i.e. the tangential component of the                 (5) ∇ × [ ∇ × A] = J 0
                                                                                   ν                  in   Ωn
magnetic field intensity and the normal component of the
magnetic flux density are continuous on Γnc . On the
                                                                           can be obtained. Here ν = 1 / μ is called the magnetic
boundary      ∂Ω n = ΓH n ∪ ΓB , and on            ∂Ω c = ΓH c ∪ ΓE        reluctivity.
boundary conditions are prescribed as H × n = 0 on ΓH n ,                      The tangential component of magnetic field intensity is
                                                                           prescribed     by   a    Neumann       boundary     condition
and    B ⋅ n = 0 on ΓB , and H × n = 0 on ΓH c , and                        [ν∇ × A]× n = 0 on ΓH n , and a Dirichlet boundary condition
E × n = 0 on ΓE .                                                          n × A = 0 on ΓB prescribes B ⋅ n = 0 .

194                                                          PRZEGLĄD ELEKTROTECHNICZNY, ISSN 0033-2097, R. 84 NR 12/2008
    Two potentials can be introduced in the eddy current              and
region, either a current vector potential T or a magnetic
vector potential A to simulate either the eddy current field or              ∇ × [ ∇ × A] − ∇[ ∇ ⋅ A] + σA + σ∇V = 0, and
                                                                                 ν           ν
                                                                                 [          ]
the magnetic field. The first one must be coupled with a              (11)
reduced magnetic scalar potential Φ, while the magnetic
                                                                             − ∇ ⋅ σA + σ∇V = 0.
vector potential must be coupled with an electric scalar              Some extra boundary conditions have to be satisfied in
potential V.
                                                                      these situations, too, A ⋅ n = 0 on ΓH n , and ν∇ ⋅ A = 0 on
    Eddy currents are solenoidal, i.e. ∇ ⋅ J = 0 , from which
the current vector potential can be introduced as                     ΓB in case of (9), ρ∇ ⋅ T = 0 on ΓH c , and T ⋅ n = 0 on ΓE
 J = ∇ × T . The magnetic field intensity can be represented
                                                                      in case of (10), A ⋅ n = 0 on ΓH c , and ν∇ ⋅ Α = 0 on ΓE in
by two unknown potentials as
                                                                      case of (11).
(6) H = T0 + T − ∇Φ                                                      Coulomb gauge is satisfied automatically when applying
                                                                      edge element based FEM to represent vector potentials,
from Ampere’s law.                                                    however the right hand side of the partial differential
   The following partial differential equations can be written        equation (5) must be in the same function space, i.e. T0
from Maxwell’s equations:                                             must be used represented by edge elements [4],
      ∇ × [ρ∇ × T ] + μT − μ∇Φ = − μT0 , and                          (12) ∇ × [ ∇ × A] = ∇ × T0 .
      ∇ ⋅ [μT − μ∇Φ ] = −∇ ⋅ [μT0 ] in Ω c .
                                                                      The partial differential equations of the vector element
Here ρ = 1 / σ         is the resistance of the material. The         representation of T,Φ-formulation, and A,V-formulation are
                                                                      the same as in (4) and in (5). These formulations are called
Dirichlet boundary conditions T × n = 0 , and Φ = Φ 0 on              ungauged version of the potential formulations, however
ΓH c ,        and      the       Neumann      boundary   conditions   (9), (10) and (11) result in the gauged version, when vector
                                                                      potentials are approximated by nodal elements.
 [ρ∇ × T ]× n = 0 , and μ [T0 + T − ∇Φ ]⋅ n = 0 on ΓE
prescribe the boundary conditions of the eddy current field           Numerical simulation by FEM
problem.                                                                 Benchmark problem No. 7 of the TEAM Workshop
    The following partial differential equations are valid in         consists of an asymmetrical conductor made of aluminum
the case of eddy current field problems when applying                 ( σ = 3.526·10 S/m ) with a hole (Fig. 1) [6].
magnetic vector potential:                                                The source of magnetic field is the sinusoidal current
                                                                      (maximum ampere-turn is 2742AT, and f = 50Hz,
      ∇ × [ ∇ × A] + σA + σ∇V = 0 , and
          ν                                                           f = 200Hz) flowing in the coil placed above the plate. The
               [             ]
       − ∇ ⋅ σA + σ∇V = 0 in Ω c .                                    problem is a linear steady state eddy current field problem
                                                                      with multiply connected eddy current region. It has been
The Neumann boundary conditions [ ∇ × A]× n = 0 and
                                ν                                     solved in the frequency domain.

      [            ]
− σ A − ∇V ⋅ n = 0 on ΓH c , and the Dirichlet boundary
conditions n × A = 0 , and V = V0 on ΓE prescribe the
boundary conditions of the eddy current field problem.
    Eddy current field problem must be coupled to static
magnetic field problem when a conducting material is
placed into a non-conducting one. The coupling of the
above potentials results in the following formulations:
T,Φ−Φ, A,V−A, T,Φ−A, T,Φ−A−Φ, A,V−Φ, A,V−A−Φ. In this
case, the interface conditions have to be implemented by
the introduced potentials [1-5].
    Scalar potentials can be approximated by nodal shape
functions and the unknowns are associated to the nodes of
the FEM mesh. Nodal shape functions as well as vector
shape functions can be used to represent vector potentials.
In the first case, three unknowns are associated to the
nodes, and in the second case, unknowns are associated to
the edges.
    It is well known from the literature that these
formulations are very sensitive to gauging when
representing vector potentials by nodal shape functions. It
can be worked out by the implicit enforcement of Coulomb
gauge on the left hand side of partial differential equations
in (5), (7) and (8), i.e. the following modifications must be
applied [2]:

(9) ∇ × [ ∇ × A] − ∇[ ∇ ⋅ A] = J 0 ,
        ν           ν
                                          [      ]
                                                                      Fig.1. Geometry and FEM mesh of TEAM Problem No. 7
          ∇ × [ρ∇ × T ] − ∇[ρ∇ ⋅ T ] + μ T − ∇Φ = − μT0 , and
(10)                                                                       The reference solution is obtained by the A,V−A-
          ∇ ⋅ [μT − μ∇Φ ] = − μ∇ ⋅ T0 ,                               formulation. The eddy currents have a path around the hole
                                                                      filled with air, as it can be seen in Fig. 2. The use of the

PRZEGLĄD ELEKTROTECHNICZNY, ISSN 0033-2097, R. 84 NR 12/2008                                                                  195
reduced magnetic scalar potential Φ in the eddy current free          well the A,V-A-Φ is faster than the nodal version of these
region results in a false solution, as it is illustrated in Fig. 3.   formulations. The nodal version of the scalar potential
This is the situation when applying the A,V−Φ-formulation or          based formulation is faster than the edge based ones.
the T,Φ−Φ-formulation. It is noted that the aim of applying
reduced scalar potential in the air region is to decrease the
degree of freedom however the solution is incorrect
because the aluminum plate is a multiply connected region.

Fig.2. Correct distribution of eddy currents in the aluminum

Fig.3. False distribution of eddy currents in the aluminum

     Applying some tricks, this problem can be eliminated.
The problem of the A,V−Φ-formulation can be solved by
using the so-called A,V−A−Φ-formulation, when the
magnetic vector potential A is introduced in the hole, too,
i.e. the multiply connected region can be transformed into a
simple connected one. The incorrect solution of the T,Φ−Φ-
formulation can be corrected by employing the T,Φ−A-
formulation or the more economic T,Φ−A−Φ-formulation.
The magnetic vector potential A is introduced in the whole
air region, or only in the hole, respectively. The other
possibility is filling the hole with a conducting material with
very low conductivity, i.e. the T,Φ−Φ-formulation can be
used without any modification.
     Fig. 4 and Fig. 5 show some comparisons between
measured [6] and simulated results. Here, the z component
of the magnetic flux density and the y component of the
eddy current density have been shown along the line
x=0,…,288mm, y=72mm, and y=144mm, z=34mm, and
x=0,…,288mm,           y=72mm,      z=0mm,    and    z=19mm,
respectively. The solutions of the different potential
formulations are practically identical. Degree of freedom
(DOF), number of iteration of the linear system solver, and           Fig. 4. Bz component, f = 50Hz
the CPU time (AMD Athlon 64 X2 Dual Core Processor
4600+ 2.41GHz, 4GByte RAM) can also be studied in Table
1-Table 4. The edge element based A,V-A formulation as

196                                                          PRZEGLĄD ELEKTROTECHNICZNY, ISSN 0033-2097, R. 84 NR 12/2008
                                                                       Table 2. Data of computations
                                                                         DOF          Iteration Time
                                                                      55554/nodal        980    1450
                                                                      55554/nodal       9995    9919
                                                                      58377/vector       185     259
                                                                      58377/vector       180     255

                                                                       Table 3. Data of computations
                                                                      DOF             Iteration Time
                                                                      53204/nodal     50        319
                                                                      53204/nodal     60        345
                                                                      55949/vector 216          278
                                                                      55949/vector 277          355

                                                                       Table 4. Data of computations
                                                                      DOF             Iteration Time
                                                                      53422/nodal     47        340
                                                                      53422/nodal     44        325
                                                                      55999/vector 2353         2409
                                                                      55999/vector 1092         1171

                                                     The paper is a summary of some potential formulations and
                                                     their FEM representation. Detailed presentation is shown in
                                                     the book [5]. The advantages and disadvantages of the
                                                     presented methods have been shown by a problem
                                                     containing a multiply connected region. The aim of further
                                                     research is to solve nonlinear eddy current field problems.

                                                     Acknowledgments: This paper was supported by the János
                                                     Bolyai Research. Scholarship of the Hungarian Academy of
                                                     Sciences (BO/00064/06), by the Hungarian Scientific
                                                     Research Fund (OTKA PD 73242), and by the Széchenyi
                                                     István University (15-3210-02).

                                                     [1] Jackson J. D., Classical electrodynamics, J. Wiley, New York,
                                                     [2] Bíró O., Richter K. R., CAD in electromagnetism, Adv. in Elect.
                                                         and Elect. Physics, Vol. 82, 1991.
                                                     [3] Jin J., The finite element method in electromagnetics, John
                                                         Wiley and Sons, New York, 2002.
                                                     [4] Bíró O., Edge element formulations of eddy current problems.
                                                         Comput. Meth. Appl. Mech. Engrg., Vol. 169, 1999, pp. 391–
                                                     [5] Kuczmann M., Iványi A., The finite element method in
                                                         magnetism, Akadémiai Kiadó, Budapest, to be published in
                                                     [6] Fujiwara K., Nakata T., Results for benchmark problem 7.
                                                         COMPEL, Vol. 11, pp. 335–344, 1990.

                                                     Author: Dr. Miklós Kuczmann, PhD, Széchenyi István University,
                                                     Department of Telecommunications, Laboratory of Electromagnetic
                                                     Fields, Egyetem tér 1, H-9026 Győr, Hungary, E-mail:

Fig. 5. Jy component, f = 200Hz

                 Table 1. Data of computations
                DOF             Iteration Time
                122020/nodal    532        1924
                122020/nodal    857        2454
                132344/vector 132          616
                132344/vector 132          616

PRZEGLĄD ELEKTROTECHNICZNY, ISSN 0033-2097, R. 84 NR 12/2008                                                       197

To top