# Nonlinear Vibration Analysis of Rotor considering Cogging and Harmonic Effects

Authors: Shanle Li; Feng Liu; Hongyan Wang; Haijun Song; Kuilong Yu
Journal: Shock and Vibration
 (2021)

---

## Abstract

This paper aims to investigate nonlinear vibration characteristics of rotor system considering cogging and harmonic effects. Firstly, relative permeance with eccentric was established and then corrected by correction factor caused by the cogging effect. Based on the new formula of relative permeance, the expression of unbalanced magnetic force was obtained, and the coefficient of cogging effect was defined. Motion equations of rotor system were established, and Runge–Kutta method was used to solve the equations. Results showed that errors between finite and analytical results were smaller considering cogging and harmonic effects. When the harmonics were taken into consideration, the vibration of rotor increases sharply. When the cogging and harmonics were taken into consideration simultaneously, the vibration of rotor decreased instead, which means that stator slots have the effect of reducing vibration in rotor system. Rotor vibration was axis symmetry with static eccentricity rather than central symmetry with no eccentricity, and double, four times, and six times supply frequency always existed in the components of main frequency with eccentric.


---

## Body

## 1. Introduction

The existence of eccentric is rather common in electric machinery, and issues with unbalanced magnetic pull (UMP) caused by eccentric have been investigated for many years [1]. Generally speaking, the electromagnetic forces are determined on electronic parameters including the stator windings and permanent magnetic and mechanical parameters such as stator slots and air-gap eccentric between stator and rotor. The air-gap eccentric has been studied in many papers [2–5], and it contains static eccentric, dynamic eccentric, and the combination of both. For an ideal motor, the air gap between the rotor and stator is even, which means that the electromagnetic force is equal to zero. However, air gap is not always uniform due to manufacturing tolerance, incorrect assembly, and so on. Consequently, UMP is generated. The UMP leads the rotor to move towards the stator in the direction of minimum air gap. Meanwhile, the eccentric increases under the motion. The coupling interaction may cause vibrations for the rotor and work failure under some extreme condition [6–10].The expression of the UMP is the foundation for the study of the vibrations and dynamical response of electrical machinery. Many papers have been published on the subject. Among those methods, the finite element method [11–13] and analytical method [14–19] are widely used. Since the finite element method only gives results of the UMP, and it cannot give insight into the production of UMP, it is often used to verify the accuracy of the analytical method. Researchers are paying attention to the analytical method. In the early days, the expressions of UMP are given as linear equations. In Behrend’s paper [20], UMP is assumed to be in proportion to the air-gap length, and the linear expression of UMP is obtained. Later on, Covo [21] improved the linear expression with the saturation taken into consideration. Although the linear analytical model has solved some actual problems, but the results are not reliable when the eccentricity is lager. In 1965, Funke [22] put forward a nonlinear model of the relationship between the UMP and the eccentricity and analyzed the influence of UMP on a synchronous machine. Smith [23] introduced the winding analysis to study the influences of harmonics on UMP. Li [24] used the conformal mapping method to calculate the UMP with eccentricity. Among all the methods to calculate the UMP, air-gap permeance method is the most widely used method. In this method, the air-gap permeance is expressed as Fourier series, and the analytical expression of air-gap flux density is obtained. Then, the expressions of UMP are obtained by using Maxwell tensor method. In Guo’s paper [25], expressions of UMP with different pole number were obtained, and the influence of pole number and rotor angular speed on rotor shaft orbits was investigated. Afterwards, many papers have been published based on the expression of UMP in Guo’s paper. Gustavsson [26] investigated the influence of UMP on the stability of a hydrogenerator. Wu [27] studied the influence of UMP and mass eccentric on the stability of a synchronous generator. Zhang [28] studied the nonlinear dynamic characteristics of a rotor-bearing system with rub-impact for hydraulic generating set under the UMP. Chen [29] obtained the approximate solution to the equation of nonlinear vibration under UMP and discussed the stability of the steady response. Xu [30] investigated the vibration of a generator rotor considering both dynamic and static eccentricities. Xiang [31] analyzed the stiffness characteristics of the rotor system and the stability of equilibrium points in free vibration based on a Jeffcott rotor model. Zhang [32] established an 8-dof model with the ball bearing taken into consideration, and the effects of static eccentricity, rotor offset, and radical clearance of bearings were investigated. Liu [33] analyzed the localized oscillations and nonlocalized oscillations of the rotor-bearing system with the nonlinear restoring force taken into consideration. Chai [34] studied the effects of stator oval deformation, rotor centrifugal distortion on the air-gap magnetic flux field, radial and tangential magnetic force, and electromagnetic torque.Although a lot of papers have been published on the nonlinear vibration of rotor system. Most of the works are focusing on the influence of eccentricity on rotor vibration. There are hardly papers concentrating on cogging and harmonic effects. The purpose of this paper is to investigate the nonlinear vibration considering cogging and harmonic effects. In Section2, the relative permeance with eccentric was established and then was corrected by correction factor of cogging effects. Based on the formula, the expression of unbalanced magnetic force considering cogging and harmonic effects was obtained and compared with finite results under different static eccentrics. In Section 3, motion equations of rotor system were established and solved by Runge–Kutta. In Section 4, eccentric, harmonic, and cogging effects on vibration were investigated in detail. Finally, the conclusions were drawn in Section 5.


## 2. Analysis of Unbalanced Magnetic Pull

### 2.1. Calculation of Relative Permeance with Eccentricity

