Analytical Modeling of Static Eccentricities in Axial Flux Permanent-Magnet Machines with Concentrated Windings

Abstract

The aim of this paper is to calculate the static eccentricity (SE) of a double rotor axial flux permanent magnet (AFPM) machine by using a general analytical model. The flux density in the air gap under healthy conditions is calculated firstly, where the axial and circumferential magnetic flux densities are obtained using a coupled solution of Maxwell’s equations and Schwarz-Christoffel (SC) mapping. The magnetic flux densities under SE conditions are calculated afterwards using a novel bilinear mapping. Some important electromagnetic parameters, e.g., back electromotive force (EMF), cogging torque and electromagnetic (EM) torque, are calculated for both SE and healthy conditions, and compared with the finite element (FE) model. As for the double rotor AFPM, SE does not contribute much effect on the back EMF and EM torque, while the cogging torque is increased. At each calculated section, FE models were built to validate the analytical model. The results show that the analytical predictions agree well with the FE results. Finally, the results of analytical model are verified via experimental results.

1. Introduction

Axial flux permanent magnet machines (AFPMMs) have a number of distinct advantages over radial flux permanent magnet machines (RFPMMs). Axial flux permanent magnet (AFPM) can be designed to have a higher power-to-weight ratio, resulting in a need for less core material and reduced complexity to adjust the air gap [1]. Moreover, because of the disc shaped rotor and stator structure, an AFPM is smaller in size than its RFPMMs counterparts. This particular structural characteristic makes it easy to address the space limitation in some applications such as electric vehicles (EV) [2], hybrid electric vehicles (HEV) [3], wind turbos [4] and flywheel storage systems [5], etc.

However, due to the smaller contact surface between the rotor and the shaft, it is more difficult to design a rotor-shaft mechanical joint with a high mechanical integrity [6], which makes AFPMMs more susceptible to imperfect assembly issues such as static eccentricities (SEs), dynamic angular misalignment [7] and static/dynamic eccentricity misalignment [8]. Under angular eccentric conditions, the rotor is inclined and the air gap is asymmetric. If the rotor shaft assembly is sufficiently rigid, which means that the level of static eccentricity (SE) does not change during periodic operation, the axial force becomes unbalanced causing vibrations and unbalanced magnetic forces (UMFs). This leads to deterioration of both the machine’s performance and lifetime. The UMFs can also act on the stator core, exposing the stator windings to unnecessary and potentially harmful vibrations and in the most extreme case, causing the stator and the rotor to scrape [9]. It is therefore important to estimate SE by measuring some suitable machine quantities.

Various techniques could be used for modeling AFPMMs under SE conditions, e.g., an analytical model [10], finite element (FE) model [11], magnetic equivalent circuit (MEC) [12] or a combination of these methods. The 3D-FE model can precisely model the influence of eccentricity owing to the accurate computation of air gap region, but it is time consuming, thus, some references have proposed alternative methods to reduce the computation time. The quasi-3D method, which is widely adopted by AFPM designers, still takes a large amount of computation time for SEs [13]. Another alternative approach is the field reconstruction (FR) method [14,15]. The FR method is a mixed analytical-FE method and partly depends on the FE model. Although, it is much faster than the 3D-FE model, it still takes a lot of time for modeling the different SE geometries. The analytical approach, due to its fast and acceptable results, is regarded as an efficient and favorable method compared to the other methods mentioned above.

Analytical solutions are categorized in different groups. One group is based on the MEC [12]. The MEC is not accurate enough for some qualities like cogging torque and forces because of the nodes limit. In addition, the model needs to be rearranged at each simulation step due to the requirement of alignment between stator and rotor reluctances. Moreover, it is more complex when dealing with the non-uniform air gap. Other possible analytical solutions are appearing in the literature [16], but the methods used in both papers ignored the circumferential flux density when modeling SE conditions. The Conformal Mapping (CM) method has so far not been used for accurate magnetic field analysis for AFPMSMs. In [17], Zarko used this method to calculate the relative permeability of the air gap considering the slotting effect, but one of the main assumptions is that the motors have an infinite slot opening which implies the shape of slot is ignored. Following this paper, Zarko improved this method in [18] to handle this problem, but the computational complexity increases as well.

