A Simple Analytical Model of Static Eccentricity for PM Brushless Motors and Validation through FEM Analysis

: The paper firstly summarizes a simple analytical model of the air gap flux-density distribution for isotropic permanent magnet (PM) synchronous machines, in the presence of static eccentricity. The model was proposed by the authors in a previous paper and is based on an efficacious analytical expression of the variable length of air gap magnetic field lines which occur in eccentric brushless machines with surface-mounted permanent magnets. The approximate expression of the air gap field makes it possible to achieve a mathematical model with concentrated parameters close to that of a PM machine without eccentricity. The expression of the armature voltages and electromagnetic torque are found, also with reference to steady-state operating conditions at fixed rotor speed and impressed currents. The differences introduced by the considered type of eccentricity are evaluated and highlighted especially with reference to the air gap inductance and to waveforms and frequency spectra of voltages and shaft torque. Numerical results in a case-study of an 8-pole, 110 kW PM motor are compared to those obtained by using finite element analysis.


Introduction
Eccentricity is one of the widely diffused unexpected working conditions in electrical machines. When an eccentricity anomaly occurs in a machine, the length of the air gap between rotor and stator has a non-uniform distribution. As is well known, if the rotation axis of the rotor is fixed in time and parallel to that of the stator, the eccentricity is called "static"; otherwise it is named "dynamic". Many technical papers in literature deal with the study of the effects of eccentricity in induction machines' performance [1,2].
In the last two decades, brushless drives have largely been used due to the necessity to increase the power density and the overall efficiency of electric drives. Also, for this type of electric machine, eccentricity is one of the most common faults [3,4]. Regarding the study of the behavior of brushless electrical machines in the presence of eccentricity, some researchers are oriented to identifying anomalies during operation using mechanical vibrations and related noise [5]. Other authors proposed novel approaches based on stator current signatures [6].
The study of eccentricity effects using finite element analysis (FEA) is proposed in [7,8]; these techniques are characterized by accurate results, but the high computational cost is not always compatible for use in a real-time diagnostic process.
An analytical approach appears more useful for the diagnosis of eccentricity; different methodologies are applied in the literature [9][10][11][12] in order to develop an accurate mathematical model suitable for correctly taking eccentricity into account.
With the aim to reduce the complexity of the models, it is possible to introduce some simplifying hypotheses on the length of magnetic flux-density in the air gap and on the magnetization law of the permanent magnets [13,14]. Using these assumptions, in the paper the authors deepen a simplified mathematical model, already introduced in [15], capable of simulating the effects of static eccentricity. Referring to PM motors with surface-mounted permanent magnets, the paper is focused on the comparison between waveforms and spectra of some relevant quantities evaluated with and without static eccentricity. The results of the model are verified by means of finite element method (FEM) analysis in terms of no-load electromotive force (e.m.f.) and armature voltages and electromagnetic torque in an assigned steady state condition with impressed armature currents.

Type of Motor Under Analysis
The analysis of the effects of static eccentricity is carried out on isotropic radial-flux brushless motors, with surface-mounted magnets (Figure 1), which cover about 75-90% of the pole pitch. The symmetrical three-phase stator winding is of distributed type (more than one slot per pole and per phase); the coils of each phase are all considered in series.

Basics of the Analytical Approach to the Study of Static Eccentricity
The stator and rotor iron surfaces facing the air gap are assumed to be cylindrical. Under the usual simplifying hypotheses in motor modeling, in ordinary (non eccentric) conditions, the magnetic field lines in the air gap are radial and with constant length In the presence of static eccentricity, the stator and rotor axes are no longer coincident, but parallel with a distance d ( Figure 2b) and remain fixed over time during rotor rotation. As is known, eccentricity is the quantity usually defined as Considering a generic point P on the internal surface of the stator, from Figure 2b we can easily deduce that the segments PP' (along the stator radius OsP) and PP" (along the rotor radius OrP) are practically coincident. So, according to [15], we can consider the segments of PP' type as approximated paths of the flux-density lines in the air gap, also called the magnetic air gap. Obviously, they are perpendicular only to the stator surface.
The length δ of the air gap is then variable along a stator angular abscissa α, with a minimum where the size of the air gap is magnified for greater clarity, the segment PP' can be calculated as: where αγ is the angular position of γ, as shown in Figure 2b.
With the hypothesis that d 2 << Rr 2 , Equation (1) becomes: Using the above definition of eccentricity, the variable length of the air gap can be written as: In ordinary conditions the air gap magnetic permeance is constant and equal to 0 when eccentricity is present, the permeance is a function of α and does not depend on the angular position of the rotor: The expansion in Fourier series of λ(α) is: as shown in [15], the permeance in p.u. becomes: from which: Usually few terms of the sum (5) are enough for an acceptable approximation. In order to highlight how the air gap permeance is modified by the static eccentricity compared to the ordinary condition, Figure 3 shows a qualitative behavior of λ(α)/g0 for assigned values of ε and αλ (ε = 0.25, αλ = π/2). In the figure the higher value of λ with respect to g0 is also represented.