When the air gap between the stator and the rotor is not uniform, the magnetic flux density becomes asymmetric, and eccentric appears. The cross section of a motor with eccentricity is shown in Figure1. O0O and Oc are the geometry of the stator, the initial geometrical center of the rotor, and geometrical center of the rotor. For simplicity, a coordinate with the origin located at O is established. O0O⟶ is the static eccentricity. OOc⟶ is the dynamic eccentricity. O0Oc⟶ is the composite eccentricity. In this paper, we assume that the coordinate of the static displacement in stator coordinate system is x0,y0, and the coordinate of the dynamic displacement in rotor coordinate system is x,y. β is the eccentricity angel. Rs, Rr are the inner radius of the stator and the outer radius of the rotor.Figure 1
Rotor system and air gap of PMSM with composite eccentricity.
(a)(b)Figure1(a) shows the model of single-disk rotor with hinged ends. The rotor disc is located at the center between the two hinges, which means that the rotor has two degrees of freedom. The coordinate of mass center of rotor in XO0Y coordinate system can be used to represent the position of the rotor disc. Figure 1(b) shows the composite eccentricity of rotor. B and A are the point on the inner surface on the stator and the point on the outer surface on the rotor. The air-gap length between the rotor and stator can be calculated as(1)δα,t=OcB¯−OcA¯.TheOcB¯ and OcA¯ are as follows:(2)OcB¯=Rs2−e2sin2α−β−ecosα−β,OcA¯=Rr,where e is the composite eccentricity. It can be expressed as follows:(3)e=x+x02+y+y02.Sincee is much smaller than Rs, Rs2−e2sin2α−β can be simplified as Rs2, which makes the air gap can be expressed as follows:(4)δα,t=Rs−Rr−ecosα−β=δ0−ecosα−β,where δ0 is the average air gap with no eccentricity, and the position angel can be expressed as follows:(5)β=arcosx+x0e,y+y0>0,2π−arcosx+x0e,y+y0<0.If we take the magnet thickness into consideration, the average air gap with no eccentricity should be added withhm/ur. Then, the air gap can be expressed as(6)δα,t=δ0+hmur−ecosα−β.The relative permeance can be calculated as(7)Λeα,t=μ0δα,t.where μ0 is the air permeability.As shown in Figure2, relative permeance is not constant due to cogging effect. It is more like the shape of the stator tooth, which is more like the actual relative permeance between the stator and rotor. Figure 2(b) shows the relative permeance with static eccentricity. It is nonuniform due to the static eccentricity. And the minimum value is in the static eccentric direction, while the maximum is in the direction away from the rotor shifting.Figure 2
Comparison of relative permeance. (a) Relative permeance with no eccentricity. (b) Relative permeance with 10% static eccentricity relative to the air gap between stator and rotor.
(a)(b)


### 2.2. Correction Factor of Relative Permeance

The relative permeance with stator slots can be obtained by introducing correction factor of stator slots. And it can be expressed as follows:(8)Λα,t=λα,tΛeα,t,where λα,t is the correction factor caused by stator slots.The correction factorλα,t can be calculated as(9)λα,t=1Kc1+∑k=1∞Akcosks1α,where Kc is Carter’s coefficient, and Ak is the harmonic permeance coefficient.Carter’s coefficientKc can be calculated as(10)kC=t1t1−γ1δ,where t1 is the slot pitch, and γ1 can be expressed as follows:(11)γ1=4π0.5κarctan0.5κ−ln1+0.5κ2.According to paper [35], the coefficient in equation (4) can be expressed as follows:(12)Ak=−βκ1k4π0.5+kΓ20.78−2kΓ2sin1.6πkΓ,where(13)βκ=0.51−11+κ2,Γ=b14t1.For permanent magnet synchronous motors,(14)κ≈b14δ+hM/μ0.


### 2.3. Expression of Unbalanced Magnetic Pull

