New Scientific Contribution on the 2-D Subdomain Technique in Polar Coordinates : Taking into Account of Iron Parts

Abstract: This paper presents a new scientific contribution to the two-dimensional red(2-D) subdomain technique in polar coordinates taking into account the finite relative permeability of the ferromagnetic material. The constant relative permeability corresponds to the linear part of the nonlinear B(H) curve. As in the conventional technique, the separation of variables method and the Fourier series are used for the resolution of magnetostatic Maxwell equations in each region. The general solutions of the magnetic field in subdomains, as well as the boundary conditions (BCs) between regions are different from the conventional method. In the proposed method, the magnetic field solution in each subdomain is a superposition of two magnetic quantities in the two directions (i.e., rand Θ-axis), and the BCs between two regions are also in both directions. For example, the scientific contribution has been applied to an airor iron-cored coil supplied by a constant current. The distribution of local quantities (i.e., the magnetic vector potential and flux density) has been validated by a corresponding 2-D finite-element analysis (FEA). The obtained semi-analytical results are in very good agreement with those of the numerical method.

Currently, the works on design are based on (semi-)analytical models (i.e., equivalent circuit, conformal transformation and Maxwell-Fourier methods).This type of model consists of a (non)linear system of N analytical equations solved analytically or numerically.In comparison with the other methods, under certain geometrical and physical assumptions, these models permit obtaining accurate analytical expressions of the magnetic field and are known as fast for the local/global electromagnetic performances prediction.At present, Maxwell-Fourier methods are one of the most used semi-analytic approaches with very accurate results (i.e., error less than 5%) on the electromagnetic performances calculation.These models are based on the formal resolution of Maxwell's equations in Cartesian, cylindrical or spherical coordinates by using the separation of variables method and the Fourier's series.Taking into account iron parts and/or the effect of local/global saturation is still a scientific challenge in Maxwell-Fourier methods, which is rarely explored in the literature [16][17][18].Recently, Dubas et al. (2017) [1] realized an overview of the existing (semi-)analytical models in Maxwell-Fourier methods with the effect of local/global saturation, which can thus be classified as follows: • Multi-layer models (i.e., Carter's coefficient [19,20], saturation coefficient [21,22], concept wave impedance [23][24][25][26] and convolution theorem [27][28][29][30]); • Eigenvalues model, viz., the method of truncation region eigenfunction expansions (TREE) [31,32]; • Subdomain technique [1,33,34]; • Hybrid models, viz., the analytical solution combined with numerical methods [35,36] or (non)linear magnetic equivalent circuit [37][38][39].
The consideration of the effect of local/global saturation appears in hybrid models, where the solution is established analytically in concentric regions of very low permeability (e.g., air-gap and magnets), and other methods (e.g., numerical or magnetic equivalent circuit) are sought in regions where the saturation effect cannot be neglected.The other models (i.e., multi-layers models, TREE method and subdomain technique) are more focused on the global saturation.Some details and (dis)advantages of these techniques can be found in [1].In most semi-analytical models based on the subdomain technique, the iron parts are considered to be infinite permeable due to the variation of material proprieties in the various directions, so that the saturation effect is neglected [16][17][18].The first paper introducing the iron parts in the magnetic field calculation by using the subdomain technique is [1], where the authors solve partial differential equations (PDEs) of the magnetic potential vector in Cartesian coordinates in which the subdomains connection is performed directly in both directions (i.e., xand y-edges).The 2-D magnetostatic model has been applied to an air-or iron-cored coil supplied by a constant current.In [33], the authors propose a 2-D semi-analytical model in spoke-type magnet synchronous machines based on the subdomain technique in polar coordinates with the Taylor polynomial of degree three by focusing on the consideration of iron.The iron magnetic permeability is supposed constant corresponding to the linear zone of the nonlinear B(H) curve.The subdomains' connection is carried out in both directions (i.e., rand Θ-edges).The general solution of the magnetic field is obtained by using the traditional boundary condition (BCs), in addition to new radial BCs (e.g., between the magnets and the rotor teeth, between the teeth and the slots of the stator), which are traduced into a system of linear equations according to Taylor series expansion.In [34], this semi-analytical model has been extended taking into account the initial magnetization curve in each soft-magnetic subdomain by an iterative procedure.
In the literature, to the authors' knowledge, there exists no exact 2-D subdomain technique in polar coordinates taking into account iron parts with(out) the nonlinear B(H) curve and not using the Taylor series expansion to satisfy the r-edges BCs.Thus, in this paper, the research work contributes to the continuous improvement of the 2-D subdomain technique.Moreover, it is an extension of [1] in polar coordinates (r, Θ).Section 2 presents this new scientific contribution.By applying the principle of superposition on the magnetic quantities in order to respect the BCs on the various edges, the general solution of the magnetic field is decomposed in Fourier's series into two general solutions in both directions (i.e., rand Θ-edges).It allows the evaluation of the local distribution of flux densities in the iron parts with a global saturation, does not have numerical convergence problems contrary to others models and would easily introduce the current penetration effect in the conductive materials.The semi-analytical solution is exact as in [1] and does not use the Taylor polynomial to satisfy the r-edges BCs contrary to [33,34].For example, it was applied to an air-or iron-cored coil supplied by a constant current.The iron magnetic permeability is constant corresponding to the linear zone of the nonlinear B(H) curve [1,33].Nevertheless, as in [29,30,34], the saturation effect could be taken into account by an iterative calculation considering, at each iteration, a constant relative magnetic permeability according to the nonlinear B(H) curve.However, this is beyond the scope of the paper.In Section 3, in order to confirm the effectiveness of the proposed technique, all semi-analytical results are then compared to those found by 2-D finite-element analysis (FEA) [40].The comparisons are very satisfying in amplitudes and waveforms.

Model Description and Assumptions
Figure 1 represents the physical and geometrical parameters of an air-or iron-cored coil with N t turns of copper wire supplied by a constant current I.The electromagnetic device is surrounded by an infinite box with a null value of magnetic vector potential at it boundaries.
Physical and geometrical parameters (see Table 1) of air-or iron-cored coil where ⊗ and are respectively the forward and return conductor.The variables are: N t the number of turns; I the supply current; S c the conductor surface; u v the vacuum magnetic permeability; u c the copper magnetic permeability; and u iron the iron magnetic permeability.
The analytical prediction of the magnetic field based on the 2-D subdomain technique is done by solving magnetostatic Maxwell equations in polar coordinates (r, Θ) with the following assumptions: • The magnetic vector potential has only one component along the z-axis (i.e., A = {0; 0; A z }), and then, the end-effects are not considered; • All materials are isotropic, and the permeabilities are supposed as constants in both directions (i.e., rand Θ-axis); • All electrical conductivities of materials are supposed as nulls (i.e., the eddy-currents induced in the copper/iron are neglected).

Problem Discretization in Regions
In Figure 2, we present the studied electromagnetic device, which is divided into seven regions with µ = C st , viz., • Region 5 (i.e., the air or iron in the middle of the coil for the air or µ 5 = µ iron for the iron; • Region 6 (i.e., the forward conductor

Governing PDEs in Polar Coordinates: Laplace's and Poisson's Equations
According to Equation (A1) (see Appendix A), the distribution of the magnetic vector potential in polar coordinates (r, Θ) is governed by: ) where J zk is the current density of the coil defined by: in which S c is the conductor surface and C k (with C 6 = 1 and C 7 = −1) is the coefficient that represents the current direction in the conductor.According to Appendix A, the resolution of Laplace's and Poisson's equations by using the separation of variables method and Fourier's series permit obtaining two potentials in both directions, viz., A Θ z• for the Θ-edges (in Equation (A2b)) and A r z• for the r-edges (in Equation (A2c)).The spatial frequency (or periodicity) of A Θ z• and A r z• is respectively defined by β• h• and λ• n• with h• and n• the spatial harmonic orders.

Definition of Boundary Conditions
In electromagnetics, the general solutions of various regions depend on the BCs at the interface of two surfaces, which are defined by the continuity of the normal flux density B ⊥ and parallel field intensity H [1]. On the outer BCs for (Θ 1 ∧ Θ 6 , ∀r) and (∀Θ, r 1 ∧ r 4 ), A z satisfies the Dirichlet BC (see Figure 2), viz., A z = 0.
Figure 3 represents the respective BCs at the interface between the various regions in both directions (i.e., rand Θ-edges).

Region 2
The same method as Region 1 is used to define the general solution in Region 2. By posing d Θ h = 0 in Equation (A11) (see Appendix B), A z2 satisfying the BCs of Figure 3b and the solution of Equation (1a) is given by: the components of B 2 = {B r2 ; B Θ2 ; 0} by: where h2 is the spatial harmonic orders in Region 2, c2 Θ h2 the integration constant and   Using a Fourier series expansion of G 2 (Θ) (see Figure 3b) over the interval

Region 3
The solutions of A z3 , B r3 and B Θ3 are determined by the Case-Study No. 1 (i.e., A z imposed on all edges of a region) in Appendix B. The BCs on the Θ-edges of the region (see Figure 3c) are met by posing e r n = 0 in Equations (A6)-(A8).Therefore, A z3 satisfying the BCs of Figure 3c and the solution of Equation (1a) is given by: the r-component of B 3 by: the Θ-component of B 3 by: where h3 and n3 are the spatial harmonic orders in Region 3; c3 Θ h3 , d3 Θ h3 and f 3 r n3 the integration constants; Using Fourier series expansion of A z1 | ∀Θ∧r=r 2 and A z2 | ∀Θ∧r=r 3 (see Figure 3c) over the interval With a weighting function g (r) = r −1 and using a Fourier series expansion of A z6 | Θ=Θ 2 ∧∀r (see Figure 3c) over the interval r = [r 2 , r 3 ], the integration constant f 3 r n3 is determined in Appendix C with: 2.5.4.Region 4 The solution in Region 4 is obtained using the same development as Region 3. By posing f r n = 0 in Equations (A6)-(A8) (see Appendix B), A z4 satisfying the BCs of Figure 3d and the solution of Equation (1a) is given by: the r-component of B 4 by: the Θ-component of B 4 by: where h4 and n4 are the spatial harmonic orders in Region 4; c4 Θ h4 , d4 Θ h4 and e4 r n4 the integration constants; β4 h4 = h4 • π τ Θ4 with τ Θ4 = Θ 6 − Θ 5 ; and λ4 n4 = n4 • π/τ r4 with τ r4 = ln (r 3 /r 2 ).

Introduction
The objective of this section is to validate the 2-D subdomain technique in polar coordinates (r, Θ) on the magnetic field distribution in relation to the numerical method.The physical and geometrical parameters of studied electromagnetic device are given in Table 1.
Table 1.Physical and geometrical parameters of the air-or iron-cored coil.

Parameters, Symbols (Units) Values
Number of turns of the coil, N t (-) 60 Supply current, I (A) 20 Conductor surface, S c (mm 2 ) 120 Current density of the coil, J zk (A/mm 2 ) ±10 Effective axial length, L z (mm) 60 Geometrical parameters in the Θ-axis, {Θ 1 ; Θ 2 ; Θ 3 ; Θ 4 ; Θ 5 ; Θ 6 } (deg.){0; 17; 21; 29; 33; 50}  Geometrical parameters in the r-axis, {r 1 ; r 2 ; r 3 ; r 4 } (mm) {21; 81; 100; 160} Relative magnetic permeability of the iron, µ iron (-) 1,500 Number of harmonics for Region 1, H1 max (-) 260 Number of harmonics for Region 2, H2 max (-) 260 Number of harmonics for Region 3, {H3 max ; N3 max } (-)  {88; 124} Number of harmonics for Region 4, {H4 max ; N4 max } (-)  {88; 124} Number of harmonics for Region 5, {H5 max ; N5 max } (-)  {42; 124} Number of harmonics for Region 6, {H6 max ; N6 max } (-)  {21; 124} Number of harmonics for Region 7, {H7 max ; N7 max } (-)  {21; 124} For this validation, the air-or iron-cored coil has been modeled using Cedrat's Flux2D (Version 10.2.1, Altair Engineering, Meylan Cedex, France) software package (i.e., an advanced finite-element method-based numeric field analysis program) [40].FEA is done with the same assumptions as the semi-analytical model (see Section 2.1).The linear system is given in Appendix C and has been implemented in MATLAB R (R2015a, MathWorks, Natick, MA, USA) by using the sparse matrix/vectors.A discussion of the numerical problems (viz., harmonics and ill-conditioned systems) of such semi-analytical models has been clarified in [1].The Maxwell-Fourier methods exhibit a similar problem to the numerical methods due to the periodicity of Fourier series and, consequently, to the finite number of harmonics.Hence, A z and B = {B r ; B Θ ; 0} in the various regions (see Section 2.5) have been computed with a finite number of spatial harmonics terms H1 max -H7 max (for the Θ-edges) and N3 max -N7 max (for the r-edges).As indicated in [41,42], these spatial harmonics terms, given in Table 1, have been imposed according to an optimal ratio, i.e., for H1 max given, The linear system size depends on the number of: (i) regions; (ii) BCs; and (iii) harmonics of each subdomain.In our study, the linear system named Equation (A17) consists of 2036 elements, which is much smaller than the 2-D FEA mesh having 3,081 surfaces elements of second order (viz., the triangles number of system) with the number of excellent quality elements equal to 100%.For information, the 2-D FEA mesh for an air-or iron-cored coil is illustrated in Figure 4.The personal computer used for this comparison has the following characteristics: HP Z800 Intel(R) Xeon(R) CPU @ 2.4 GHz (with two processors) RAM 16 Go 64 bits.The computation time of 2-D subdomain model is divided by two (viz., 0.5 s for 2-D subdomain model and 1 s for the 2-D FEA).

Results Discussion
The validation paths of A z and B = {B r ; B Θ ; 0} for the semi-analytic and numeric comparison are given in Figure 5.The waveforms of global quantities are shown on different paths in Figure 6 for A z and in Figures 7-11 for the components of B. The solid lines represent the global quantities computed by the 2-D FEA, and the circles correspond to the 2-D subdomain model.Comparing those results with 2-D FEA, it can be shown that a very good evaluation is obtained for A z and for the components of B, whatever the paths, for both the air-and iron-core.This confirms that the effect of global saturation can be taken into account accurately.According to the concept of symmetry, it can be seen that in polar coordinates, there is only one symmetry of A z on Path 5 (viz., B r = 0 and B Θ = 0 in Figure 11) unlike the same electromagnetic device in Cartesian coordinates.Indeed, in Cartesian coordinates, there exist two symmetry axes of A z on Path 2 and Path 5 [1].In polar coordinates, Path 2 does not correspond to a symmetry axis of A z ; consequently, B r = 0 and B Θ = 0 in Figure 8.It will be noted that the r-component of B in Region 5 is more intense with the magnetic core (see Figure 8a) and that the magnetic leakages in the middle of the device in the Θ-axis are equivalent for an air-and iron-cored coil (see Figure 8b).It is interesting to note that numerical peaks appear in the FEA results (see Figures 6e, 7, 8b and 11b), which are mainly due to the mesh.The relative error is less than 1.5% for the various global quantities (see Figure 6a,c for the maximum error).

Conclusions
It has been demonstrated that there exists no exact semi-analytical model based on the 2-D subdomain technique in polar coordinates taking into account iron parts with(out) the nonlinear B(H) curve.An improved 2-D subdomain method in polar coordinates (r, Θ) to study the magnetic field distribution in the iron parts with a finite relative permeability has been presented in this paper.Nevertheless, the research work is an extension of [1] in polar coordinates (r, Θ).
The proposed new subdomain model is applied to an air-or iron-cored coil supplied by a constant current.The magnetic field solutions in the subdomains and the BCs between regions are carried out in the two directions (i.e., rand Θ-axis).The iron relative permeability used in this model is constant and corresponds to the linear part of the nonlinear B(H) curve.However, the whole B(H) curve of the magnetic material can be applied with an iterative algorithm as in [29,30,34].The proposed subdomain method in polar coordinates (r, Θ) takes less computing time than the FEA (approximately two-fold versus to FEA).It is very suitable for the design and optimization of the electromechanical systems in general and electrical machines in particular.The semi-analytical results have been validated with FEA, and good agreement has been obtained in both amplitudes and waveforms.
The major scientific contribution could be applied to rotating electrical machines (e.g., radial-flux machines) in polar coordinates with(out) magnets supplied by a direct or alternate current (with any waveforms).Moreover, one advantage of this technique would be their exploitation in three-dimensional studies in order to reduce the computation time, which remains a major problem in this numerical method.
where A zP is the particular solution of A z respecting the second member ES in Equation (A1), C Θ 0 -F Θ h and C r 0 -F r h the integration constants, β h and λ n the spatial frequency (or periodicity) of A Θ z and A r z and h and n the spatial harmonic orders.
Using B = ∇ × A, the components of magnetic flux density B = {B r ; B Θ ; 0} can be deduced by: which leads to: and: Figure A1a shows a region (for Θ ∈ [Θ r , Θ l ] and r ∈ [r l , r t ]) whose A z is imposed on all edges.By respecting the BCs and applying the principle of superposition on the magnetic quantities, Figure A1a is redefined by Figure A1b.In Case-Study No. 1, A z = A Θ z + A r z , i.e., Equation (A2), is redefined by: the component B r = B Θ r + B r r of B, i.e., Equation (A4), by: and the component B Θ = B Θ Θ + B r Θ of B, i.e., Equation (A5), by: where c Θ h , d Θ h , e r n and f r n are new integration constants; with τ r = ln (r t/ r l ); and E ¡ σ (w, x, y) and P ¡ σ (w, x, y) are [44]: When A z = 0 on Θ-edges and A z imposed on r-edges (see Figure A2), A z with A r z = 0 in Equation (A6) is expressed by: the r-component of B with B r r = 0 in Equation (A7) by: the Θ-component of B with B r Θ = 0 in Equation (A8) by: Figure A3a shows a region (for Θ ∈ [Θ r , Θ l ] and r ∈ [r l , r t ]) whose B r and A z are respectively imposed on rand Θ-edges.By respecting the BCs and applying the principle of superposition on the magnetic quantities, Figure A3a is redefined by Figure A3b.In Case-Study No. 2, A z = A Θ z + A r z , i.e., Equation (A2), is redefined by: the component B r = B Θ r + B r r of B, i.e., Equation (A4), by: and the component B Θ = B Θ Θ + B r Θ of B, i.e., Equation (A5), by: where c Θ 0 , d Θ 0 , c Θ h , d Θ h , e r n and f r n are new integration constants.

5 Figure 2 .
Figure 2. Definition of regions in the air-or iron-cored coil.

Figure 5 .
Figure 5. Validation paths for the semi-analytic and numeric comparison.

)
Appendix B. Simplification of Laplace's Equations According to Imposed BCs Appendix B.1.Case-Study No. 1: A z Imposed on all Edges of a Region

Figure A1 .
Figure A1.A z imposed on all edges of a region: (a) general and (b) principle of superposition.

Figure A2 .
Figure A2.Particular case:A z = 0 on Θ-edges and A z imposed on r-edges of a region.

Figure A3 .
Figure A3.B r imposed on r-edges and A z imposed on Θ-edges of a region: (a) general and (b) principle of superposition.