Mathematical Models and Methods in Modern Science

Finite Element Method for calculation of magnetic field produced from a helical turn in linear and nonlinear medium IULIA CATA, DUMITRU TOADER Department of Physical Foundation of Engineering "Politehnica" University of Timisoara Romania [email protected] [email protected] Abstract: - The paper presents a numerical model using finite element method to calculate the magnetic field produced by current through a helical turn, whose thickness is comparable to the winding radius. The magnetic field was determined for helical turn considering homogeneous-linear and inhomogeneous-nonlinear medium, using the finite element software package Vector Field Opera. It analyzes the numerical model, in order to optimize it, in relation to the number of finite elements spatial, the extension of domain to calculate the field, current carrying conductor length and material characteristics of medium of the magnetic field existence.

Key-Words: - magnetic field, helical turn, numerical models, finite element method.

1 Introduction

imposed current density; σ ⋅

Calculation of the magnetic field coils is usually, considering filamentary and circular turns [1], [2], [3], [4]. In reality turns are filamentary having thickness and wound helically. For this reason it is necessary design a model for calculating the magnetic field to eliminate these shortcomings. Since the analytical models presented in the literature considers filamentary conductors is appropriate to use numerical models. Using of ferromagnetic cores (nonlinear medium) in the construction of coils, requires analysis of the magnetic field in nonlinear medium. The paper presents a numerical model to calculate the magnetic field in both linear and nonlinear medium. Numerical model has been implementing as in software package Vector Fields Opera.

density. The functional is,  B 1   A ∂ A  F = ∫ ∫ (B − BR )d B − ∫ (J a − σ ⋅ )d Adv − V 0 ∂t   0  µ(B)  (2.2) − ∫ A × H ⋅ n dΣ − ∫ ( J S ⋅ A)ds Σ

)

Σ

S

B - the magnetic flux density, B R = µ 0 ⋅ M p - the magnetic flux density for

permanent magnets, µ0 – permeability of free space,

M p - the permanent magnetization , µ (B) - is the permeability tensor, Σ - is the boundary domain, dΣ - is the surface element of the boundary, n Σ - is the normal unit vector of domain boundary, A × H - the density of magnetic energy transferred through boundary Σ, J S - the current density of surface S contained in domain V, A - magnetic vector potential. In the case of linear medium, homogeneous and without permanent magnetization (missing permanent magnets) relation (2.2) becomes,  B 1     A ∂A F = ∫  ∫ Bd B  − ∫ ( J a − σ ⋅ )d A dv − V ∂t  0 µ  0  

Variational model of the magnetic field in cvasistationar regime for conductive environments fixed, nonlinear, with magnetic anisotropy, inhomogeneous and with permanent magnets admits extremize in relation with magnetic vector potential. In this case the current density is [5, 6], ∂A ∂A (2.1) = J a −σ ⋅ J = σ ⋅ E = −σ ⋅ gradV − σ ⋅ ∂t ∂t where: E - electric field strength; σ – conductivity;

(

)

− ∫ A × H ⋅ n dΣ − ∫ ( J S ⋅ A)ds Σ

Σ

(2.3)

S

Numerical solution of mathematical model variation for electromagnetic field involves

A - magnetic vector potential; J a -

ISBN: 978-1-61804-055-8