To overcome the aforementioned problems, this paper adopts an accurate and improved CM method combined with the Schwarz-Christoffel (SC) method for the no-load and on-load conditions under SEs. This method depends on the Matlab SC Toolbox instead of the FE model [19].

This paper is organized as follows: in Section 2, a description of double rotor prototype AFPM with concentrated coil is provided. Section 3 introduces the electromagnetic calculations and SC transformation used in this paper. SE and its effects on magnetic field are also described in this section. Afterward, the models of both the healthy condition and the static eccentric magnetic distribution are introduced in Section 4. In Section 5, the simulation results, e.g., the back electromotive force (EMF), cogging torque and EM torque, are compared and validated through the simulations done by the FE model developed by the commercial software JMAG (https://www.jmag-international.com/). In Section 6, the experimental results verify the accuracy of this method, which is followed by a brief summary and conclusion of the current study.

2. Description of the Prototype Axial Flux Permanent Magnet Machine

Figure 1a shows the construction of the test machine used in this paper. The AFPM machine consists of two rotors and one segmented stator. Soft magnetic composite (SMC) is used for the stator core and axially magnetized NdFeB permanent magnets (PMs) are mounted on the surface of each rotor, which is made of #45 steel, the SMC modules are fastened by resin and holders. Three phase concentrated windings are distributed on each stator segment to avoid overlapping. This machine topology benefits from advantages including short end windings with high winding fill factor, small stator iron loss under high-speed operation and a convenient manufacturing and assembly process, plus reduced stator weight due to the absence of the stator yoke. The dimensions and specifications of the investigated machine are shown in

Table 1. Parameters of the 6-Pole, 9-Slot axial flux permanent magnet machine (AFPMM). SMC: Soft magnetic composite; PM: Permanent magnet.

Parameter	Symbol	Value	Unit
Rated power	P	10	kW
Rated voltage	U	200	V
Rated speed	np	20,000	rpm
Number of poles/slots	p/Qs	6/9	-
Stator outer radius	Ro	70	mm
Stator inner radius	Ri	45	mm
Air gap length	g	3	mm
Remnant flux density of PM	Br	1.03	T
Stator core	-	SMC	-
Permanent magnet	-	NdFeB	-

.

3. Calculation of the Magnetic Field

In this paper, the Quasi-3D-method is adopted in the analytical model. The machine is divided into a certain number of layers. Thus, an axial machine could be considered as a combination of several individual linear machines, which allows separate analysis of each plane. One 2D model developed in 2D polar coordinates is shown in Figure 1d and the average radius of a particular ith layer is given by:

$$R_{ave,i}=R_i+\frac{R_o-R_i}{2n_s}(2i-1)$$

(1)

where n_s is the number of slices.

To improve the accuracy of the results, six slices are considered to calculate the healthy and fault conditions in this paper. The processes of quasi 3D method and two poles of a 2D model in Cartesian coordinates are illustrated in Figure 1b–d.

As for the magnetic flux calculation of air gap, the same approach in [20] is used to calculate the magnetic induction in the air gap caused by the magnets and the armature current. Moreover, to take the slot effect into account, the SC mapping [21] is adopted.

By applying Maxwell’s equations on selected slice and by summing the air gap flux caused by the armature current and magnets, the flux distribution can be obtained. The following assumptions are made during calculations in order to reduce the complexity of the computations:

(1)

The magnetic material has a uniform magnetization and the relative recoil permeability μ_r is constant and has a value close to unity such as in NdFeB materials.

(2)

For the computation of armature reaction field, the magnet regions are regarded as free space.

(3)

Magnetic saturation is absent and the rotor iron cores have infinite magnetic permeability.

(4)

Eddy current effects are neglected, which avoids the need for the complex eddy current field formulation.

3.1. Model of Permanent Magnets (PMs)

For NdFeB, the magnetic induction in the PMs is written as:

$$B=B_r+{\mathsf{\mu }}_{rec}H={\mathsf{\mu }}_0{M}_r+{\mathsf{\mu }}_{rec}H$$

(2)

where μ_rec = μ₀ μ_r is the recoil permeability.

As shown in Figure 2, the magnetization vector is assumed to be along the axial direction and may be described by Fourier series containing only cosine terms:

$$M_y(x)={\displaystyle \sum _{n=1,3,\mathrm{5...}}^{\infty }{M}_n\cos (\frac{n\mathsf{\pi }x}{{\mathsf{\tau }}_p})}$$

(3)

where τ_p is pole pitch in circumferential direction and α_p is the ratio of magnet pole arc to pole pitch, h_m is axial thickness of magnets, L is the axial distance between rotor and back plates and M_n is shown as follows:

$$M_n=\frac{4B_r}{n{\mathsf{\pi }\mathsf{\mu }}_0}\sin (\frac{n{\mathsf{\pi }\mathsf{\alpha }}_p}{2})$$

(4)

For a PM machine with linear demagnetization characteristic, the scalar magnetic potentials ($\mathsf{\phi }$) in both the air space and the PMs are governed by Laplace’s equation when the rectangular coordinate system is adopted [21].

$$\frac{{\partial }^2\mathsf{\phi }}{\partial x^2}+\frac{{\partial }^2\mathsf{\phi }}{\partial y^2}=0$$

(5)

The components of magnetic strength are related to $\mathsf{\phi }$ by:

$$H_x=-\frac{\partial \mathsf{\phi }}{\partial x}\begin{array}{c}{H}_y=-\frac{\partial \mathsf{\phi }}{\partial y}\end{array}$$

(6)

The following boundary conditions are applied to solution of Equation (5):

$$\begin{array}{l}\{\begin{array}{l}{H}_{xI}(x,y){|}_{y=L}=0\\ H_{xII}(x,y){|}_{y=0}=0\end{array}\\ \{\begin{array}{l}{B}_{yI}(x,y){|}_{y=h_m}=B_{yII}(x,y){|}_{y=h_m}\\ H_{xI}(x,y){|}_{y=h_m}=H_{xII}(x,y){|}_{y=h_m}\end{array}\end{array}$$

(7)

By applying the previous boundary conditions, the axial and circumferential flux density (B_y and B_x respectively) in region I and II are solved. For region I (the air space):

$$B_{xI}={\mathsf{\mu }}_0{H}_{xI}={\mathsf{\mu }}_0{M}_n{\displaystyle \sum _{n=1,3,5,\mathrm{...}}^{\infty }\frac{\sinh (n\mathsf{\pi }{h}_m/{\mathsf{\tau }}_p)}{\Delta }}\times \sinh \left\{\frac{n\mathsf{\pi }(L-y)}{{\mathsf{\tau }}_p}\right\}\sin (\frac{n\mathsf{\pi }x}{{\tau }_p})$$

(8)

$$B_{yI}={\mathsf{\mu }}_0{H}_{yI}={\mathsf{\mu }}_0{M}_n{\displaystyle \sum _{n=1,3,5,\mathrm{...}}^{\infty }\frac{\sinh (n\mathsf{\pi }{h}_m/{\mathsf{\tau }}_p)}{\Delta }}\times \cosh \left\{\frac{n\mathsf{\pi }(L-y)}{{\mathsf{\tau }}_p}\right\}\cos (\frac{n\mathsf{\pi }x}{{\mathsf{\tau }}_p})$$

(9)

where:

$$\Delta ={\mathsf{\mu }}_r\cosh (\frac{n\mathsf{\pi }{h}_m}{{\mathsf{\tau }}_p})\sin \left\{\frac{n\mathsf{\pi }(L-h_m)}{{\mathsf{\tau }}_p}\right\}+\cos \left\{\frac{n\mathsf{\pi }(L-h_m)}{{\mathsf{\tau }}_p}\right\}\sinh (\frac{n\mathsf{\pi }{h}_m}{{\mathsf{\tau }}_p})$$

(10)

3.2. Model of Armature Reaction Current

Taking into account the armature winding current effect, the fractional slot windings are regarded as thin wires. The current density J can be obtained by multiplying the value of current sheet distribution of each phase shown in Figure 3, and dividing by the thickness of slot opening t_so as follows:

$$J=(J_A\cdot I_A\sin (\mathsf{\omega }t)+J_B\cdot I_B\sin (\mathsf{\omega }t-\frac{2\mathsf{\pi }}{3})+J_C\cdot I_C\sin (\mathsf{\omega }t+\frac{2\mathsf{\pi }}{3}))/t_{so}$$

(11)

The current sheet spatial distribution of slots in the stator in per unit (J_A, J_B and J_C) of each phase is plotted in Figure 3 [21].

Afterwards, the current density J is expressed by Fourier series to obtain the spatial distribution of the current, which could be derived as:

$$J={\displaystyle \sum _{n=1,2,3}^{\infty }(A_n\cos (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p})+B_n\sin (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p}))}$$

(12)

where A_n and B_n are the coefficients obtained from Fourier series.

As shown in Figure 4, the scalar magnetic potentials $\mathsf{\phi }$ due to the current also can be described by Laplace’s equation. To solve the equation, the following boundary conditions are applied:

$$\begin{array}{l}\{\begin{array}{l}{H}_{xIII}(x,y){|}_{y=0}=0\\ H_{xIV}(x,y){|}_{y=L}=0\end{array}\\ \{\begin{array}{l}{H}_{xIII}(x,y){|}_{y=L^{\prime }}-H_{xIII}(x,y){|}_{y=L^{\prime }}=J\\ H_{yIII}(x,y){|}_{y=L^{\prime }}=H_{yIV}(x,y){|}_{y=L^{\prime }}\hfill \end{array}\end{array}$$

(13)

where L′ is the axial position of a typical armature winding current sheet as shown in Figure 4.

The PMs are physically located in the lower region III and the armature current sheet is very close to the stator back plate. The result in region III will be:

$$B_{xIII}=-{\displaystyle \sum _{n=1,3,5,\mathrm{...}}^{\infty }{\Delta }^{\prime }\sinh (\frac{n\mathsf{\pi }y}{p{\mathsf{\tau }}_p})(A_n\cos (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p})}+B_n\sin (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p}))$$

(14)

$$B_{yIII}=-{\displaystyle \sum _{n=1,3,5,\mathrm{...}}^{\infty }{\Delta }^{\prime }\cosh (\frac{n\mathsf{\pi }y}{p{\mathsf{\tau }}_p})(A_n\sin (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p})}-B_n\cos (\frac{n\mathsf{\pi }x}{p{\mathsf{\tau }}_p}))$$

