Investigation of Volumic Permanent-Magnet Eddy-Current Losses in Multi-Phase Synchronous Machines from Hybrid Multi-Layer Model

: This paper investigates the permanent-magnet (PM) eddy-current losses in multi-phase PM synchronous machines (PMSM) with concentric winding and surface-mounted PMs. A hybrid multi-layer model, combining a two-dimensional (2-D) generic magnetic equivalent circuit (MEC) with a 2-D analytical model based on the Maxwell–Fourier method (i.e., the formal resolution of Maxwell’s equations by using the separation of variables method and the Fourier’s series), performs the eddy-current loss calculations. First, the magnetic ﬂux density was obtained from the 2-D generic MEC and then subjected to the Fast Fourier Transform (FFT). The semi-analytical model includes the automatic mesh of static/moving zones, the saturation effect and zones connection in accordance with rotor motion based on a new approach called "Air-gap sliding line technic". The results of the hybrid multi-layer model were compared with those obtained by three-dimensional (3-D) nonlinear ﬁnite-element analysis (FEA). The PM eddy-current losses were estimated on different paths for different segmentations as follow: (i) one segment (no segmentation), (ii) ﬁve axial segments , and (iii) two circumferential segments, where the non-uniformity loss distribution is shown. The top of PMs presents a higher quantity of losses compared to the bottom.


Introduction
Currently, there is trend of the electrification of transport applications, where PMSM are increasingly being used. The reason is their good torque density, lower losses and higher power-to-mass ratio [1]. In addition to advantages associated to this topology, a multi-phase system allows continuous operation in degraded mode, which is a powerful asset in traction application [2,3].
One of the major concerns in PMSM is the PM demagnetization risk. This occur by temperature rise caused by the eddy-currents [4]. The high-order magnetic flux density harmonics, which are not synchronous with the rotor, cause those parasitic eddy-current losses in the PMs [5]. Eddy-current losses may be separated into two cases according to the operating point: • No-load (without stator current): the reluctance variation due to the stator slot-opening [5,6]; • Load (with stator current): the spatio-temporal magnetomotive force (MMF) harmonics [7,8].
In design process, the PM eddy-current losses must be calculated. Several methods have been developed and categorized into three families viz.: (i) (Semi-)analytical methods based the formal resolution of Maxwell's equations [9,10], (ii) numerical methods [11,12], and (iii) hybrid methods [13,14]. An overview on the eddy-current loss calculation has been realized in [15], where several investigations have been chronologically listed. In [16], sources, calculation methods, reduction techniques, and thermal analysis of PM eddy-current losses were revised. A review of (semi-)analytical models based on the subdomain technique for PM eddy-current calculation was made by [17].
The fact that the eddy-currents flow is along the loops inside the PM volume leads to the classification of their calculation as a 3-D problem [18]. The 2-D calculation ignores end-effects which cause a large error in the PM eddy-current computation. To consider end-effects and obtain a more accurate calculation, a 3-D model must be used [19].
This paper investigates the PM eddy-current losses in massive conducting parts of a multi-phase PMSM from hybrid multi-layer model. The latter consists of combining two models, viz.: (i) a 2-D generic MEC, and (ii) a 2-D analytical model based on the Maxwell-Fourier method (i.e., the formal resolution of Maxwell's equations by using the separation of variables method and the Fourier's series). The 2-D generic MEC determines the magnetic flux density distribution in the massive conductive parts (e.g., the PMs,...) without the skin effect (i.e., without the eddy-current reaction field). Moreover, this model can give the distribution of other local/integral physical quantities in the electrical machine with the nonlinear of B(H) curve. The 2-D analytical model based on the Maxwell-Fourier method calculates the magnetic field distribution in massive conductive parts considering the skin effect as well as the resultant eddy-current density. The boundary conditions (BCs) imposed on the 2-D analytical model are equivalent to the magnetic field obtained from the MEC. These BCs are considered homogeneous on the edges and are only equal to normal magnetic field in the top or middle of massive conductive parts. This is a strong hypothesis, because the magnetic field is naturally non-uniform. In our case, the model gives the ability to use both conditions (i.e., in the top or middle of massive conductive parts), allowing to see the impact on the eddy-current loss calculation.
The proposed approach is introduced in section 2 with a short focus of the 2-D generic MEC and its validation. Section 3 shows the obtained PM eddy-current losses on different paths for different segmentations. In the last section, the results of the hybrid multi-layer model are compared with those obtained by 3-D nonlinear FEA [21].