Armature Flux-Density in the Air Gap
As already explained in [15], by taking into account the permeance expression (6) and neglecting iron losses and magnetomotive force (m.m.f.) drops in all iron paths, the instantaneous distribution Bs(α, t) along α of the flux-density generated in the air gap by a stator 3-ph system of currents can be expressed by: with: • p pole-pair number; • a wye-connected 3-ph winding, whose p coils of each phase are all in series.
Moreover, in (6) the term , ( ), B s n t which is dimensionally a flux-density, is given by: where: • n is an index equal to 1 and (6ν±1) with ν ∈(1,∞); i. e. n = 1, 5, 7, 11, 13, ....; • N is the number of turns in series per phase; •  n is the winding-factor for the n th spatial harmonic; • qn is a coefficient equal to: • i (t) is the instantaneous currents space-vector equal to: is the instantaneous current value in the k th phase.
From Equation (9) we can deduce that the various terms depend on the quantities p n h ± and, for some of them, only the λ coefficients of order pn exist. As a sample, Figure 4 shows the behavior along α of the flux-density Bs(α,t) generated by one stator phase current (e.g., i1(t)) at its rated value. The case of a PM-motor with 8 poles and 48 armature slots is considered: magnitude and position of the eccentricity are respectively ε = 0.25 and αλ = π/2.
More detailed parameters of the motor are reported in Section III.B, Table 1. The variability of the magnetic permeance can be deduced easily from the upper and lower envelopes of the flux-density curve which replicates the permeance behavior of Figure 3. The quantity BM * in Figure 4 is the maximum flux-density value of one stator phase in the healthy non eccentric condition.

Rotor Flux-Density in the Air Gap
We assume that the remanence in rotor magnets is a generic BRe(β) symmetrical and periodic function of the rotor angular abscissa β, which has clockwise orientation and origin 0 β = on a polar axis. The air gap field-lines generated by the rotor magnets are assumed to coincide with those of the stator as in Figure 2b. Consequently, the rotor magnets work on circuits with variable permeance not only along β for a fixed position ϑ, but also over time for each fixed point at the abscissa β. Obviously, at constant speed ωr, the ϑ position of the rotor is a linear function of time: where ϑ0 is the angle between the α and β axes at instant 0. According to the hypotheses introduced in [15], the flux-density in the air gap is assumed constant along each field-line in Figure 2b. Since in rotor coordinates the length of these lines is δ (α) = δ (β +ϑ), denoting with lm the radial length of the magnets, the Br(β,ϑ) flux-density generated by the rotor magnets in the air gap is linked to the remanence BRe(β) through the expression: where Cϑ is constant for each ϑ rotor position; it can be evaluated imposing the magnetic field to be solenoidal.
Considering a generic magnetization function BRe(β), whose expansion in Fourier series is: and referring to the expression of the air gap permeance (5) where the coefficients pu h λ are given by (7) At first, we can consider the ideal case of magnets which generate perfectly sinusoidal air gap flux-density in the absence of eccentricity, i.e., only the coefficient Br,1(ϑ) ≠ 0 in (15). With reference to an 8-pole motor (see Table 1 for its rated values and main data) in which a static eccentricity is arisen with ε = 0.25, Figure 5 shows the flux-density curves for two arbitrary ϑ positions of the rotor, assuming αλ = π/2.  As a second example, we can refer to the same motor of Table 1, but with different magnetization compared to the previous case. Figure 6 shows a pole pitch τp of the rotor; the magnets cover the τa polar arc and can be divided into some parts (three parts in the example); each of the latter is uniformly magnetized along radial directions. The value of the residual induction uniformly obtained in the magnets is assumed to be equal to 1.2 T at the temperature of 120 °C.  Figure 7 shows the field lines generated by rotor magnets in one pole of the considered motor of Table 1 in no-load operating conditions, neglecting the presence of the armature slots and assuming the iron permeability to be constant. The behavior along α of the flux-density Br(α,ϑ) is drawn in Figure 8, for an assigned ϑ position (ϑ = π/8 rad); the diagram a) refers to the healthy case without eccentricity, while the diagram b) shows the variations introduced by static eccentricity with ε = 0.25.
The diagrams of Figure 8 show a good match between the curves obtained theoretically or using FEM, in both cases without (a) and with (b) eccentricity. As can be seen from Figure 5 and Figure 8, the static eccentricity significantly influences the air gap distribution of Br(α,ϑ) flux-density. In both diagrams the curves obtained by FEM analysis (Appendix A) are drawn with solid blue lines, while brown dashed lines represent the results of Equations (15) and (16), with reference to an uniform remanence BRe = 1.2 T in the magnets, that cover polar arcs of 88% of the pole pitch (see Table  1). As it is possible to note, the FEM analysis are carried out for a 2-D geometry, usually adopted for a radial flux motor. The 2-D CAD software neglects the behavior of end-windings and the edge effects of stator and rotor, but the end-windings are considered no-active magnetic parts of the motor (they affect only the leakage inductance value and the resistance calculation, and can be easily compensated introducing analytical corrections).