(

Where

2 The variational approach to finite elements of the magnetic field

t – time;

∂A - induced current ∂t

100

Mathematical Models and Methods in Modern Science

Partitioning the domain into finite elements impose ensure the continuity conditions of boundary.

extremize of functional (2.3) in conditions boundary specified [7, 8, 9, 10]. Such an approximation can be made using Rayleigh-Ritz-Galerkin's procedure, by which the extremize of functional transforms in determination of extreme of a function which has more variable, becoming a problem of classical mathematical analysis. Considering in the set of definition of functional F[Ψ], a system of n linear independent functions

3 Three-dimensional finite element The three-dimensional finite element, at which interpolation polynomial is linear, has 4 nodes (p=4), so i = 1,4 . In the cartesian coordinate system, energy functional (2.3), for linear environments without permanent magnetization, domain does not contain current density surface and Dirichlet conditions on the boundary becomes, F = Fx + Fy + Fz (2.7)

{φi}, i = 1, n called coordinate functions. Function Ψ which extremize the functional F[Ψ] is approximated by a linear combination of coordinate functions obtaining, n

Ψ = ∑ wi ⋅ ϕ i

(2.4)

i =1

where

where wi is the variational parameters. Introducing Eq. (2.4) in functional F[Ψ] this become function of coefficients wi. So its extremize leads to the next system (2.5),

∂F ( wi ) =0 ∂wi

 1  ∂A ∂A  2  ∂A  y  − J ax ⋅ Ax + σ x Ax  dxdydz  z − Fx = ∫    V 2µ ∂z  ∂t  ∂y     2  1  ∂A ∂A    ∂A x Fy = ∫  − z  − J ay ⋅ Ay + σ y Ay  dxdydz   V 2µ ∂x  ∂t   ∂z   2  1  ∂A ∂A    ∂A  y − x  − J az ⋅ Az + σ z Az  dxdydz Fz = ∫  V 2µ ∂y  ∂t  ∂x   

(2.5)

Accuracy of process depends on way to choose coordinate functions {φi} and their number. In finite element method (FEM) the functions of coordinate are approximated on portions through interpolation polynomials. The domain where the existence of electromagnetic field to partition into disjunctive sub domains and for each sub domain (finite element) specifies the functions of coordinates (interpolation polynomials). In these conditions the energy functional is,

Expression of potential vector A for finite element “e” accepting linear interpolation polynomial becomes,

A e = Aex ( x, y , z , t ) ⋅ i + Aey ( x, y , z , t ) ⋅ j + 4

+ Aez ( x, y, z , t ) ⋅ k = ∑ N ei ( x, y, z ) ⋅ A ei ( x, y, z, t ) . i =1

m

F [Ψ ] = ∑ Fe (Ψe )

So Fx, Fy, Fz for finite element “e” results the expressions,

(2.6)

e =1

where "m" - the number of finite element of domain. 2  1  ∂ 4 4  ∂ 4 Fx = ∫   (∑Nei(x, y, z) ⋅ Aeiz(x, y, z, t)) − (∑Nei(x, y, z) ⋅ Aeiy(x, y, z, t)) − Jax ⋅ (∑Nei(x, y, z) ⋅ Aeix(x, y, z, t)) + Ve 2µ ∂z i=1 i =1   e  ∂y i=1 4 ∂ 4  +σ ( ∑ N ei ( x , y , z ) ⋅ A eix ( x , y , z , t ) ) ⋅ ( ∑ N ei ( x , y , z ) ⋅ A eix ( x , y , z , t ) )   dxdydz (2.8) ∂ t i =1 i =1 

2 4  1  ∂ 4 ∂ 4  Fy = ∫   (∑Nei(x, y, z)⋅ Aeix(x, y, z, t)) − (∑Nei(x, y, z)⋅ Aeiz(x, y, z, t)) − Jay ⋅ (∑Nei(x, y, z) ⋅ Aeiy(x, y, z, t)) + Ve 2µ ∂x i=1 i=1   e  ∂z i=1 4 ∂ 4  +σ ( ∑ N ei ( x , y , z ) ⋅ A eiy ( x , y , z , t ) ) ⋅ ( ∑ N ei ( x , y , z ) ⋅ A eiy ( x , y , z , t ) )   dxdydz (2.9) ∂ t i =1 i =1 

ISBN: 978-1-61804-055-8

101

Mathematical Models and Methods in Modern Science

2  1  ∂ 4 4  ∂ 4 Fx = ∫   (∑Nei(x, y, z) ⋅ Aeiy(x, y, z,t)) − (∑Nei(x, y, z) ⋅ Aeix(x, y, z,t)) − Jaz ⋅ (∑Nei(x, y, z) ⋅ Aeiz(x, y, z,t)) + Ve 2µ ∂y i=1 i=1   e  ∂x i=1 4 ∂ 4  +σ ( ∑ N ei ( x , y , z ) ⋅ A eiz ( x , y , z , t ) ) ⋅ ( ∑ N ei ( x , y , z ) ⋅ A eiz ( x , y , z , t ) )   dxdydz (2.10) ∂ t i =1 i =1 

Similarly, the Post-Processor provides facilities necessary for the calculation of electromagnetic fields, with the possibility of displaying graphic and numerical form. Establishing the optimal numerical model using finite element method implemented in Vector Fields Opera program is to calculate the magnetic field created by current through the helical turn wrapped around a cylinder (Fig.1) of radius a. Helical turn has circular section of radius a1 and length h. Numerical calculation was made considering a = a1 = 1.35 mm, h = 81mm. Helical turns modeling was done using PreProcessor program, being made from a predetermined number of segments so that this wrap exactly over the cylinder. Helical turn is then used in the two programs to create the whole model. Domain where is calculated the magnetic field is cylindrical with length h1 + h + h2 and radius b (Fig.2), where h1 = h2 = 10h.

where Nei is trial functions. Considering the equations (2.8), (2.9), (2.10) results a system of linear equations of the form,





[M ] ⋅ [Ai ] + [C ] ⋅  ∂ Ai  + [F i ] = 0  ∂t 

(2.11)

where [M] – square matrix of linear system, [C] – column matrix of current density induced coefficients, [ F i ] – column matrix of imposed current density, [ A i ] – column matrix of potential magnetic vector. Equation (2.11) allows calculation of magnetic potential vector in all nodes from domain.

4 Finite Element in Vector Fields To calculate the electromagnetic field was used software package Vector Fields Opera [5] where the finite element method is implemented. It includes programs for analyzing the electromagnetic field in plane and space. Opera-3D includes two programs of modeling: Geometric Modeller and PreProsessor, a computer program Post-Processor and eight analysis programs. All programs use finite element method for solving partial differential equations. Modeling programs, Geometric Modeler and Pre-Processor provides facilities for creating finite element models with complex geometry conductors; define the characteristics of material (including nonlinear, anisotropic and hysteretic materials). Manner of realization the mesh is different in the modeling programs; Geometric Modeller uses finite element tetrahedron type and Pre-Processor finite element parallelepiped type. To achieve a volume discretization in Geometric Modeller finite element size is inserted and in Pre-Processor variables as the number of elements chosen coordinate system.

Fig.2 Domain of magnetic field analysis In calculations the current value was considered 1A, and the results are presented as graphics or table. Optimizing the numerical model requires consideration of the possibility that the existence of magnetic field to be linear or nonlinear, for which the paper analyzes the two cases.

4.1 Homogeneous and linear medium In Table 1 shows comparative values for magnetic field strength on the cylinder radius, inside and the middle of turn (rotation angle τ = 180 and h/2), using the two modeling programs (Geometric Modeller and Pre-Processor), when the outside domain is extended to b = 10a, b = 20a, b = 40a, b=40a, b=80a and 21h (Fig. 2) considering a mesh step higher. Analyzing Table 1 it can see that extending the magnetic field values change in small limits, so that it can be considered the extension to 10 times the cylinder radius as sufficient.

Fig.1 Helical turn and cylinder

ISBN: 978-1-61804-055-8

102

Mathematical Models and Methods in Modern Science

b[mm] a[mm] 0 0.135 0.27 0.405 0.54 0.675 0.81 0.945 1.08 1.215 1.35

13.5 65.5263 67.8484 71.4382 75.0081 78.8110 82.9103 87.3467 92.1501 97.3465 102.9823 109.2121

Table 1. Magnetic field intensity 27 54 108 13.5 27 54 Pre-Processor Geometric Modeller 65.6795 67.9871 71.5436 75.0708 78.8203 82.8560 87.2192 91.9401 97.0451 102.5831 108.7243

65.8725 68.1623 71.6806 75.1604 78.8514 82.8174 87.1005 91.7311 96.7359 102.1651 108.1984

66.0336 68.3068 71.7911 75.2298 78.8708 82.7779 86.9932 91.5471 96.4663 101.8019 107.7380

61.4967 64.1050 66.6504 68.6194 71.7046 74.7908 77.8751 80.0020 82.1113 85.8394 89.6525

61.4782 64.0852 66.6284 68.5908 71.6890 74.7883 77.8856 80.0128 82.1227 85.8625 89.6878

61.4871 64.0974 66.6385 68.5931 71.6732 74.7542 77.8334 79.9852 82.1202 85.8640 89.6924

108 61.4819 64.0903 66.6308 68.5841 71.6722 74.7614 77.8488 79.9848 82.1041 85.8314 89.6430

3 - the mesh step of radius is a/10, respectively (b-a)/10, of the length domain analyzed is h/20, h2/20 and h1/20 and of the angle from the centre τ/10; 4 - the mesh step of the angle from the centre is τ/20, of radius is a/10, respectively (b-a)/10 and of the length domain analyzed is h/20, h2/20 and h1/20; 5 - the mesh step of radius is a/20, respectively (b-a)/20, of the length domain analyzed is h/20, h2/20 and h1/20 and of the angle from the centre τ/10; 6 - the size of finite element (introduced in Geometric Modeller) from the exterior domain is twice than that from inside of cylinder radius a (to enter value 2, respectively 1); 7 - the size of finite element (introduced in Geometric Modeller) from the exterior domain is twice than that from inside of cylinder radius a (to enter value 1, respectively 0.5). Analyzing Fig. 3 it is noticed that using finite elements of different shapes resulting differences between magnetic field intensity values obtained in the two types of modeling. It can be concluded that variant 7 is the best because Geometric Modeller uses the tetrahedron type finite elements, so modifying the finite element size it changes the three spatial variables (r, z, τ) simultaneously. In Table 2 are shown the number of finite elements and the number of nodes used in each of the following analysis.

Next is considered the field extended up to 10 times the cylinder radius and analyze how the mesh step influences the magnetic field intensity values. Fig.3 presents the magnetic field intensity resulted from the two programs, where notations have the following meanings: 1 - the mesh step of the radius is a/10, respectively (b-a)/10, of the length domain analyzed is h/10, h2/10 and h1/10 and of the angle from the centre τ/10; 2 - the mesh step of radius is a/20, respectively (b-a)/20, of the length domain analyzed is h/10, h2/10 and h1/10 and of the angle from the centre τ/10;

Fig. 3. Magnetic field intensity on radius

Table 2. Number of elements and number of nodes No. of element No. of nodes

1 60000 62031

ISBN: 978-1-61804-055-8

2 240000 248031

Pre-Processor 3 4 120000 240000 122061 244061

103

5 480000 488061

Geometric Modeller 6 7 918755 9458719 169075 1640791

Mathematical Models and Methods in Modern Science

variants obtained with modeling program PreProcessor. Analyzing Fig.5 results that extending the domain in which it analyzes the magnetic field the values of intensity increase. From a certain value increasing become insignificant and can be considered the extension to 54 mm to be sufficient.

Comparing the variant 7 with all the other is noticed a number of finite elements and nodes greater than at least 10 times the other variants. For this reason we considered that variant 7 is the most accurate, reason the other results were compared with those obtained in variant 7. Further is analyzed how conductor length affects the magnetic field intensity values. It is considered 1 turn, 3 turns and 5 turns(Fig.4).

Fig. 5 Magnetic field intensity on radius Next is considered the field extended up to 40 times the cylinder radius and analyze how the mesh step influences the magnetic field intensity values. Fig.6 presents the magnetic field from the two programs, where notations have the following meanings: 1 - the mesh step of the radius is a/10, respectively (b-a)/10, of the length domain analyzed is h/10, h2/10 and h1/10 and of the angle from the centre τ/10; 2 - the mesh step of radius is a/20, respectively (b-a)/20, of the length domain analyzed is h/10, h2/10 and h1/10 and of the angle from the centre τ/10; 3 - the mesh step of radius is a/10, respectively (b-a)/10, of the length domain analyzed is h/20, h2/20 and h1/20 and of the angle from the centre τ/10; 4 - the mesh step of radius is a/20, respectively (b-a)/20, of the length domain analyzed is h/20, h2/20 and h1/20 and of the angle from the centre τ/10; 5 - the size of finite element (introduced in Geometric Modeller) from the exterior domain is ten times greater than that from inside of cylinder radius a (to enter value 10, respectively 1); 6 - the size of finite element (introduced in Geometric Modeller) from the exterior domain is twice than that from inside of cylinder radius a (to enter value 1, respectively 0.5).

Fig.4 Magnetic field intensity on radius where notations have the following meanings: 1M, 3M and 5M represent 1 turn, 3 turns and 5 turns obtained with modeling program Geometric Modeller respectively 1P, 3P and 5P represent 1 turn, 3 turns and 5 turns obtained with modeling program Pre-processor. From Fig. 4 results that starting with 3 turns the magnetic field value no longer changes, reason that the analysis of a long conductor, with a large number of turns, it can be made considering only 3 turns. 4.2. Inhomogeneous and nonlinear medium If the environment was considered nonlinear model analysis is done from the same points of view, introducing the magnetization curve of steel over which is wrapped the helical turn. Fig. 5 shows values for magnetic field strength on to the cylinder radius, inside and the middle of turn (rotation angle τ = 180 and h/2), using the two modeling programs (Geometric Modeller and Pre-Processor), when the outside domain is extended the b = 10a (13.5), b = 20a (27), b = 40a (54), b=80a (108) and 21h (Fig. 2) considering a mesh step higher. In Fig.5 notations have the following meanings: iM, i = 13.5, 27, 54, 108 represent variants obtained with modeling program Geometric Modeller respectively iP, i = 13.5, 27, 54, 108 represent

ISBN: 978-1-61804-055-8

104

Mathematical Models and Methods in Modern Science

requires appropriate choice of field sizes considered, the mesh step, namely the finite element size and conditions of the border. Magnetic field strength decreases to the ferromagnetic nonlinear medium comparative with air due to phenomenon magnetization. Because the ferromagnetic nonlinear medium, µ isn't constant, it is necessary introducing into program the magnetization curve B(H). ACKNOWLEDGMENT This work was partially supported by the strategic grant POSDRU/88/1.5/S/50783, Project ID 50783 (2009), co-financed by the European Social Fund – Investing in People, within the Sectored Operational Programmed Human Resources Development 2007-2013.

Fig. 6 Magnetic field intensity on radius From Fig. 6 it can be concluded that variant 6 is the best because Geometric Modeller uses the tetrahedron type finite elements. Further is analyzed how conductor length affects the magnetic field intensity values. It is considered 1 turns, 3 turns, and 5 turns.

References: [1]

[2]

[3]

[4] [5]

[6]

Fig.7 Magnetic field intensity on radius Figure 7 shows that using Geometric Modeller five turns are not enough in the analysis of a long conductor and for Pre-Processor considering conductor with three turns is enough.

[7]

4 Conclusion

[8]

From the study results the following important conclusions: -Numerical model that uses finite elements form a tetrahedron leads to the best results but require the longest calculation time (about two hours and variant a only 1 minute and the 24 seconds); -The numerically analyzed and implemented in the software package Vector Field Opera allows precise calculation of the magnetic field, but this

ISBN: 978-1-61804-055-8

[9]

[10]

105

Căta, D. Păunescu, D. Toader, Calculation of magnetic field intensity vector in the axis of a wire cable with helically wound wires, Scientific Bulletin of the Politehnica, University of Timisoara, Romania, Transactions on Mathematics and Physics, Tom 56(70), Fascicule 2, 2010, p.73-85X2. T. Tominaka, Calculations using the helical filamentary structure for current distributions of a six around one superconducting strand cable and a multifilamentary composite, J. Appl. Phys. 96 5069– 80, 2004 T. Tominaka, Analytical field calculations for various helical conductors, IEEE Trans. Appl. Supercond. 14 1838–41, 2004 Vector Field Opera 13.0 User Guide. Şora C., De Sabata I., Bogoevici N., Heler A., Daba D., Vetreş I., ş.a., The Foundations of Electrotechnics, Publisher Polytechnic Timisoara 2008. Mîndru G., Rădulescu M., Numerical analysis of electromagnetic field, Lithography Institutului Polytechnic Cluj Napoca, 1983. S. Thongmak, W. Kunpasuruang, Ch. Rattanakul, A Delay-Differential Equations Model of Bone Remodeling Process, Proceedings of the 5th International Conference on Applied Mathematics, Simulation, Modelling (ASM '11) Corfu Island, Greece, July 14-16, 2011, pp. 64-68 J. Donea , A. Huerta, Finite element methods for flow problems, John Wiley & Sons, Chichester, 2003. E. Madenci, I. Guven, The finite element method and application in engineering using Ansys, Springer Science Business media LLC, ISBN 0387-28289-0, New York, S.U.A., 2006 O. C. Zienkiewicz, R. L. Taylor, The Finite Element Method, Fifth Edition, Butterworth – Heinemann, ISBN 978-0-7506-5049-6, Oxford, UK, 2000