Hybrid Multi-Layer Model Description
The performed approach is composed of two main modules, viz.: (i) a 2-D generic MEC, and (ii) a 2-D analytical model based on the Maxwell-Fourier method. The first one allows the calculation of magnetic flux density B without the skin effect over all the PMSM. In the PM zones, the magnetic flux density is levied on different paths parallel in the tangential direction. Indeed, the magnetic flux density varies depending on the PM depth direction. These magnetic flux densities will be subjected to the FFT giving the amplitude as well as the frequency of different harmonic numbers. Furthermore, these outputs constitute the main input parameters to estimate the PM eddy-current losses by using the developed formal Maxwell model. Note that only the harmonics are going to be taken into account when calculating the eddy-current losses. The eddy-current losses due to the average value are null because the resultant eddy-current density in the PMs does not exist. The developed approach for the eddy-current loss calculation is explained in the flowchart shown in Figure 1.

Geometrical and Physical Parameters of PMSM
The performed model was applied to a five-phase surface-mounted PMSM, where the geometrical and physical parameters are given, respectively, in Tables 1 and 2. Figure 2 gives the topology of the machine and the spatial distribution of various phases. The stator has a double-layer concentrated winding distribution (viz., the non-overlapping winding with all teeth wound), supplied by sinusoidal current waveform. The five-phases windings are star-connected. The PMs, radially magnetized, and the sheet metal are, respectively, NdFeB and M270-35A. The dashed lines represent the different paths used to proof this approach. Based on this path, the local quantities will be levied and furthermore compared to those of FEA [21].
Attention should be paid to Figure 3 which presents the different tangential paths located in both South and North PM. This path will be used to get the different values of the magnetic field along the symmetric axis of the PM.

2-D Generic MEC
In this work, the approach proposed in [22] (viz., the generalized nonlinear adaptive MEC) was applied onto the 5-phases PMSM. Indeed, the 2-D MEC is based on the discretization principle, which allows bidirectional (BD) automatic mesh reluctance generation. In a general manner, the electrical machine can be discretized in a set of mesh elements that can be, in turn, discretized in a high number of BD blocks. As an example, the mesh element spotted by the grey color in Figure 4 can be discretized in a set of BD blocks, as shown in Figure 5.
The main advantages of such a model include of the flexibility, generic, compromise between computation time and precision. The semi-analytical model includes the automatic mesh of static/moving zones, the saturation effect and zones connection in accordance with rotor motion based on a new approach called "Air-gap sliding line technic". This technique was applied in [23][24][25][26] on different electrical machine configurations, viz., axial-flux interior PM machine [23,24], coaxial magnetic gear equipped with surface-mounted PMs [25], wound-rotor synchronous machine [27], and radial-flux interior PM machine [26]. This novel technique of connection between static/moving zones as well as the approach proposed in [22] is applied for the first time on a radial-flux surface-mounted PMSM having multi-phases for automotive application. Nevertheless, it should be noted that there are other techniques permitting to connect static/moving zones in electrical machines which have been overviewed in [22].
Comparing to [26], in our case, the stator tooth-pitch will be constituted by five radial zones, and four radial zones for the rotor pole pitch. For the tangential zones, the rotor pole pitch is consists of three zones, while the stator pitch consists of five zones. The discretization number in each zones is given by Table 3.
The developed model respects the flowchart given in Figure 6, where the different steps of the nonlinear system solving are presented. To account for the saturation effect, an iterative process is considered by using the fixed point iteration method. The nonlinear B(H) curve of the iron is then introduced by using the Marrocco's interpolation function. In Figure 6, µ ri represents the initial relative permeability in the magnetic circuit and N sat is the total iterations number for taking into account the saturation effect. The 2-D nonlinear adaptive MEC (where the loop fluxes ψ are the unknowns) can be expressed by (1) where [ ] is the diagonal matrix of reluctances, [ψ] is the loop fluxes vector, [χ] is the topological (or incidence) matrix, [F] is the loop MMFs vector, and [MMF] is the branch MMFs vector. These matrixes and vectors are dependent upon the discretization of the various zones. Their dimensions are invariable with the time (or the mechanical angular position between the rotor and the stator) [22]. Finally, dividing this later by the corresponding reluctance surfaces leads to the magnetic flux density in each magnetic branch.