According to the electric machine theory, the magnetomotive force (MMF) for a general motor can be expressed as(15)Fα,t=Fsα,t+Frα,t,where Fsα,t and Frα,t are the magnetomotive force of the stator windings and the permanent magnet. In this paper, the harmonics caused by current harmonics are ignored. The magnetomotive force of stator windings and the permeant magnetic can be expressed as follows:(16)Fsα,t=∑v=6k∓1∞Fsvcosvpα±ωt,k=0,1,2,…,(17)Frα,t=∑u=2n−1∞Frucosupα−uωt,n=1,2,…,where v is the harmonic order of the stator magnetomotive force, u is the harmonic order of the rotor magnetomotive force, and Fsv and Fru are the amplitude of the magnetomotive force of the stator windings and the permanent magnet. And they can be expressed as follows:(18)Fsv=1v1.35N1kdpvI1p,Fru=4πuBrhmμ0sinuαpπ2,where I1 is amplitude of the current of the stator windings, p is the number of pole-pair numbers, N1 is the number of the stator turns in series per phase, kdpv is the winding factor for the v−th harmonic, Br is the magnet remanence, hm is the magnet thickness, αp is the pole-arc/pole-pitch ratio, and μ0is the vacuum permeability.As shown in Figures3 and 4, armature MMF (stator MMF) and permeant MMF (magnetic MMF) are a longer standard trigonometric function due to the harmonic effect. It is the composition of trigonometric functions.Figure 3
Comparison of Armature MMF. (a) Arm MMF at a given time. (b) Armature MMF at a given point in time history.
(a)(b)Figure 4
Comparison of PM MMF. (a) PM MMF at a given time. (b) PM MMF at a given point in time history.
(a)(b)According to the air-gap permeance approach, the magnetic flux density of air gap can be expressed as follows, and the waveform is shown in Figure5. It shows that the air-gap flux density is nonuniform due to cogging and harmonic effects.(19)B=Λα,tfα,t.Figure 5
Comparison of air-gap flux density. (a) Air-gap flux density at a given time. (b) Air-gap flux density at a given point in time history.
(a)(b)Since the tangential magnetic flux density is much smaller than radical magnetic flux density, the tangential magnetic flux density in this paper is ignored. The radical component of Maxwell stress is expressed as(20)σr=Br22μ0.The analytical expression of UMP can be obtained by integrating the Maxwell stress on rotor surface:(21)Fux=RrL∫02πσrcosαdα,Fuy=RrL∫02πσrsinαdα.Taking equations (8), (16), and (17) into (21) and only considering the first several harmonics, the expression of UMP considering harmonics and stator slots when the pole-pair number is 1 can be expressed as(22)Fux=2+∑k=1NAk2kc2F0cosγ+F1cos5γ+F2cos2ωt−γ+F3cos2ωt−3γ+F4cos4ωt−3γ+F5cos4ωt−5γ+F6cos6ωt−γ+F7cos6ωt−5γ+F8cos2ωt+3γ+F9cos2ωt+5γ+F10cos4ωt+γ+F11cos4ωt+3γ+F12cos6ωt+γ,Fuy=2+∑k=1NAk2kc2F0sinγ+F1sin5γ+F2sin2ωt−γ+F3sin2ωt−3γ+F4sin4ωt−3γ+F5sin4ωt−5γ+F6sin6ωt−γ+F7sin6ωt−5γ+F8sin2ωt+3γ+F9sin2ωt+5γ+F10sin4ωt+γ+F11sin4ωt+3γ+F12sin6ωt+γ,where(23)F0=14Fr12+Fr32+Fr52+Fs12+Fs52+2Fr1Fs12Λ0Λ1+Λ1Λ2+Λ2Λ3,F1=2Fs1Fs512Λ2Λ3,F2=14Fr12+Fs12+2Fr1Fr3+2Fr1Fs1+2Fr3Fr5+2Fr3Fs1Λ0Λ1+12Λ1Λ2+12Λ2Λ3,F3=14Fr12+Fs12+2Fr1Fr3+2Fr1Fs1+2Fr3Fr5+2Fr3Fs1Λ0Λ3+12Λ1Λ2,F4=142Fs1Fr3+2Fs1Fr5+2Fr1Fr3+2Fr1Fr5Λ2Λ3+12Λ1Λ2,F5=142Fs1Fr3+2Fs1Fr5+2Fr1Fr3+2Fr1Fr512Λ2Λ3,F6=142Fs5Fr52Λ0Λ1+Λ1Λ2+Λ2Λ3,F7=14Fr32+2Fs1Fr5+2Fr1Fr512Λ2Λ3,F8=142Fs1Fs5+2Fr1Fs5Λ0Λ3+12Λ1Λ2,F9=142Fs5Fr1+2Fs1Fs512Λ2Λ3,F10=142Fs5Fr32Λ0Λ1+Λ1Λ2+Λ2Λ3,F11=142Fs5Fr3Λ0Λ3+12Λ1Λ2,F12=142Fs5Fr52Λ0Λ1+Λ1Λ2+Λ2Λ3.For simplicity, the coefficient of cogging effect is defined, and it can be expressed as(24)λs=∑k=1∞Ak+22kc.


### 2.4. Verification and Discussion

In order to verify the accuracy of the analysis, a finite element model (FEM) of the same parameters is established. Figure6 shows a comparison between the analytical results and finite element analytical (FEA) results of the air-gap flux density. It shows that the analytical results with cogging and harmonic effects taken into consideration are in great consistency with the FEA results. Table 1 shows the results of the analytical method and FEM method. It shows that the harmonics have the effects of increasing the UMP, while the stator slots have the effects of decreasing the UMP. In addition, the errors between the analytical results taking the harmonics and stator slots into consideration and FEM results are the smallest, which means that the analytical method is accurate, and it can be used to analyze the dynamical response of the rotor system.Figure 6
Comparison of air-gap flux density. (a) Air-gap flux density at a given time. (b) Air-gap flux density at a given point in time history.
(a)(b)Table 1
Comparison of unbalanced magnetic pull.
Static eccentricity0%10%20%30%40%50%FEM results03.34 × 1056.92 × 1051.15 × 1061.76 × 1062.65 × 106Analytical results (no harmonics and stator slots)03.04 × 1056.40 × 1051.05 × 1061.60 × 1062.40 × 106Errors0%8.98%7.51%8.70%9.09%9.43%Analytical results (with harmonics)03.83 × 1057.82 × 1051.27 × 1062.05 × 1062.85 × 106Errors0%14.67%13.01%10.43%16.48%7.55%Analytical results (with harmonics and stator slots)03.52 × 1057.31 × 1051.22 × 1061.86 × 1062.78 × 106Errors0%5.39%5.64%6.09%5.68%4.91%


## 2.1. Calculation of Relative Permeance with Eccentricity

