Application of Mixed Finite Element and Natural Element Method in anti-periodic electromagnetic devices 1

Abstract

In this paper, the Mixed Finite Element and Natural Element (Mixed FEM-NEM) Method is applied to determine the flux distribution in magnetic devices using its reduced magnetic circuit models. In this solution, the NEM is applied exclusively in the gap of the devices where deformations of the mesh are expected in case of movement and the FEM is applied in the remainder magnetic circuit. In order to reduce the computational cost, the use of anti-periodic boundary conditions is proposed allowing the reduction of the entire domain model. To show the effectiveness of the proposed methodology, the problem is solved by using FEM and NEM in separate, and by the mixed FEM-NEM. The obtained results are compared using an error estimator and the computational cost.

Keywords

Magnetic flux reduced magnetic circuit model Mixed FEM-NEM Method anti-periodic boundary condition

1. Introduction

The Finite Element Method (FEM) is still the most popular numerical method for solving many problems in applied physics and engineering. However, difficulties may appear when solving problems involving the need of mesh reconstruction such as mobile parts or deformations [1]. In these cases, the Finite Element Method requires extensive computational and human effort [2]. In order to avoid these drawbacks, the Mesh-Free Methods have been proposed [3]. In the Meshless Methods a cloud of nodes without connectivity relations that covers the entire domain is used to solve the problem. Due to its characteristics, these methods are increasingly being used to solve electromagnetic problems where remeshing procedures are needed [4, 5]. It is well known that Meshless Methods provide high accuracy solutions but present some difficulties to handle boundary and interface conditions [6, 7]. So, to eliminate these issues, the Natural Element Method (NEM) was proposed [8].

The NEM is a Meshless Method based on the Voronoï diagram and the natural neighbor concept [9]. It provides highly accurate solution and the main interest in this method lies in its interpolation property, which allows enforcing directly the essential and periodic boundary conditions [9, 10]. In addition, due to its meshless characteristic, it does not present problems of mesh deformation in mobile interfaces, which makes it an interesting choice to solve problems involving moving parts [10, 5].

However, the computational cost could be high if the NEM is used alone to solve the entire domain of the device [9]. Therefore, in this kind of problem it is important to consider the symmetry in order to reduce the number of unknowns and the use of periodicity might be very useful to explore it. In this way, the anti-periodic boundary condition allows to work with part of device instead of the entire magnetic circuit. In case of the electrical machine model only $1/4$ of the domain needs be modeled [9, 4].

A strategy that combines the traditional FEM with the Meshless NEM is yet another way to reduce computational cost while maintaining versatility of both methods [10]. In the Mixed FEM-NEM approach, the NEM is used in the gap region, where the largest deformations of the mesh are expected, while the FEM is used to treat the remainder magnetic circuit of the device [10].

In this paper, the mixed FEM-NEM approach is applied to solve the magnetic circuits of two devices: a simplified model of a linear induction motor and a rotate electrical machine. Both use the anti-periodic conditions to reduce the evaluation domain. To show the effectiveness of the proposed methodology, the problems are also solved by using the FEM and the NEM alone. The obtained results are compared in terms of accuracy and computational cost.

2. Problem formulation

For the purpose of analysis, induction motors are simulated. The aim here is to obtain the magnetic flux distribution inside the induction motors through the solution of the magnetic vector potential in whole domain. So, the magnetic vector potential can be determined by the solution of Poison equation [9, 4]:

$\displaystyle\left\{{{\begin{array}[]{*{20}c}\nabla.\left({\frac{1}{\mu}\nabla A% }\right)=-{J}\\ {A}={A}_{d}\text{ on }\Gamma_{d}\\ \end{array}}}\right.$ (1)

Equation (1) is known as the strong formulation of the problem where ${A}$ is the magnetic potential vector, ${\mu}$ is the permeability, ${J}$ is the current density and ${A}_{d}$ is the magnetic potential imposed on the Dirichlet boundary ${\Gamma}_{d}$ . Note that it is possible to consider ${A}$ and ${J}$ as scalars since the current is imposed in only one direction. Based on Eq. (1), the following weak formulation can be written [9]:

$\displaystyle\int_{\Omega}\frac{1}{{\mu}}\nabla{A}\nabla wd\Omega=\int_{\Omega% }-Jwd\Omega,$ (2)

In Eq. (2) ${\Omega}$ is the problem domain surrounded by the surface ${\Gamma}={\Gamma}_{d}\cup{\Gamma}_{{ap}}$ where ${\Gamma}_{{ap}}$ is the anti-periodic boundary and ${w}$ is the scalar weight function.

The discrete equation is obtained by using the Galerkin Method. Many numerical techniques can be used in the evaluation of Eq. (2). However, this paper is focused on the mixed FEM-NEM.

Figure 1.

(a) The ${n}_{1}$ node, its Voronoï cell ${V}_{1}$ and its six natural neighbors. (b) NEM shape function computation where ${S}_{1}\left({x}\right)$ represents the subarea associated with the node ${n}_{1}$ .

3. Natural Element Method (NEM)

3.1 Voronoï diagram

The Natural Element Method uses the concept of natural neighbors which is based on the construction of Voronoï diagram on a cloud of nodes [9]. For instance, consider a set of nodes ${N}=\left\{{{n}_{1},{n}_{2},{n}_{3},\cdots,{n}_{M}}\right\}$ distributed in the whole domain ${\Omega}$ in 2D space as depicted by Fig. 1a. The Voronoï diagram is a subdivision of the domain into cells, where each cell ${V}_{i}$ is associated with a node ${n}_{i}$ , such that any point in ${V}_{i}$ is closer to ${n}_{i}$ than to any other node ${n}_{j}$ with ${j}\neq{i}$ . The regions ${V}_{i}$ are the Voronoï cells of ${n}_{i}$ and could be defined in mathematical terms as [9]:

$\displaystyle{V}_{i}=\left\{{{\rm{\bf x}}\in{\rm R}^{2}:{l}\left({{\rm{\bf x}}% ,{\rm{\bf x}}_{i}}\right)<l\left({{\rm{\bf x}},{\rm{\bf x}}_{j}}\right)\forall% {j}\neq{i}}\right\}.$ (3)

In Eq. (3) ${\rm{\bf x}}_{i}$ are the coordinates of the node ${n}_{i}$ and ${l}\left({{\rm{\bf x}},{\rm{\bf x}}_{j}}\right)$ is the distance between node ${n}_{j}$ and any point ${\rm{\bf x}}$ . Thus, ${V}_{i}$ is the region of the plane that contains the closest ${\rm{\bf x}}$ points to node ${n}_{i}$ in comparison to any other node in ${N}$ . The construction of each Voronoï cell can be obtained by the intersection of the segments joining the node ${n}_{i}$ to its neighbor nodes, and a straight line, normal to each one of these segments, traced at the central point of those segments [8]. This diagram subdivides the studied domain into a set of polygons which defines the natural neighbors of the node in its center. The Delaunay triangulation, which is the dual of the Voronoï diagram, is constructed by connecting the nodes which the Voronoï cells have common boundaries [9]. The interpolation of a point can be obtained by introducing a fictitious point in the Voronoï diagram and calculating the value of the corresponding shape function.

3.2 The shape function

Based on the Voronoï diagram, the natural element shape function can be calculated. In literature, several formulas are used to calculate this shape function [9]. Among the most used, are the Sibson functions which may be determined in analogy with classical FEM shape functions as the ratio of surfaces in the case of triangles [8]. The same principle is applied to Voronoï cells to achieve Sibson shape functions. At a point ${\rm{\bf x}}$ shown by Fig. 1b the Sibson shape function is given by Eq. (4) where each ${S}_{i}\left({\rm{\bf x}}\right)$ represents the subarea of Voronoï cell centered on ${\rm{\bf x}}$ and linked to the natural neighbor ${n}_{i}$ and ${S}\left({\rm{\bf x}}\right)$ is the total area of Voronoï cell related to that point [8].

$\displaystyle\phi_{i}\left({\rm{\bf x}}\right)=\frac{{S}_{i}\left({\rm{\bf x}}% \right)}{{S}\left({\rm{\bf x}}\right)},$ (4)

where ${i}$ ranges from 1 to ${NN}$ (number of natural neighbors of point ${\rm{\bf x}}$ ), and:

$\displaystyle{S}\left({\rm{\bf x}}\right)=\mathop{\sum}\limits_{{i}=1}^{{NN}}{% S}_{i}\left({\rm{\bf x}}\right).$ (5)

More details about the construction of NEM shape functions can be found in [8]. Equation (4) verifies the same properties of FEM shape functions. The Kronecker delta, interpolation, partition of unity and linear completeness properties are verified [5, 8, 9]. Thus, the magnetic potential vector can be written as follows:

$\displaystyle{A}^{h}\left({\rm{\bf x}}\right)=\mathop{\sum}\limits_{{i}=1}^{{% NN}}{A}_{i}\phi_{i}\left({\rm{\bf x}}\right).$ (6)

It is important to note that, due to the existence of subdomains with different materials, as well as other Meshless Methods, it is necessary to use the visibility criterion for the treatment of material discontinuity [6, 7]. This criterion modifies the original form of Voronoï diagram on the boundary of different materials that becomes a constrained Voronoï diagram. Moreover, if the visibility criterion is introduced in the NEM, the natural neighbors become constrained natural neighbors. The set of natural neighbors will be restricted by applying the visibility criterion. This implies that the method used is the so-called Constrained Natural Element Method (C-NEM) [6]. Thus, the approximation function for C-NEM takes account only the natural neighbors nodes visible from point ${x}$ .

In the C-NEM, the support of the constrained natural neighbor shape function $\phi_{i}^{C}$ is similar to $\phi_{i}$ for NEM interpolation, that consists of the union of the circumcircles which contain the node ${n}_{i}$ . However, these circumcircles can extend across the boundaries, but cannot contain non-visible neighbors [6]. Figure 2 shows the support of $\phi_{i}^{C}$ referring to a node ${n}_{i}$ for a regular grid of nodes.

Figure 2.

Support for shape function in C-NEM. (a) Grid of nodes – empty nodes $\phi_{i}^{C}=0$ and nodes fulfilled $\phi_{i}^{C}=1$ , (b) Three-dimensional view [6].

The C-NEM maintains the properties of the NEM as interpolation, partition of unity and linearity of the shape functions in both convex and non-convex domains [6, 2]. These properties facilitate the imposition of the essential boundary conditions as well as the periodic boundary conditions [9].

As the NEM and the FEM shape functions share the same proprieties, the unknowns could be interpolate conveniently by one of them. So, in the Mixed FEM-NEM approach, the NEM is used in the gap region, where the largest deformations of the mesh are expected, while the FEM is used to treat the remainder magnetic circuit of the device.

4. Results

Consider the structure showed at Fig. 3 to validate the proposed approach. This structure can be considered a simplified model of a linear induction motor. The circles are aluminum conductors carrying a current density in the alternating directions of 1 MA/m ${}^{2}$ and involved by an iron material with $\mu_{r}=$ 7000. The dots in the centers of the circles represent current going out of the sheet plane while the cross represents current entering the sheet plane.

Figure 3.

Anti-periodic structure.

Figure 4.

The resulting induction flux distribution for an anti-periodic structure using Mixed FEM-NEM.

The anti-periodic boundary conditions are set at the left and right boundaries while zero Dirichlet ones are set at the top and bottom boundaries [9]. The potential is evaluated on a line orthogonal to the gap showed by a dashed line in Fig. 3.

Figure 4 shows the resulting induction flux distribution for a 3,658 nodes discretization solved by the Mixed FEM-NEM method.

The problem is also solved using only FEM and only C-NEM for the same number of nodes. Figure 5 shows the potential distribution along the orthogonal direction in the gap (dashed line in Fig. 3) obtained by the Mixed FEM-NEM method and the individual ones.

The scalar error ${\varepsilon}$ for each method is obtained by the error estimator proposed by [10, 11]. This error is based on the verification of the constitutive relationship that links the magnetic flux density ${\rm{\bf B}}$ and the magnetic field ${\rm{\bf H}}$ , given by [10]:

$\displaystyle{\rm{\bf B}}=\mu{\rm{\bf H}}$ (7)

where ${\mu}$ is the magnetic permeability. So, the error estimator can be compute by [11]:

$\displaystyle{\varepsilon}^{2}=\frac{\|{\rm{\bf B}}_{A}-{\mu}\left({{\rm{\bf H% }}_{A}+{\rm{\bf H}}_{S}}\right)\|_{\Omega}^{2}}{\|{\rm{\bf B}}_{A}+{\mu}\left(% {{\rm{\bf H}}_{A}+{\rm{\bf H}}_{S}}\right)\|_{\Omega}^{2}}$ (8)

where ${\rm{\bf X}}_{\Omega}^{2}=\int_{\Omega}{XX}/{\mu}$ , ${\rm{\bf B}}_{A}$ is the magnetic flux density admissible, ${\rm{\bf H}}_{A}$ is the magnetic field admissible and the ${\rm{\bf H}}_{S}$ is the coercive magnetic field, while ${\rm{\bf B}}_{A}$ and ${\rm{\bf H}}_{A}$ are obtained through the potential vector ${\rm{\bf A}}$ of the Eq. (1), the ${\rm{\bf H}}_{S}$ is obtained by the solution of the scalar potential $\nabla{w}$ [11]:

$\displaystyle-\nabla\left({{\mu}\nabla{w}}\right)=-\nabla\left({\mu{\rm{\bf H}% }_{S}}\right)\text{in }\Omega$ (9)

Figure 5.

Potential distribution along the orthogonal direction in the gap (dashed line in Fig. 3).

The error is evaluated along the orthogonal direction in the gap and showed in Fig. 6. Note that the Mixed Method is more accurate than FEM. However, the appropriate choice of the method should take into account the accuracy and the computational cost. The biggest drawback of NEM lies with the computational cost. Using as reference the Natural Element Method, the use of the Mixed Method resulted in 86% reduction in computational time with loss of precision of only 25%. Therefore the results show that the use of the Mixed Method is qualitatively and computationally advantageous.

Figure 6.

Error estimator along the parallel line to anti-periodic boundary.

Figure 7.

Geometry of a reduced electric induction machine. The dimensions are given in meters and $\Omega$ represents the entire domain to be studied including the stator, rotor and air gap regions [12].

Figure 8.

Magnetic flux distribution using mixed FEM-NEM.

Figure 9.

Components of the magnetic potential vector obtained through the Mixed FEM-NEM and evaluated at the lower anti-periodic boundary.

Figure 10.

Error estimator along the anti-periodic boundary.

The same procedure is applied to solve a rotate induction motor. The C-NEM is used only in the gap of the machine while the remaining domain is solved by FEM. The motor to be analyzed is a 2HP three-phase induction type, four-pole, current-fed and squirrel cage-type rotor. There are a total of 36 slots on the stator and 28 slots on the rotor. The material of the rotor bars is aluminum, the stator coils are copper wire and both are encased by iron silicon which is a linearly responsive material [4]. So, using anti-periodic boundary conditions and due to symmetry considerations, only $1/4$ of it needs to be modeled [4]. The Fig. 7 shows the reduced 2D computational domain model for this motor.

For simulation, the effects of the induced current are ignored and it is considered that each slot of the stator is covered by 44 ampere-turns. Assuming the flow is concentrated only inside the machine, $A=$ 0 is imposed on ${\Gamma}_{d}$ (the Dirichlet boundary). The data about the geometry and discretization of the machine are obtained using the software FEMM that uses the Finite Element Method for the calculation of electromagnetic fields [12].

The choice of anti-periodic boundary conditions for the induction motor is due to the way the stator coils are distributed. These coils have the following phase sequence ( ${\cal B};-{\cal C};{\cal A}$ ), ( $-{\cal B};{\cal C};-{\cal A}$ ), ( ${\cal B};-{\cal C};{\cal A}$ ) e ( $-{\cal B};{\cal C};-{\cal A}$ ), such that [4]:

$\displaystyle\left\{{{\begin{array}[]{*{20}l}{{\cal A}\left({t}\right)={\cal A% }_{\text{max}}\cos\left({\omega t}\right)}\\ {{\cal B}\left({t}\right)={\cal B}_{\text{max}}\cos\left({\omega t+120^{\circ}% }\right)}\\ {{\cal C}\left({t}\right)={\cal C}_{\text{max}}\cos\left({\omega t-120^{\circ}% }\right)}\\ \end{array}}}\right.$ (10)

where ${t}$ represents the instant of observation given in the second and, $\omega$ is the angular frequency given in radians/second.

Figure 8 shows the obtained magnetic flux distribution evaluated using mixed FEM-NEM with 12,837 nodes and considering the relative permeability $\mu_{r}=$ 7,000 in the iron domain.

By evaluating the magnitude of the magnetic vector potential at the lower anti-periodic border, it is possible to obtain the absolute value and the real and imaginary parts of vector ${A}$ . These data are also obtained by Mixed FEM-NEM Method as shown in Fig. 9.

The same problem is also solved using only FEM and only C-NEM and the accuracy of the results are compared with Mixed Method.

The error is evaluated at the lower periodic boundary and presented in Fig. 10. Note that the Mixed Method is more accurate than FEM.

The time spent by the Finite Element Method algorithm for solving the problem is taken as the reference. Based on it, the C-NEM solved the problem spending 8.33 times more while the Mixed FEM-NEM Method spent 4.39 times more than the reference. An analysis of the computational cost allows verifying that the use of the Mixed Method reduces up to 50% the time needed to obtain the result using the C-NEM in whole domain.

5. Conclusion

The use of C-NEM guarantees greater accuracy when compared to FEM. The Mixed Method manages to combine part of the precision of the C-NEM and the speed of the FEM.

The computational cost can be reduced by using the periodic boundary conditions which allow the reduction of evaluation domain.

The Mixed Method works with the application of the Meshless Method only in the region of the gap of the device where the mesh could suffer deformations.

The analysis of results has demonstrated the applicability and the accuracy of the Mixed FEM-NEM Method to solve efficiently anti-periodic devices.

References

Liu

et al., A direct coupling method of meshless local petrov-galerkin (MLPG) and finite element method (FEM), International Journal of Applied Electromagnetics and Mechanics51(1) (2016), 51.

Gonçalves

B.M.F.

et al., Comparison of the Natural Element and Finite Element Interpolation Functions Behavior, Int. J. Numer. Model. e2230 (2017), 1–7.

Liu

G.R.

, Mesh Free Methods – Moving Beyond the Finite Element Method, CRC Press, 2003.

Coppoli

E.H.R.

Silva

R.S.

and Mesquita

R.C.

, Periodic Boundary Conditions in Element Free Galerkin Method, COMPEL: Int. J. Comput. Math. Electr. Electron. Eng.28(4) (2009), 922–934.

Illoul

Yvonet

Chinesta

and Clénet

, Application of the Natural-Element Method to Model Moving Electromagnetic Devices, IEEE Trans. Magn.42(4) (2006), 727–730.

Yvonnet

Ryckelynck

Lorong

and Chinesta

, A new extension of the natural element method for non-convex and discontinuous domains: the constrained natural element method (C-NEM), Internat. J. Numer. Methods Engrg.60 (2004), 1451–1474.

Gonçalves

B.M.F.

et al., The Constrained Natural Element Method Applied in Electromagnetic Heterogeneous Models, Sensors & Transducers198(3) (March 2016), 55–60.

Sukumar

, The Natural Element Method in Solid Mechanics, Ph.D. dissertation, Northwestern Univ., USA, 1998.

Gonçalves

B.M.F.

et al., Periodic Boundary Conditions in Natural Element Method, IEEE Trans. Magn.52(3) (2015), 1–4.

10.

Illoul

Le Menach

Clenet

and Chinesta

, A mixed finite element/meshless natural element method for simulating rotative electromagnetic machines, Eur. Phys. J. Appl. Phys.43 (2008), 197–208.

11.

Bruyère

Illoul

Chinesta

and Clénet

, Comparison Between NEM and FEM in 2-D Magnetostatics Using an Error Estimator, IEEE Trans. Magn.46 (2008), 1342–1345.

12.

Meeker

, Finite Element Method Magnetics. Available in: http://www.femm.info/wiki/HomePage. Access: 09/23/ 2015.