Validation of 2-D generic MEC
The tangential and radial components of B at t = 0 sec are given in Figures 7-12 for the different path crossing , viz.: the stator yoke, stator teeth, stator air-gap, rotor air-gap, rotor PMs and rotor yoke. The components of B in the stator teeth [see Figure 8] are very saturated. The results are in agreement with 2-D nonlinear FEA, which confirms that the nonlinear solving method is accurate. The errors are due to the discretization of static and moving zones in both tangential and radial directions. For example, a sensibility study versus stator/rotor discretization on the computation time as well as on the iron losses determined from the magnetic flux density in the stator/rotor iron have been performed in [27]. The difference is mainly linked to the discretization for the rotor yoke where Nd r 3 = 4 , which can be observed on the Figure 12a of the tangential component of B in the rotor yoke. The computation time is reduced by twice compared to the 2-D nonlinear FEA.
The simulation in time stepping leads to the results represented in Figures 13-15. The electromagnetic torque evolution at locked rotor is given in Figure 13; this help us to find the shift angle to obtain the maximal electromagnetic torque at synchronism shown in Figure 14. The results give a decent satisfaction. The average electromagnetic torque given by 2-D nonlinear FEA is equal to 32 Nm while the 2-D nonlinear adaptive MEC is slightly above 30 Nm. The error is equal to 6.25 %. This small difference can be linked to the approach used to take into account the magnetic saturation, the rotor motion as well as the air-gap, stator and rotor discretization.
As mentioned in Section 2, the performed model will be used to extract the local physical quantities in PMs in order to estimate the eddy-current losses. Moreover, it should be noted that B θ in the PMs is considered negligible in relation to B r (see Figure 11). The losses are then determined from the normal magnetic field in the PMs. The model can give the waveform of B according to the time or the rotor position. The simulations are done for the path crossing the middle of PMs for one periodic pattern as shown in Figure 15a. The results are in agreement with the 2-D nonlinear FEA. Figure 15b gives the magnetic flux density evolution in the middle of North and South PMs.
Note that the average value of B is assumed to be not considered for the eddy-current loss calculation. The evolution of ∆B r for both middle of North and South PMs is shown in Figure 15b. It can be seen that the variation of B is practically identical. Thus, it can be concluded that the losses in North and South PMs are identical. Furthermore, only the results given for the South PM are given.

General Solution of Magnetic Field with the Skin Effect
The magnetic field with the skin is determined from an analytical model in 2-D Cartesian coordinates (i.e., the PMs are assumed rectangular) based on the Maxwell-Fourier method (i.e., the formal resolution of Maxwell's equations by using the separation of variables method and the Fourier's series). This permits to estimate the resultant eddy-current density in PMs in both directions (i.e., axial and circumferential axis). The BCs imposed on the 2-D analytical model are equivalent to the magnetic field obtained from the MEC. These BCs are considered homogeneous on the edges of massive conductive parts, which are equal to H s . This is a strong hypothesis because the magnetic field is naturally non-uniform. The value of magnetic field is defined as the normal magnetic field in the top or middle of PMs determined by 2-D generic MEC i.e., without the skin effect, viz., H s = B r (0, 0)/µ mp . In our case, the model gives ability to use both conditions (i.e., in the top or middle of massive conductive parts), allowing to see the impact on the eddy-current loss calculation. Furthermore, the results will be given for different paths parallel to the tangential direction. The various BCs are shown in Figure 16.
In quasi-stationary approximation, inside a linear (non)magnetic material of electrical conductivity without electromagnetic sources, the partial differential equation in magnetodynamic in term of H mp can be defined by [28] Using the complex notation, the magnetic field H mp = 0; H mp σy ; 0 in massive conductive parts can be written as where j = √ −1 and ω = 2π f is the electrical pulse. Therefore, (3) becomes in 2-D Cartesian coordinates the complex Helmholtz's equation, viz., By using the separation of variables method and by applying the BCs, the 2-D general solution of H mp σy in both directions (i.e., xand z-edges) can be written as Fourier's series [28] where λ k = kπ/d is the periodicity of H mp σy in the z-axis, δ k = α 2 + λ 2 k , and k are the spatial harmonic orders.
It should be noted that H mp σy is assumed to be invariant in the y-axis according to the study path. Moreover, when α = 0 (viz., σ mp = 0 S/m and f ∼ = 0 + Hz) then f σ (x, z) = 1 thus giving H mp σy = H s .