When the air gap between the stator and the rotor is not uniform, the magnetic flux density becomes asymmetric, and eccentric appears. The cross section of a motor with eccentricity is shown in Figure1. O0O and Oc are the geometry of the stator, the initial geometrical center of the rotor, and geometrical center of the rotor. For simplicity, a coordinate with the origin located at O is established. O0O⟶ is the static eccentricity. OOc⟶ is the dynamic eccentricity. O0Oc⟶ is the composite eccentricity. In this paper, we assume that the coordinate of the static displacement in stator coordinate system is x0,y0, and the coordinate of the dynamic displacement in rotor coordinate system is x,y. β is the eccentricity angel. Rs, Rr are the inner radius of the stator and the outer radius of the rotor.Figure 1
Rotor system and air gap of PMSM with composite eccentricity.
(a)(b)Figure1(a) shows the model of single-disk rotor with hinged ends. The rotor disc is located at the center between the two hinges, which means that the rotor has two degrees of freedom. The coordinate of mass center of rotor in XO0Y coordinate system can be used to represent the position of the rotor disc. Figure 1(b) shows the composite eccentricity of rotor. B and A are the point on the inner surface on the stator and the point on the outer surface on the rotor. The air-gap length between the rotor and stator can be calculated as(1)δα,t=OcB¯−OcA¯.TheOcB¯ and OcA¯ are as follows:(2)OcB¯=Rs2−e2sin2α−β−ecosα−β,OcA¯=Rr,where e is the composite eccentricity. It can be expressed as follows:(3)e=x+x02+y+y02.Sincee is much smaller than Rs, Rs2−e2sin2α−β can be simplified as Rs2, which makes the air gap can be expressed as follows:(4)δα,t=Rs−Rr−ecosα−β=δ0−ecosα−β,where δ0 is the average air gap with no eccentricity, and the position angel can be expressed as follows:(5)β=arcosx+x0e,y+y0>0,2π−arcosx+x0e,y+y0<0.If we take the magnet thickness into consideration, the average air gap with no eccentricity should be added withhm/ur. Then, the air gap can be expressed as(6)δα,t=δ0+hmur−ecosα−β.The relative permeance can be calculated as(7)Λeα,t=μ0δα,t.where μ0 is the air permeability.As shown in Figure2, relative permeance is not constant due to cogging effect. It is more like the shape of the stator tooth, which is more like the actual relative permeance between the stator and rotor. Figure 2(b) shows the relative permeance with static eccentricity. It is nonuniform due to the static eccentricity. And the minimum value is in the static eccentric direction, while the maximum is in the direction away from the rotor shifting.Figure 2
Comparison of relative permeance. (a) Relative permeance with no eccentricity. (b) Relative permeance with 10% static eccentricity relative to the air gap between stator and rotor.
(a)(b)


## 2.2. Correction Factor of Relative Permeance

The relative permeance with stator slots can be obtained by introducing correction factor of stator slots. And it can be expressed as follows:(8)Λα,t=λα,tΛeα,t,where λα,t is the correction factor caused by stator slots.The correction factorλα,t can be calculated as(9)λα,t=1Kc1+∑k=1∞Akcosks1α,where Kc is Carter’s coefficient, and Ak is the harmonic permeance coefficient.Carter’s coefficientKc can be calculated as(10)kC=t1t1−γ1δ,where t1 is the slot pitch, and γ1 can be expressed as follows:(11)γ1=4π0.5κarctan0.5κ−ln1+0.5κ2.According to paper [35], the coefficient in equation (4) can be expressed as follows:(12)Ak=−βκ1k4π0.5+kΓ20.78−2kΓ2sin1.6πkΓ,where(13)βκ=0.51−11+κ2,Γ=b14t1.For permanent magnet synchronous motors,(14)κ≈b14δ+hM/μ0.


## 2.3. Expression of Unbalanced Magnetic Pull