(15)

where:

$${\Delta }^{\prime }=\frac{{\mathsf{\mu }}_0\cosh ((n\mathsf{\pi }/p{\mathsf{\tau }}_p)\cdot (L-L^{\prime }))}{\sinh (n\mathsf{\pi }\cdot L/p{\mathsf{\tau }}_p)}$$

(16)

3.3. Model of Stator Slotting

The effect of stator slotting can be linked to the previous flux density by defining a vector potential (λ) in each slot, the axial and circumferential components of flux density can be written in the form [18]:

$$B_Z=B_T\cdot {\mathsf{\lambda }}^{*}=(B_y+j\cdot B_x)\cdot (Re(\mathsf{\lambda })+j\cdot Im(\mathsf{\lambda }))$$

(17)

To this aim, the CM method is used to calculate the λ by considering the slot effect in electrical machines, but this method has a defect mentioned above. The numerical SC mapping and its Matlab SC Toolbox could draw the real slot shape [22]:

$$\mathsf{\lambda }={\mathsf{\lambda }}_y+j{\mathsf{\lambda }}_x$$

(18)

where λ₀ is the slotless air gap complex permeance in the T-plane, λ_x and λ_y are the circumferential and axial components of the complex relative air gap permeance in the original Z-plane. Furthermore, K, T, Z and S represent the K-plane, T-plane, Z-plane and S-plane, respectively.

It should be noticed that logarithmic mapping is used for the axial flux electric machine to convert the circular geometry in the S-plane to the Z-plane. However, the AFPM machines in the quasi 3D model (Figure 1) can be regarded as linear permanent magnet machines presented by the Z-plane. Hence, the SC mapping only transforms Z-plane to T-plane. However, in order to calculate the relative air gap permeance via Hague’s solution, the complex K plane is introduced in this paper. Thus, the vector potential can be deduced as:

$$\mathsf{\lambda }=\frac{\partial K}{\partial T}\cdot \frac{\partial T}{\partial Z}\cdot \frac{1}{{\mathsf{\lambda }}_0}$$

(19)

The SC transformation in this step maps one canonical domain (e.g., a rectangle, disk, bi-infinite strip, or the upper (lower) half-plane) into the interior (exterior) of the corresponding polygon. This transformation is defined as follows:

$$Z=f(T)=A{\displaystyle \int {\displaystyle \prod _{k=1}^{n-1}{(t-t_k)}^{-\frac{{\mathsf{\alpha }}_t}{\mathsf{\pi }}}\mathrm{d}t+C}}$$

(20)

where A and C are unknown complex constants, n is the number of polygon corners with interior angle α_t, t₁, …, t_k are the points in the canonical domain (in the T-plane) corresponding to the polygon corners.

The analytical solution of Equation (20) is very difficult for geometries with more than three vertices [22], therefore, the SC Toolbox provides a numerical solution by a library of command-line functions which could be seen in [19]. In this paper, the rectangle domain (T-plane) is created by calling f = correctmap (p, alpha), p is the point and alpha is the corresponding angle as shown in Figure 5.

The next step aims to map the interior of annular domain in K-plane to rectangular domain in T-plane. This mapping can be calculated analytically:

$$T(t)=j(\log (K(k))\frac{\Delta x}{2\mathsf{\pi }}+\frac{\Delta y}{2}-j\frac{\Delta x}{2})$$

(21)

With the determinate rectangular domain vertices by calling v_t = evalinv (f, v_z), the length (Δx) and width (Δy) of the T-plane could be calculated.

Finally, Hague solution for magnetic field in an annular domain could be used, and the final equation which maps the field solution from the K-plane to Z-plane is given as Equation (19), where $\partial T/\partial Z$ can be obtained in SC Toolbox by calling evaldiff (f, Z) and:

$$\frac{\partial K}{\partial T}=j\frac{\Delta x}{2\mathsf{\pi }}\frac{1}{T}$$

(22)

Afterwards, the magnetic flux density at slotted air gap field could be calculated by Equation (17).

3.4. Model of Static Eccentricity

SE is a misalignment condition of stator and rotor axis, thus, the rotor axis is deflected from that of the stator which causes a non-uniform air gap. In a healthy state (Figure 6a), the air gap length g₀ is uniform along the circumferential direction. However, with the occurrence of SE (Figure 6b,c), the rotor shaft experiences a deflection and the symmetry of the rotor deviates from that of the stator by an angle β. The air gap varies from small to large around the circumferential of the stator but the relative position does not vary with time. In other words, in the case of SE the position of the minimum air gap length is fixed in space.

The first step in modeling the eccentricity is to evaluate the air gap length in SE condition [23], the SE factor can be defined as:

$$SEF=\frac{r}{g_0}\times 100\%$$

(23)

The air-gap length at mean radius R_mid can be written as follows:

$$g_m(\mathsf{\theta })=g_0(1-SEF\cos (\mathsf{\theta }-{\mathsf{\gamma }}_0))$$

(24)

where θ is the stator position measured from a reference point (γ₀) of the minimal air gap. According to Figure 6b, the air gap length in all positions can be understood as:

$$g(R,\mathsf{\theta })=g_m(\mathsf{\theta })-g^{\prime }(R,\mathsf{\theta })$$

(25)

where g’ is the air gap deviation at the R_mid:

$$g^{\prime }(R,\mathsf{\theta })=(R-R_{mid})\sin \mathsf{\beta }\cos (\mathsf{\theta }-{\mathsf{\gamma }}_0)$$

(26)

Therefore, the air gap length in all positions can be written as follows:

$$g(R,\mathsf{\theta })=g_0(1-\frac{R}{R_{mid}}SEF\cdot \cos (\mathsf{\theta }-{\mathsf{\gamma }}_0))$$

(27)

Because of the high ratio of machine diameter to length in the AFPMs, 40% SEF means that the maximum declined distance of rotor plate is 1.2 mm in this case, which is already a significant value compared to the air gap length of 3 mm.

In this paper, as quasi-2D model processed in Figure 1, the minimal air gap happens at γ₀ = 0 rad and the maximum air gap at θ= π rad. Finally, the air gap length can be deduced as follows:

$$g(R,\mathsf{\theta })=g_0(1-\frac{R}{R_{mid}}SEF\cdot \cos (\frac{\mathsf{\pi }x}{p{\mathsf{\tau }}_p}))$$

(28)

3.5. Bilinear Mapping

Bilinear mapping [24] is a particular mapping method which is used for the treatment of the boundaries which are circular but eccentric or intersecting. The purpose of this process is to map the pair of eccentric circles (E-plane) into a pair of concentric circles (S-plane). After mapping E-plane to the S-plane, the following mapping could be done by previous study. The general form of this CM can be defined as:

$$E(s)=\frac{1}{s+C}-K$$

(29)

$$H_S=H_E{(\frac{\partial E}{\partial S})}^{*}$$

(30)

where parameters C and K are determined by the rotor and stator radii and the magnitude of eccentricity.

Markovic used this method to calculate a non-salient motor under SE. Bilinear maps preserve the scalar magnetic potential but they cannot conserve the vector field potential (Ω_v) because flux density in concentric plane is coupled with the coordinate system. Since vector field potential calculation is complex and it is very difficult to estimate the scale factor required to map magnetic field into desired eccentric domain, in this paper, the following form is used to convert the air gap flux density of concentric AFPM motor in Z-plane to the eccentric model in E-plane [25]:

$$B_E=B_Z\cdot {\mathsf{\lambda }}_e$$

(31)

where λ_e is the relative permeance function due to the effect of SE and can be defined as:

$${\mathsf{\lambda }}_e(r,\mathsf{\theta })=B_{ecc}(r,\mathsf{\theta })/B_{no-ecc}(r,\mathsf{\theta })$$

(32)

The flux density in the slotless air gap with eccentricity (B_ecc) and without eccentricity (B_no-ecc) can be estimated by modifying the air gap length in Equation (9).

The calculation processes can be seen in Figure 7.

4. Simulation Verification

To validate the analytical model, an AFPMM with nine slots and six poles has been studied. The basic parameters are shown in Table 1. The air gap flux distribution is one of the most important characteristics in electrical machines because any change of this feature would affect the performance of motor. In this section, as the first step of investigation, the air gap flux density under healthy condition is presented. The results of the analytical model are compared with the results calculated by FE model. As the second step, the flux density under SE at the same conditions is also compared with the FE results. Moreover, the effect of eccentricity on air gap flux of AFPM will be discussed in the final step.

4.1. Air Gap Flux under Healthy Condition

The axial and circumferential components of complex relative permeance in the middle of air gap are shown in Figure 8. It shows that, under healthy condition, the complex relative air gap permeance is periodic for each slot pitch.

With Equation (17), the axial and circumferential components of the flux density in the slotted air gap are given by:

$$B_{Zy}=B_y{\mathsf{\lambda }}_y+B_x{\mathsf{\lambda }}_x$$

(33)

$$B_{Zx}=B_x{\mathsf{\lambda }}_y-B_y{\mathsf{\lambda }}_x$$

(34)

The axial and circumferential flux density at the center of air gap under no-load condition are shown in Figure 9. Figure 10 shows the flux density at the rated current (10 A), 20 krpm. The figures show that the analytical model is in a very good agreement with the FE model.

4.2. Air Gap Flux Density under SE Condition

In order to verify the SE model, the rotor of AFPM is assumed with 40% eccentricity. It is shown in previous section that the slotless flux density caused by PMs under SE and healthy condition could be solved analytically, thus, the relative permeance function of SE can be determined. The FE model is represented by two rotor disc to ensure the real flux distribution. The 3D FE model and its grid could be seen in Figure 11. The rotor axis is set as shown in Figure 6a.

The axial and circumferential components of one side rotor disk with 40% eccentricity under no-load and rated load condition conditions are shown in Figure 12 and Figure 13, respectively. It can be seen that the axial and circumferential components of the flux density between the FE model and analytical method do not agree well due to the interaction between the circumferential and axial flux density, but the error still remains in a reasonable range. Overall, the proposed method in this paper is an effective way to investigate the eccentric influence of AFPM.

5. Results and Discussion

In this section, some important machine characteristics, e.g., back-EMF and torque, will be presented and compared. Since these characteristics under different SE conditions only differ in magnitude, the 40% eccentricity is chosen to carry out the analysis.

As shown in Figure 14, the flux passes from one rotor disc and enters the other rotor disc. Therefore, the coils are assumed to be symmetrical with respect to the stator teeth axis. Hence, both rotors contribute with the same amount of magnetic flux (ϕ_Lk,ϕ_Rk). Therefore, the flux linkage of the kth tooth coil with N_c turns is simply given by:

$${\mathsf{\psi }}_k=N_c({\mathsf{\phi }}_{Lk}+{\mathsf{\phi }}_{Rk})/2$$

(35)

5.1. Machine Main Characteristics under No Load Condition

The no-load air gap total flux density under healthy condition compared with 40% eccentricity condition, the analytical results are shown in Figure 15.

Comparing with the healthy conditions, it can be seen that the flux density near minimum air gap has increased but the flux density near maximum air gap has decreased. This is owing to the partly increased yet partly decreased air gap length around circumferential direction.

The back EMF and cogging torque are also important for PM machines. For this purpose, the rotor movement should be modeled to calculate these parameters. The back EMF waveform of the AFPM can be calculated from the axial component of the no-load flux density distribution (B_y) with the knowledge of the armature winding distribution. The voltage induced in a single coil in each layers of quasi-3D method can be calculated as:

$$E_{c,i}=-N_c\frac{\mathrm{d}{\mathsf{\psi }}_{c,i}}{\mathrm{d}t}$$

(36)

where Ψ_c,i is the flux linkage of the coil in ith layer and is equal to the integral of the air gap flux density distribution across one coil pitch, as given by:

$${\mathsf{\psi }}_{c,i}=\frac{R_{ave,i}(R_o-R_i)}{n_s}{\displaystyle \underset{{\mathsf{\theta }}_0}{\overset{{\mathsf{\theta }}_0+{\mathsf{\theta }}_c}{\int }}{B}_{y,i}\mathrm{d}\mathsf{\theta }}$$

(37)

where θ₀ is the coil starting side angle from origin and θ_c is the angle of coil pitch, and B_y,i is the axial component of flux density in ith layer.

The back-EMF per phase is calculated in all coils with the phase winding connected in parallel. For phase A, according to the coils arrangement shown in Figure 3, the back-EMF in ith layer E_A,i is given by:

$$E_{A,i}=-N_c\frac{R_{ave,i}(R_o-R_i)}{n_s}\frac{\mathrm{d}{\mathsf{\psi }}_{c,i}}{\mathrm{d}t}$$

(38)

And based on the quasi-3-D consumption, the back EMF for the whole machine is calculated as:

$$E_A={\displaystyle \sum _i^{n_s}{E}_{A,i}}$$

(39)

The back-EMF of the prototype AFPM with and without SE condition obtained from both FE and SC methods are compared in Figure 16.

It can be seen that in healthy condition, the induced voltage of three coils have the same shape and amplitude. However, monitoring back EMF is not an effective approach to detect eccentric fault in double rotor AFPMs since the back EMF is not affected by eccentricity sufficiently in SE condition. Eccentricity increases the air gap, but on the other hand, the air gap length of corresponding position of the other rotor is decreased. These two effects seem to compensate each other, given the same results.

The back EMFs and its fast Fourier transform (FFT) under SE condition varies from 0% to 65% is shown in Figure 17 and Figure 18 respectively. When the SE factor is increasing, the amplitude of back EMFs is slightly increasing. Moreover, the SE has a strong impact on the first order of FFT.

In the case of SE, the following differential flux linkages act on the phase loops [26], in this paper, phase A is disconnected the parallel path (shown in Figure 16) and:

$${\mathsf{\psi }}_{dl}={\mathsf{\psi }}_{lI}-{\mathsf{\psi }}_{lII}$$

(40)

where ${\mathsf{\psi }}_{lI}$ and ${\mathsf{\psi }}_{lII}$ are magnetic flux linkages of two disconnected parallel paths of phase A.

By time deriving, the no-load loop EMF follows:

$$E_{dl}=\mathrm{d}{\mathsf{\psi }}_{dl}/\mathrm{d}t$$

(41)

Thus, monitoring loop EMF could be a practical approach to detect defects of AFPMs. Figure 19 shows the distribution of phase A loop EMF during one complete rotor revolution at rated speed. It can be seen that the loop EMF fluctuates at each position.

The cogging torque is another important parameter for PM machines which occurs in slotted AFPM machine due to the interaction between the PM magnetic field and the varying air gap permeance.

By assuming negligible saturation, the cogging torque is independent of armature currents and is created only by the PM magnetic field. The period of the cogging torque waveform is calculated as follows:

Cogging torque period = 2π/LCM (Q_s, 2p)

(42)

where LCM (Q_s, 2p) is the least common multiple of the stator slot number and the rotor pole number.

The cogging torque is calculated by the changing rate of total air gap co-energy including the region of PMs [27]:

$$T_{cog}=\frac{\partial W}{\partial \mathsf{\theta }}=\frac{\partial }{\partial \mathsf{\theta }}\left({\displaystyle \underset{\tau }{\iiint }{B}_{t-total}^2\mathrm{d}\mathsf{\tau }}\right)$$

(43)

It should be noted that the AFPM machine in this paper has two rotor disks, which means that both rotors contribute to the cogging torque. Thus, the cogging torque of each disk should be superimposed. The cogging torque waveforms obtained under healthy and SE condition are compared in Figure 20. Comparing with the results obtained from FE model, we can see that the results obtained by analytical model do not agree well, this is mainly due to the subtraction in Equation (40), but the trend of analytical model is still in accordance with that of FE model.

In order to validate the results, 3D FEM calculation is also depicted, which shows a good agreement with analytical simulations. From the simulation results and comparisons, the cogging torque of each disk has the same shape and amplitude, but under SE condition, the cogging torque is increased.

5.2. Machine Main Characteristics under on Load Condition

The output torque is comprised of EM and cogging torque. Neglecting saturation, currents in each phase are updated according to rotor position. The electromagnetic torque for each layer could be calculated by integrating Maxwell stress tensor near the center of the air gap:

$$T_{em,i}=\frac{R_{ave,i}^2(R_o-R_i)}{{\mathsf{\mu }}_0{n}_s}{\displaystyle \underset{0}{\overset{2\mathsf{\pi }}{\int }}{B}_{\mathrm{arm}x,i}{B}_{army,i}d\mathsf{\theta }}$$

(44)

The electromagnetic torque of machine is obtained as:

$$T_{em}={\displaystyle \sum _{i=1}^{n_s}{T}_{em,i}}$$

(45)

Also, the average electromagnetic torque can be simply calculated from

$$T_{em}=\frac{3}{\sqrt{2}}p{N}_c{k}_w{\mathsf{\varphi }}_{f}I$$

(46)

where k_w is the winding factor,φ_f is the magnetic flux excited by the permanent magnet per pole and I is the stator current. From pervious study, the magnetic flux of each phase can be compensated by the right and left magnetic flux, thus, SE does not have significant effects on the EM torque. This also can be proved in Equation (46).

As shown in Figure 21, the results of the proposed method are in good agreement with the simulation results of 3D-FE model.

6. Experimental Results

The double rotor AFPM in Table 1 is the prototype machine investigated in this paper, as shown in Figure 22a,b. The experimental set-up and devices are presented in Figure 22c. The AFPM motor is driven by a frequency converter shown in Figure 22c to test its motor performance, and the no-load parameters.

The back EMF is measured by a driving motor via a belt and the drive speed is 300 rpm, the measured back EMF is shown in Figure 23. The analytical model provides good correspondence with the experimental results shown in Figure 16.

Table 2. Comparison of Back EMF with Experiment Data.

Item	FE Model	Analytical Model	Experiment
Back EMF Coefficient	5.30	5.30	5.28
Error in Back-EMF	0.37%	0.37%

shows the back EMF coefficient obtained by analytical model, FE model and experiment, respectively. However, there is a slight error between the analytical results and the measured value which is mainly caused by the assembly process.

In terms of the computational time, the analytical model requires only 206 s to obtain the cogging torque and electromagnetic torque with one mechanical cycle for SE condition. The 3D FE model, on the other hand, requires almost 27 h to compute the performance in one electrical cycle. Although the analytical method shows a slight error compared to the experimental result, it is still acceptable and thus, can be regarded as a meaningful approach which could save time and achieve a satisfying result. Moreover, the analytical model could also be used for further calculation, e.g., optimization.

7. Conclusions

In this paper an approach for modeling double rotor AFPM with SE has been presented. To take into account the intrinsic 3D nature of the AFPM machine, a quasi-3D method is used. The magnetic field, cogging torque and output torque of AFPM machine under healthy and 40% eccentricity conditions are quickly calculated with great accuracy compared to the 3D FE model and experimental results. In addition, it can be found that the eccentricity plays no effects on the back EMF and the electromagnetic torque. However, the cogging torque increases with SE. In addition, this paper makes a remarkable contribution to the computational time savings while maintaining as much the high accuracy the FE model can provide. Finally, an experimental validation has been carried out to verify the accuracy of the simulation results obtained by the proposed method.