Eddy-Current Loss Formulation
The instantaneous density of power flow ∏ ∏ ∏ at a point is defined by the Poynting's vector [28] ∏ ∏ ∏ = E × H mp (7) Using J = σ · E = ∇ × H mp , this power across a closed surface, in terms of complex vectors, is given by the instantaneous apparent power The real part of the surface integral of the complex Poynting's vector gives the instantaneous ohmic losses and the imaginary part the instantaneous magnetic energy.
From BCs at edges of PMs (see Figure 16), the average of s app over an electrical cycle T = 2π/ω can be defined by where h l = h a / (2 · Nd r 2 ) is the layers thickness in the r-axis of PMs with l the path number. After development, (9) is then given by It should be noted that P = {S app } represents the average eddy-current losses in PMs.

Magnetic Flux Density vs PM Thickness
Once the 2-D nonlinear adaptive MEC was validated on different stator and rotor paths. The radial component of B in the middle of North and South PMs according to the different tangential paths was evaluated. It should be noted that the paths number is imposed by the radial discretization number applied in the PMs region. The waveform of B in the PMs is done by stepwise simulation to get the radial component evolution according to the time or the rotor position.
Under the assumptions that the continuous magnitude B m of B mp σy does not contribute to the eddy-current losses and that only the waveform ∆B r is subjected to the FFT. Hence, the most influential harmonics can be used to estimate the PMs eddy-current losses. Knowing that:=

Eddy-Current Loss Evolution in different Paths
Figures 19 and 20 show the frequency spectrum of the radial component of B for different speeds. The spectral components are due to: (i) the supply (i.e., the current waveform), (ii) the slotting effect, and (iii) the spatial distribution of winding. Using the results given by the FFT for different paths in the PM, the PMs eddy-current losses can be estimated for each path as shown in Figure 21a and Figure 21b.
By applying the sum, the total eddy-current losses due to the different paths are obtained as shown in Figure 21c. It should be noted that the eddy-current loss repartition in the PMs is greater at the top of PMs (whatever the PM magnetization) and decrease with to the PMs height to achieve the lower losses in the bottom of PMs. This distribution can be explained by the variation of the magnetic flux density amplitude that is more important in the region which is close to the air-gap, contrary to the bottom of PMs where the amplitude is relatively less important.

3-D nonlinear FEA Validation With(out) PM Segmentation
In this section, the results of the hybrid multi-layer model are compared with those obtained by 3-D nonlinear FEA [21]. The PMs eddy-current losses were computed at different speed ranges for different segmentations as follow: (i) one segment (no segmentation), (ii) five axial segments, and (iii) two circumferential segments.
To prove the validation of this approach, the results are compared with those of 3-D nonlinear FEA in transient state [21]. Due to the boundary conditions (i.e., periodicity and Dirichlet conditions), the electrical machine can be reduced into 5-slots/4-poles. The mesh generation for 3-D nonlinear FEA is equal to 2,263,818 second order elements. It is given the eddy-current loss evolution in all PMs by using different segmentations with two different approaches (see Figure 22): • 1 st approach: based on all paths; • 2 nd approach: based on the middle path of PMs. Figure 23 represents the relative error between each approach and 3-D nonlinear FEA. These errors are due to medium mesh in the 3-D nonlinear FEA as well as the magnetic field evaluation on all paths of the 2-D nonlinear adaptive MEC which is a data input the analytical model based on the Maxwell-Fourier method permitting the PM eddy-current loss calculation. In both cases, the results provide a good agreement with a maximal error of 16% that is reached for the 2 nd approach based on the middle path with two circumferential segments. While the maximal relative error for the 1st approach is estimated at 11%, which is achieved for one segment (no-segmentation) at 3, 000 rpm.
The developed model returns both results based on to the two different approaches. It is better to use the 1 st approach to estimate the eddy-current losses in the circumferential segmentation case.
Concerning the computation time, the 3-D nonlinear FEA takes 16 days for one simulation while the hybrid multi-layer model takes only 12 minutes. For an operating point, the computation time is then divided by 1920 with respect to 3-D nonlinear FEA.

Conclusion
In this work, a hybrid multi-layer model, combining a 2-D generic MEC with a 2-D analytical model based on the Maxwell-Fourier method was performed. This developed model was aimed to compute the volumic PM eddy-current losses in multi-phase PMSM. It is clearly seen that the PM eddy-current losses distribution present no-uniformity according to the PMs height. In this study case, these losses are greater at the top of PMs. This can be justified by the amplitude of the radial component of B, which is the highest in the top of PMs. Then, they decrease to achieve the minimal in the PM bottom. As a conclusion, the temperature repartition can be affected by the non-uniformity of the PM eddy-current loss distribution over the PM volume.