According to the electric machine theory, the magnetomotive force (MMF) for a general motor can be expressed as(15)Fα,t=Fsα,t+Frα,t,where Fsα,t and Frα,t are the magnetomotive force of the stator windings and the permanent magnet. In this paper, the harmonics caused by current harmonics are ignored. The magnetomotive force of stator windings and the permeant magnetic can be expressed as follows:(16)Fsα,t=∑v=6k∓1∞Fsvcosvpα±ωt,k=0,1,2,…,(17)Frα,t=∑u=2n−1∞Frucosupα−uωt,n=1,2,…,where v is the harmonic order of the stator magnetomotive force, u is the harmonic order of the rotor magnetomotive force, and Fsv and Fru are the amplitude of the magnetomotive force of the stator windings and the permanent magnet. And they can be expressed as follows:(18)Fsv=1v1.35N1kdpvI1p,Fru=4πuBrhmμ0sinuαpπ2,where I1 is amplitude of the current of the stator windings, p is the number of pole-pair numbers, N1 is the number of the stator turns in series per phase, kdpv is the winding factor for the v−th harmonic, Br is the magnet remanence, hm is the magnet thickness, αp is the pole-arc/pole-pitch ratio, and μ0is the vacuum permeability.As shown in Figures3 and 4, armature MMF (stator MMF) and permeant MMF (magnetic MMF) are a longer standard trigonometric function due to the harmonic effect. It is the composition of trigonometric functions.Figure 3
Comparison of Armature MMF. (a) Arm MMF at a given time. (b) Armature MMF at a given point in time history.
(a)(b)Figure 4
Comparison of PM MMF. (a) PM MMF at a given time. (b) PM MMF at a given point in time history.
(a)(b)According to the air-gap permeance approach, the magnetic flux density of air gap can be expressed as follows, and the waveform is shown in Figure5. It shows that the air-gap flux density is nonuniform due to cogging and harmonic effects.(19)B=Λα,tfα,t.Figure 5
Comparison of air-gap flux density. (a) Air-gap flux density at a given time. (b) Air-gap flux density at a given point in time history.
(a)(b)Since the tangential magnetic flux density is much smaller than radical magnetic flux density, the tangential magnetic flux density in this paper is ignored. The radical component of Maxwell stress is expressed as(20)σr=Br22μ0.The analytical expression of UMP can be obtained by integrating the Maxwell stress on rotor surface:(21)Fux=RrL∫02πσrcosαdα,Fuy=RrL∫02πσrsinαdα.Taking equations (8), (16), and (17) into (21) and only considering the first several harmonics, the expression of UMP considering harmonics and stator slots when the pole-pair number is 1 can be expressed as(22)Fux=2+∑k=1NAk2kc2F0cosγ+F1cos5γ+F2cos2ωt−γ+F3cos2ωt−3γ+F4cos4ωt−3γ+F5cos4ωt−5γ+F6cos6ωt−γ+F7cos6ωt−5γ+F8cos2ωt+3γ+F9cos2ωt+5γ+F10cos4ωt+γ+F11cos4ωt+3γ+F12cos6ωt+γ,Fuy=2+∑k=1NAk2kc2F0sinγ+F1sin5γ+F2sin2ωt−γ+F3sin2ωt−3γ+F4sin4ωt−3γ+F5sin4ωt−5γ+F6sin6ωt−γ+F7sin6ωt−5γ+F8sin2ωt+3γ+F9sin2ωt+5γ+F10sin4ωt+γ+F11sin4ωt+3γ+F12sin6ωt+γ,where(23)F0=14Fr12+Fr32+Fr52+Fs12+Fs52+2Fr1Fs12Λ0Λ1+Λ1Λ2+Λ2Λ3,F1=2Fs1Fs512Λ2Λ3,F2=14Fr12+Fs12+2Fr1Fr3+2Fr1Fs1+2Fr3Fr5+2Fr3Fs1Λ0Λ1+12Λ1Λ2+12Λ2Λ3,F3=14Fr12+Fs12+2Fr1Fr3+2Fr1Fs1+2Fr3Fr5+2Fr3Fs1Λ0Λ3+12Λ1Λ2,F4=142Fs1Fr3+2Fs1Fr5+2Fr1Fr3+2Fr1Fr5Λ2Λ3+12Λ1Λ2,F5=142Fs1Fr3+2Fs1Fr5+2Fr1Fr3+2Fr1Fr512Λ2Λ3,F6=142Fs5Fr52Λ0Λ1+Λ1Λ2+Λ2Λ3,F7=14Fr32+2Fs1Fr5+2Fr1Fr512Λ2Λ3,F8=142Fs1Fs5+2Fr1Fs5Λ0Λ3+12Λ1Λ2,F9=142Fs5Fr1+2Fs1Fs512Λ2Λ3,F10=142Fs5Fr32Λ0Λ1+Λ1Λ2+Λ2Λ3,F11=142Fs5Fr3Λ0Λ3+12Λ1Λ2,F12=142Fs5Fr52Λ0Λ1+Λ1Λ2+Λ2Λ3.For simplicity, the coefficient of cogging effect is defined, and it can be expressed as(24)λs=∑k=1∞Ak+22kc.


## 2.4. Verification and Discussion

In order to verify the accuracy of the analysis, a finite element model (FEM) of the same parameters is established. Figure6 shows a comparison between the analytical results and finite element analytical (FEA) results of the air-gap flux density. It shows that the analytical results with cogging and harmonic effects taken into consideration are in great consistency with the FEA results. Table 1 shows the results of the analytical method and FEM method. It shows that the harmonics have the effects of increasing the UMP, while the stator slots have the effects of decreasing the UMP. In addition, the errors between the analytical results taking the harmonics and stator slots into consideration and FEM results are the smallest, which means that the analytical method is accurate, and it can be used to analyze the dynamical response of the rotor system.Figure 6
Comparison of air-gap flux density. (a) Air-gap flux density at a given time. (b) Air-gap flux density at a given point in time history.
(a)(b)Table 1
Comparison of unbalanced magnetic pull.
Static eccentricity0%10%20%30%40%50%FEM results03.34 × 1056.92 × 1051.15 × 1061.76 × 1062.65 × 106Analytical results (no harmonics and stator slots)03.04 × 1056.40 × 1051.05 × 1061.60 × 1062.40 × 106Errors0%8.98%7.51%8.70%9.09%9.43%Analytical results (with harmonics)03.83 × 1057.82 × 1051.27 × 1062.05 × 1062.85 × 106Errors0%14.67%13.01%10.43%16.48%7.55%Analytical results (with harmonics and stator slots)03.52 × 1057.31 × 1051.22 × 1061.86 × 1062.78 × 106Errors0%5.39%5.64%6.09%5.68%4.91%


## 3. Motion Equations

In this paper, the bearings are assumed to be rigid, and the gyroscopic effect is ignored. The transverse vibration differential equations can be expressed as(25)mx¨+cx˙+kx−Fux=meω2cosωt,my¨+cy˙+ky−Fuy=meω2sinωt,where m is the rotor mass, and c and k are the damping coefficient and stiffness coefficient of the rotor shaft.For the convenience of calculation, the dimensionless transformations are introduced.(26)T=ωnt,ωn=km,ξ=c2km,β=eδ0,X=xδ0,Y=yδ0,λ=ωωn.The dimensionless transformation of equation (25) is(27)X¨+2ξX˙+X−Fuxmδ0ωn2=βλ2cosλT,Y¨+2ξY˙+Y−Fuymδ0ωn2=βλ2sinλT.


## 4. Simulation and Analysis

The parameters used in this paper are shown in Table2. The simulation is conducted by MATLAB, and the stand dynamic equations are solved by fourth-order Runge–Kutta method. Results at the beginning are distorted to remove the influence of the initials.Table 2
Parameters of PMSM.
ParametersValueSystem dampingc1000 Ns/mRotor massm150 kgCurrent of the stator windingsI1940 aAverage air gapδ4 × 10−3 mShaft stiffnessk6.06 × 109 N/mMass eccentrice5 × 10−4 mStator radiusR0.17 mAir-gap lengthL0.22 m


### 4.1. Effects of the Stator Geometry Parameters on the Coefficient of Cogging Effect

In this section, the effects of ratio of slot width to average air gap, ratio of slot width to slot pitch, and ratio of air gap to slot pitch are investigated. The results are shown in Figures7–9. It shows that the coefficient of cogging effect is almost equal to 1 when the ratio of slot width to pitch is extremely small in Figure 7. That is because when the ratio is extremely small, the slots can be ignored, and the stator can be regarded as a smooth cylinder. When the ratio is much larger, the coefficient is smaller than 1 obviously, which means that the effects of slots cannot be ignored when the slot width to pitch is much larger. Figure 8 displays the effect of the ratio of slot width to pitch on the coefficient of cogging effect. It shows that the coefficient is much larger when the ratio of slot width to average air gap is equal to 4. The peak value is obtained when the ratios of slot width to pitch are equal to 0.6245, 0.3122, 0.2082, 0.1561, and 0.1249. Those ratios are gained when 0.78−2kΓ2=0. And another phenomenon can be seen that the peak value becomes larger when k is equal to a small value. Figure 9 shows that the coefficient of stator slots is almost equal to 1 when the ratio of the average air gap to slot pitch is large enough, which means that the effects of cogging can be ignored when the air gap is equal to a larger value for a certain stator. When the air gap is small, the influence cannot be ignored. Since the air gap between the stator and rotor is relatively small, the effects of cogging should be taken into consideration to invest the dynamic behavior of the rotor.Figure 7
Relationship between the coefficient and ratio of stator slot and air gap.Figure 8
Relationship between the coefficient and ratio of stator slot and pitch.Figure 9
Relationship between the coefficient and ratio of air gap and pitch.


### 4.2. Influences of Harmonics and Cogging on Dynamic Response

Figure10 shows the rotor nonlinear vibration under different rotating frequencies. It shows that the vibration taking harmonics and cogging into consideration is consistent with the orbits without taking the harmonics and cogging into consideration. However, the amplitudes of the vibration are quite different. When the harmonics is taken into consideration, the amplitudes become larger than the amplitudes without harmonics and become smaller when the cogging effects are taken into consideration at the same time. Both amplitudes are a little larger than the one without harmonics and cogging. The same conclusion can be drawn when the rotor is in high rotating speed. Figure 11(a) shows the frequency spectral of the displacement in x direction when the rotating speed is 5 Hz. It shows that the main frequencies include the components of 5 Hz and 95 Hz when the harmonics and stator slots are ignored, and the main frequencies are composed of 5 Hz, 95 Hz, and 195 Hz when the harmonics and cogging are taken into consideration. Those main frequencies are rotor rotating frequency (5 Hz), combination of the double-power frequency minus rotating frequency (95 Hz), and combination of four-power frequency minus rotating speed (195 Hz). Figure 11(b) shows the frequency spectral of the displacement in x direction when the rotating speed is 20 Hz. The main frequencies are composed of 20 Hz, 80 Hz, and 180 Hz, which are actually the rotor rotating frequency (20 Hz), combination of the double-power frequency minus rotating frequency (80 Hz), and combination of four-power frequency minus rotating frequency (180 Hz).Figure 10
p = 1, rotor shaft orbit for different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz.
(a)(b)Figure 11
p = 1, frequency spectral of displacement under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz.
(a)(b)


### 4.3. Dynamic Response of Rotor with Only Dynamic Eccentricity

In this section, static eccentricity is assumed to be zero. Rotor shaft orbits under different rotor rotating frequencies are shown in Figures12 and 13 showing the frequency waterfall chart. It shows that the amplitude of vibration increases nonlinearly with the increase in the rotating frequency. The rotor shaft orbit changes from a circle with petals around circumference to a ring, then to an ellipse, and then back to a ring. Figure 12 shows the main frequency for different rotating speeds. It shows that the main frequency includes f, 2fe−f, 4fe−f, 6fe−f and 6fe+f in all rotating frequencies. Compared with the results without considering harmonics and stator slots in paper [30], the harmonics and stator slots make the multiples of the rotating frequency and electrical frequency appear even when the rotating frequency is low.Figure 12
p = 1, Rotor orbits under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz; (c) f = 15 Hz; (d) f = 20 Hz; (e) f = 21 Hz; (f) f = 24 Hz; (g) f = 25 Hz; (h) f = 30 Hz.
(a)(b)(c)(d)(e)(f)(g)(h)Figure 13
p = 1, frequency waterfall chart of vibration under various rotating frequencies.


### 4.4. Dynamic Response of Rotor with Static Eccentricity and Dynamic Eccentricity

Figure14 shows the static eccentricity on rotor shaft orbit when the rotating frequency is equal to 10 Hz. As Figure 13 shows, the rotor shaft orbit changes from a circle with petals around circumference in Figure 12(b) to a ring, and the ring is like a gourd. Meanwhile, vibration changes from centrosymmetric to x-axis symmetric, and the displacement in y direction is much smaller than the vibration in x direction due to the static eccentric. In addition, the vibration increases, and the width of the ring becomes narrow with the increase of static eccentric. It can also be seen that the static eccentricity increases the vibration in all directions.Figure 14
p = 1, rotor vibration under different static eccentrics. (a) r0/δ = 2%; (b) r0/δ = 4%; (c) r0/δ = 6%; (d) r0/δ = 8%; (e) r0/δ = 10%; (f) r0/δ = 12%.
(a)(b)(c)(d)(e)(f)Figure15 shows the influence of static eccentric on frequency spectral. It shows that the main components of the frequency for different static eccentrics are similar due to having the same rotating frequency. The main frequencies include 2fe, 4fe, 6fe. The 4fe frequency is not obvious. The reason for those phenomena is that the static eccentric excites more harmonic frequency. For the motor in this paper, those harmonics are 2fe, 4fe, 6fe, which is different from the results in paper [30].Figure 15
p = 1, f = 10 Hz, frequency waterfall chart of vibration under various rotating frequencies.Figure16 shows the rotor vibration under different rotating frequencies with the same static eccentricity r0/δ = 10%. Results show that the vibration changes little with the increase of the rotating frequency. However, the amplitude of the vibration increases with the increase of rotating frequency, and it becomes a ring since the rotating frequency is larger than 10 Hz.Figure 16
p = 1, r0/δ = 10%, rotor orbits under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz; (c) f = 15 Hz; (d) f = 20 Hz; (e) f = 21 Hz; (f) f = 25 Hz.
(a)(b)(c)(d)(e)(f)Figure17 shows the influence of the rotating frequency on frequency spectral with the same static eccentricity r0/δ = 10%. It shows that the main components of the frequency for different rotating frequencies are similar. 2fe, 4fe, 6fe are included. But the component of 4fe is not obvious.Figure 17
p = 1, r0/δ = 10%, frequency waterfall chart of vibration under various rotating frequencies.


## 4.1. Effects of the Stator Geometry Parameters on the Coefficient of Cogging Effect

In this section, the effects of ratio of slot width to average air gap, ratio of slot width to slot pitch, and ratio of air gap to slot pitch are investigated. The results are shown in Figures7–9. It shows that the coefficient of cogging effect is almost equal to 1 when the ratio of slot width to pitch is extremely small in Figure 7. That is because when the ratio is extremely small, the slots can be ignored, and the stator can be regarded as a smooth cylinder. When the ratio is much larger, the coefficient is smaller than 1 obviously, which means that the effects of slots cannot be ignored when the slot width to pitch is much larger. Figure 8 displays the effect of the ratio of slot width to pitch on the coefficient of cogging effect. It shows that the coefficient is much larger when the ratio of slot width to average air gap is equal to 4. The peak value is obtained when the ratios of slot width to pitch are equal to 0.6245, 0.3122, 0.2082, 0.1561, and 0.1249. Those ratios are gained when 0.78−2kΓ2=0. And another phenomenon can be seen that the peak value becomes larger when k is equal to a small value. Figure 9 shows that the coefficient of stator slots is almost equal to 1 when the ratio of the average air gap to slot pitch is large enough, which means that the effects of cogging can be ignored when the air gap is equal to a larger value for a certain stator. When the air gap is small, the influence cannot be ignored. Since the air gap between the stator and rotor is relatively small, the effects of cogging should be taken into consideration to invest the dynamic behavior of the rotor.Figure 7
Relationship between the coefficient and ratio of stator slot and air gap.Figure 8
Relationship between the coefficient and ratio of stator slot and pitch.Figure 9
Relationship between the coefficient and ratio of air gap and pitch.


## 4.2. Influences of Harmonics and Cogging on Dynamic Response

Figure10 shows the rotor nonlinear vibration under different rotating frequencies. It shows that the vibration taking harmonics and cogging into consideration is consistent with the orbits without taking the harmonics and cogging into consideration. However, the amplitudes of the vibration are quite different. When the harmonics is taken into consideration, the amplitudes become larger than the amplitudes without harmonics and become smaller when the cogging effects are taken into consideration at the same time. Both amplitudes are a little larger than the one without harmonics and cogging. The same conclusion can be drawn when the rotor is in high rotating speed. Figure 11(a) shows the frequency spectral of the displacement in x direction when the rotating speed is 5 Hz. It shows that the main frequencies include the components of 5 Hz and 95 Hz when the harmonics and stator slots are ignored, and the main frequencies are composed of 5 Hz, 95 Hz, and 195 Hz when the harmonics and cogging are taken into consideration. Those main frequencies are rotor rotating frequency (5 Hz), combination of the double-power frequency minus rotating frequency (95 Hz), and combination of four-power frequency minus rotating speed (195 Hz). Figure 11(b) shows the frequency spectral of the displacement in x direction when the rotating speed is 20 Hz. The main frequencies are composed of 20 Hz, 80 Hz, and 180 Hz, which are actually the rotor rotating frequency (20 Hz), combination of the double-power frequency minus rotating frequency (80 Hz), and combination of four-power frequency minus rotating frequency (180 Hz).Figure 10
p = 1, rotor shaft orbit for different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz.
(a)(b)Figure 11
p = 1, frequency spectral of displacement under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz.
(a)(b)


## 4.3. Dynamic Response of Rotor with Only Dynamic Eccentricity

In this section, static eccentricity is assumed to be zero. Rotor shaft orbits under different rotor rotating frequencies are shown in Figures12 and 13 showing the frequency waterfall chart. It shows that the amplitude of vibration increases nonlinearly with the increase in the rotating frequency. The rotor shaft orbit changes from a circle with petals around circumference to a ring, then to an ellipse, and then back to a ring. Figure 12 shows the main frequency for different rotating speeds. It shows that the main frequency includes f, 2fe−f, 4fe−f, 6fe−f and 6fe+f in all rotating frequencies. Compared with the results without considering harmonics and stator slots in paper [30], the harmonics and stator slots make the multiples of the rotating frequency and electrical frequency appear even when the rotating frequency is low.Figure 12
p = 1, Rotor orbits under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz; (c) f = 15 Hz; (d) f = 20 Hz; (e) f = 21 Hz; (f) f = 24 Hz; (g) f = 25 Hz; (h) f = 30 Hz.
(a)(b)(c)(d)(e)(f)(g)(h)Figure 13
p = 1, frequency waterfall chart of vibration under various rotating frequencies.


## 4.4. Dynamic Response of Rotor with Static Eccentricity and Dynamic Eccentricity

Figure14 shows the static eccentricity on rotor shaft orbit when the rotating frequency is equal to 10 Hz. As Figure 13 shows, the rotor shaft orbit changes from a circle with petals around circumference in Figure 12(b) to a ring, and the ring is like a gourd. Meanwhile, vibration changes from centrosymmetric to x-axis symmetric, and the displacement in y direction is much smaller than the vibration in x direction due to the static eccentric. In addition, the vibration increases, and the width of the ring becomes narrow with the increase of static eccentric. It can also be seen that the static eccentricity increases the vibration in all directions.Figure 14
p = 1, rotor vibration under different static eccentrics. (a) r0/δ = 2%; (b) r0/δ = 4%; (c) r0/δ = 6%; (d) r0/δ = 8%; (e) r0/δ = 10%; (f) r0/δ = 12%.
(a)(b)(c)(d)(e)(f)Figure15 shows the influence of static eccentric on frequency spectral. It shows that the main components of the frequency for different static eccentrics are similar due to having the same rotating frequency. The main frequencies include 2fe, 4fe, 6fe. The 4fe frequency is not obvious. The reason for those phenomena is that the static eccentric excites more harmonic frequency. For the motor in this paper, those harmonics are 2fe, 4fe, 6fe, which is different from the results in paper [30].Figure 15
p = 1, f = 10 Hz, frequency waterfall chart of vibration under various rotating frequencies.Figure16 shows the rotor vibration under different rotating frequencies with the same static eccentricity r0/δ = 10%. Results show that the vibration changes little with the increase of the rotating frequency. However, the amplitude of the vibration increases with the increase of rotating frequency, and it becomes a ring since the rotating frequency is larger than 10 Hz.Figure 16
p = 1, r0/δ = 10%, rotor orbits under different rotating frequencies. (a) f = 5 Hz; (b) f = 10 Hz; (c) f = 15 Hz; (d) f = 20 Hz; (e) f = 21 Hz; (f) f = 25 Hz.
(a)(b)(c)(d)(e)(f)Figure17 shows the influence of the rotating frequency on frequency spectral with the same static eccentricity r0/δ = 10%. It shows that the main components of the frequency for different rotating frequencies are similar. 2fe, 4fe, 6fe are included. But the component of 4fe is not obvious.Figure 17
p = 1, r0/δ = 10%, frequency waterfall chart of vibration under various rotating frequencies.


## 5. Conclusions

In this paper, the expression of UMP, taking harmonics and cogging effects into consideration, was obtained when the pole-pair number is 1. Motion equations based on the expression were established and solved. The effects of stator structure parameters on coefficient of cogging effects were analyzed. And influences of eccentricities, harmonics, and cogging on dynamical responses were investigated. The main conclusions can be drawn as follows:(1).
The stator parameters have great influences on coefficient of stator slots. When the ratio of slot width to pitch was extremely small, the slots can be ignored, and the stator can be regarded as a smooth cylinder. When the ratio is much larger, cogging effects cannot be ignored. Especially, when the ratios of slot width to pitch are equal to 0.6245, 0.3122, 0.2082, 0.1561, and 0.1249, the coefficient of cogging effects is larger than that of the others, which means that the unbalanced pull is much larger and should be avoided in the design of motors. For a certain stator, the effects of stator slots can be ignored when the air gap is large enough.(2)
Harmonics caused by permanent magnetic and stator windings increased vibration of the rotor, while the vibration decreased when the cogging effects were taken into consideration simultaneously, which means that the cogging effects have the effects of decreasing the amplitude. In addition, the harmonics and cogging effects do not change the components of the main frequency.(3)
Vibration changed under different static eccentricities and rotating frequencies. The main components of frequency included2fe, 4fe, 6fe when static eccentric existed, which is different from the results when there is no static eccentricity. The components of the main frequency under different rotating frequencies were composed of 2fe, 4fe, 6fe with static eccentricity as well.



---

*Source: 6685588-2021-03-05.xml*