Resultant Flux-Density in the Air Gap
The total air gap flux-density Bδ (α,ϑ) is obtained by adding the separate contributions of stator (6) and rotor (15), also considering that ϑ is the time-function (12); i.e.,: where the coefficients are deduced from (10) and (16): and where the approximated expression (7) of the p.u. permeance coefficients λh pu is used.

Induced Voltages and Armature Equations
With reference to the instantaneous values of the different quantities involved and neglecting the elasticity of the shaft, the mathematical model of the PM motor considered is composed by the set of equations: where the different symbols represent: • me-electromagnetic torque; The aforementioned space vectors are equal to: where vk(t) and ik(t) are instantaneous voltage and current values of the k th armature phase; all the coils of each phase winding are considered in series. Using this hypothesis, the back e.m.f. erot,k (t) induced in the k th stator phase by the rotor magnets is the sum of the back e.m.f. erot,k,ν (t) induced in the p coils of the considered phase: Obviously, if there is eccentricity, the different rotor fluxes φk,γ (t) (with γ = 1, ..., p) linked with the various γ th coils of the k th phase are not equal, because of their different positions with respect to the axis λ in Figure 2b.

Electromagnetic Torque
The instantaneous electromagnetic torque developed by the motor can be evaluated by means of the well-known expression:  (11) and (17) the non-null terms are only those with h p n pν = + or h pν = ; therefore, the previous equation can be written as: where: since from (7) (7) we derive: Synthesizing the above considerations, the torque expression (29) becomes: The air gap flux-density , ( ) with , ( ) r n t B expressed by (16). Moreover, the terms with n = ν are separated from the others; we obtain:  with , ( ) s n t B expressed by (10), i.e.,: where Lm,1 is the air gap mutual inductance relating to the 1st spatial harmonic of the armature field in presence of eccentricity: The resulting instantaneous torque is then the sum of (40) and (41).

Numerical Investigation
In this section, some numerical results of the model introduced above are compared with those obtained using analysis with finite elements methods. The results refer to the PM motor of Table 1, under steady state operating conditions at constant ωr speed. In these conditions, the erot,k(t) back e.m.f.
The angle ϑ0 is the position of the rotor at instant t = 0 (12), while Φr,n is the amplitude of the n th harmonic of the rotor flux linked with a stator phase (26). With reference to a rotor speed of 1500 rpm, Figure 9 shows the no-load induced voltage in an armature phase (diagram a) and the related frequency spectrum (b). The f1 fundamental frequency of the induced voltage is 100 Hz. In the spectrum of Figure 9b the harmonic amplitudes are expressed in p.u. of the fundamental E1,FEM obtained in the FEM analysis. Blue lines refer to the results of FEM analysis obtained with the rotor flux-density in Figure 8; brown lines are used for the theoretical curves. An acceptable approximation can be observed, especially for the most significant harmonics.  The eccentricity does not introduce any harmonic of different order than those of the symmetric case; moreover the amplitudes of fundamental and other harmonics are a little higher (a few percent) when the eccentricity is present.
A load condition is examined with impressed armature currents. The latter are assumed to be symmetrical and sinusoidal with angular frequency In (46) En is given by (44) and In = 0 for n ≠ 1, by hypothesis. With reference to the motor of Table 1 in a load condition with rated armature currents, Figure  11 shows the waveforms of the three armature phase voltages; the diagram (a) refers to the voltages evaluated using FEM analysis, while the voltages draw in the diagram (b) are evaluated by means of the above presented model. The three voltages remain symmetrical, as predicted by the model; the two sets of curves (a) and (b) are quite similar. The same occurs for the corresponding frequency spectra shown in Figure 12.   The electromagnetic torque at steady state can be derived from (39) where: • e j I ϕ = I is the phasor of the impressed sinusoidal armature current at angular frequency The term (a) in (46) is extracted from (40) and represents the contribution to the average torque due to the interaction between stator and rotor fields. The term (b) is extracted by (41) and depends only by the armature currents. Compared to the healthy case, the eccentricity slightly changes the term Φr,1 and all the various ,1 ν σ ; however, the contribution of the latter terms is very low. The air gap magnetic anisotropy introduced by eccentricity gives rise to the term (b), whose contribution to the average torque is also quite modest and can generally be neglected.
By analyzing Equation (40) at steady state, some sinusoidal components of the torque can be deduced at angular frequency 6ω, 12ω, ... They depend on the interaction between the spatial harmonics of the linear current density and the spatial harmonics of the rotor field in the air-gap. Their analytical expressions are here omitted for simplicity, also because the magnitude of these components is not significantly different from those in the absence of anisotropy. Moreover, Equation (41) at steady state yields additional sinusoidal torque components at angular frequency 2ω due to the magnetic anisotropy caused by eccentricity. However, the amplitude of these components is very small and can be neglected.
With reference to the rated steady state operating condition, the electromagnetic torque is plotted as a function of the time in Figure 13a, where the waveform evaluated with the model equations (brown line) is compared with that obtained by the FEM analysis (blue line); the figure shows also the instantaneous torque in non-eccentric case (gray line). A magnification of a section allows to better highlight the differences between the three curves. The spectrum of the torque is shown in Figure 13b; the amplitude values are expressed in p.u. of the average value Me,med evaluated by FEM in the healthy non-eccentric case. Only harmonics of order 6n are visible, even if all of them are very small with respect to the torque average value. Also in this case there is a good correspondence between the results of the model and those of FEM analysis. In broad terms, we can notice that eccentricity does not modify the shape of the instantaneous torque but slightly increases its average value.

Conclusions
By some geometrical considerations and a set of simplifying hypotheses, an analytical lumpedparameter model was derived for a surface-magnets PM motor affected by static eccentricity. The equations presented apply for a distributed 3-phase winding with only series connections between coils of a given phase. The analytical model was validated successfully through FEM analysis for a case-study in which sinusoidal currents were assumed to be impressed in the motor, as common for AC brushless motors fed with current-controlled voltage source inverters. The trends of the results show the goodness of the proposed analytical model: in fact, the difference in the magnitude of the The numerical validation confirms the analytical model as an adequate mean to gain further insight into the effects of static eccentricity on the operation of the PM motors. These effects can be broadly summarized as follows.
Static eccentricity translates into magnetic anisotropy and thus introduces infinite additional harmonic components in the airgap flux-density distributions generated by stator currents and permanent magnets; the resulting spatial distributions are modulated and lose their periodicity along a double pole pitch.
Due to the assumed full-pitch armature winding and the series connection between phase coils, these additional spatial harmonics do not introduce in the linked fluxes and back e.m.f.s different harmonics from those already observable during healthy operation. Symmetric behavior in the three phases is preserved as well.
The eccentricity affects the values of some parameters in the model. In particular, the amplitude of the first harmonic of rotor flux experiences a little increase, resulting in a correspondent small increase in the mean value of the electromagnetic torque with respect to the healthy case. The mutual inductance is slightly increased as well. However, both variations appear rather modest.
These results suggest that, for the considered layout of the machine, static eccentricity cannot be easily detected through time or frequency inspection of the main electromagnetic variables. The highlighted variation in the machine parameters appears as the only potential indicator of the considered fault, although weak and hardly detectable in real applications.

Appendix A-Finite Element Analysis
The verification of the goodness of the analytical model is carried out through the comparison with finite element analysis. The analyses are performed using the software FEMM 4.2 [16]. The following appendix summarizes the modalities of calculation of back-e.m.f., phase armature voltage, electromagnetic torque and synchronous mutual inductance: • Calculation of back-e.m.f.
The back-e.m.f. is calculated using the Faraday's law. The methodology needs the determination of the linkage flux with a single phase. For a 2-D problem, the linkage flux of a single-phase coil with massive conductor is calculated as [17,18]: where Sc is the conductor area, Az is the z-component of the potential vector of magnetic flux density, L is the stack length of the motor. The integral in (A1) are computed both on left and right side of the phase coils. The numerical derivative of (A1) gives the back-emf. •

Phase armature voltage
The phase armature voltage [17,18] is calculated adding the voltage drop on the single-phase resistance to the back-e.m.f. The resistance is calculated through the evaluation of Joule losses. •

Electromagnetic torque calculation
The electromagnetic torque is calculated using the Maxwell stress tensor [17,18]. Considering a 2-D geometry, the expression is: where L is the magnetic stack length, r is the radius of an opportune contour lr located in the airgap and around the rotor, while Bn and Bt are the normal and tangential components of magnetic flux density on the contour. •

Synchronous mutual inductance
The value of synchronous mutual inductance Lm is obtained by means of the calculation of electromagnetic energy Egap in the airgap [17,18]. Eliminating the permanent magnet magnetization, the Egap is calculated imposing a quadrature axis current of amplitude Imax (it was verified that with a d-axis current, the value of mutual inductance of the considered motor is the same). Therefore, the mutual inductance is equal to: