# Rotation Modes Stability Analysis and Phase Compensation for Magnetically Suspended Flywheel Systems with Cross Feedback Controller and Time Delay

Authors: Yuan Ren; Xiaocen Chen; Yuanwen Cai; Weijie Wang
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

This paper analyzes the effects of time delay on the stability of the rotation modes for the magnetically suspended flywheel (MSFW) with strong gyroscopic effects. A multi-input multioutput system is converted into a single-input single-output control system with complex coefficient by variable reconstruction, and the stability equivalence of the systems before and after variable reconstruction is proven. For the rotation modes, the stability limits and corresponding vibration frequencies are found as a function of nondimensional magnetic stiffness and damping and nondimensional parameters of rotor speed and time delay. Additionally, the relationship between cross feedback control system stability and time delay is investigated. And an effective phase compensation method based on cross-channel is further presented. Simulation and experimental results are presented to demonstrate the correctness of the stability analysis method and the superiority of the phase compensation strategy.


---

## Body

## 1. Introduction

With the development of magnetic bearing technology, a magnetically suspended flywheel (MSFW) represents an exciting alternative to the traditional mechanical flywheel due to its inherent superior features such as contact free operation, zero friction, adjustable damping, and stiffness as well as the potential for low vibration and super high rotational speeds [1, 2]. As a result, it has been widely employed in power industry [3, 4], transportation, spaceflight [5–8], and so on [9, 10].To maximize the energy storage capacity, a MSFW usually chooses a radially thick rotor so that the principle mass moment of inertial is larger than the transversal counterpart, that is,Jz/Jx>1. In this case, the gyroscopic effects are especially strong, which puts a strain on the stability analysis of the system.Particularly, to improve the transient response and to reduce power consumption of a magnetic bearing, the permanent-magnetic-biased hybrid magnetic bearings (PMHMBs) are applied [1]. This bias field causes the rotor to behave as it was supported on a spring with negative stiffness. There is no time delay associated with the bias field since it is caused by a permanent magnet. However, the time delay of the electromagnetic bearing always exists since each of the components involved in it and its control system have a time delay associated with them. For an analog device, the time delay can be approximated by a linear time delay low-pass filter. The amplitude of the transfer function is assumed to be a constant and the phase is linear. However, for the digital devices, this time delay is the loop time for each measurement [11]. The total time delay is the sum of the individual time delay [12]. All of these make the stability analysis of the MSFW more complex and difficult.In recent years, considerable researches have been carried out for the stability analysis of the rotor supported by magnetic bearings with time delay. Ji [13, 14] analyzed the effect of time delays on the nonlinear dynamic behavior of the magnetic bearing system. Wang and Liu [15], Li et al. [16], and Ji [17] studied the Hopf bifurcation and its stability in a magnetic bearing system with time delays. However, these works only considered a simple single-degree-of-freedom (DOF) magnetic bearing without gyroscopic effects. As for the multi-DOF cases with high rotor speed and strong gyroscopic effects, Kascak et al. [18] employed nondimensional magnetic stiffness and damping as well as nondimensional system parameters to describe the stability limits of a PD controller. However, they did not consider the cross feedback controller. Kascak et al. [19] presented the stability boundaries of a spinning rotor with parametrically excited gyroscopic system. Bauomy [20] performed the stability analysis of a rotor active magnetic bearing system with time varying stiffness. Inayat-Hussain [21] studied the nonlinear dynamics of a magnetically supported rigid rotor in auxiliary bearings. Andersen et al. [22] and Sugai et al. [23] analyzed the rotation modes stability of the magnetically suspended rotor with passive cylinder-magnet bearings. Additionally, Pei [24] presented the stability boundaries of a spinning rotor with parametrically excited gyroscopic system. Fan and Pan [25] described the stability of rotor-shaft system including electromagnetic exciters. All in all, few literatures studied the multi-DOF magnetic bearings with time delays and gyroscopic effects.To resolve the gyroscopic effects, considerable research has been conducted and various approaches have been proposed [26–32]. Among these methods, decentralized PD plus filtered cross feedback control is the most well-known and easy-to-implement approach [29]. To further improve the system nutation stability, different phase compensation methods are developed, such as unsymmetrical current sampling resistance method [11, 12] and inverse system method with phase compensation [33, 34], where the phase compensations are performed with both decentralized and cross-channels. For simplicity, it is called the combinational compensation method in this paper. In fact, there are other phase compensation strategies in addition to this, such as phase compensation with decentralized channels and phase compensation with cross-channels. Thus, the question naturally raised is what is the most effective way for the system stability among the three phase compensation methods? This has not been effectively resolved.To resolve the problems above, this paper describes a modeling method which is used to theoretically map the stability boundaries of the MSFW controller with cross feedback and time delay. As for the rotation modes, the relationships between the time delays and the critical nutation frequency are further analytically analyzed, and then an effective phase compensation strategy is proposed.This paper is organized as follows. First, in Section2, the system is modeled with complex coefficient. Then, its rotation modes stability is further studied in Section 3. Simulation and experiments are developed in Section 4. Finally, Section 5 concludes this paper.


## 2. Modeling of the MSFW System and Its Variable Reconstruction

### 2.1. Modeling of the MSFW System

Here, we review the model of a rigid rotor supported by magnetic bearings, which has been established in [7]. Figure 1 shows the coordinate system definition of a MSFW, where O is the geometric center of the magnetic bearing stator; lm refers to the distance between the central point of radial magnetic bearing and point O; Ω is the rotor speed, X-, Y-, and Z-axis form the generalized coordinate system of the rotor position; px and py are the torques in the X and Y directions, respectively; α and β are the rotor angular displacements about the X- and Y-axis. fax, fay, fbx, fby, faz, and fbz are the magnetic forces along the magnetic bearing coordinate system AX-, AY-, BX-, BY-, AZ-, and BZ-axis.Figure 1
Force and coordinate system definition of a MSFW.Then, the dynamic model of a rigid rotor supported by PMHMBs can be described as(1)mx¨=fax+fbx,Jyβ¨-JzΩα˙=py=lmfax-fbx,my¨=fay+fby,Jxα¨+JzΩβ˙=px=lmfby-fay,mz¨=faz+fbz,where m is the mass of the flywheel rotor; Jx, Jy, and Jz are the moments of inertia of the rotor about the X-, Y-, and Z-axis, respectively, and Jx=Jy=Jr; and x, y, and z are the linear displacements of the mass center of the rotor from point O in the X-, Y-, and Z-axis, respectively.Figure2 shows the schematic of the cross feedback control system, where gb(s), gcr(s), and gf(s) are the transfer functions of the basic controller (such as a PID controller), the cross feedback controller, and the antialias filter, respectively.Figure 2
Schematic diagram of the decentralized PD plus filter cross feedback controller.As for an axisymmetric rotor supported by PMHMBs, the magnetic forcefλ can be linearized as [7](2)fλ=kiλiλ+khλhλ,λ=ax,ay,bx,by,az,bz,where kiλ and khλ denote, respectively, the current stiffness and displacement stiffness of the λ channel. Note that the four radial channels of the MSFW are often symmetric; for simplicity, ki and kh denote the current and displacement stiffness of the radial magnetic bearings, respectively.For PMHMBs, there is no time delay associated with the bias field since it is produced by a permanent magnet [24], and the delay only exists in the electromagnetic field. Accordingly, the closed-loop model of the control system with the cross feedback controller can be given by [33, 35](3)mx¨t-2khxt=-2kiksgbxt-τ,Jyβ¨t-JzΩα˙t-2khlm2βt=-2lmlskiksgbβt-τ+gcrαt-τ,my¨t-2khyt=-2kiksgbyt-τ,Jxα¨t+JzΩβ˙t-2khlm2αt=-2lmlskiksgbαt-τ-gcrβt-τ,mz¨t-2khzzt=-2kizksgbzzt-τ,where ga denotes the power amplifier of the HMBs, ks is the proportional gain of the displacement sensor, ls is the distance from the geometric center of the rotor to the sensor, and τ is the control system delay.It can be seen from (3) that the cross feedback control only influences the system rotation modes stability and does not endanger the system translation modes stability. As for the translation modes of the MSFW, their stability analysis can be equivalent to that of the single-DOF magnetic bearing control system, which has been clearly analyzed. Accordingly, this paper focuses on the rotation modes stability analysis of the MSFW system with cross feedback control system.The rotation modes model can be described as(4)Jyβ¨t-JzΩα˙t-2khlm2βt=-2lmlskiksgbβt-τ+gcrαt-τ,Jxα¨t+JzΩβ˙t-2khlm2αt=-2lmlskiksgbαt-τ-gcrβt-τ.


### 2.2. Variable Reconstruction and Its Stability Equivalence Analysis

As for the rotation modes model (4), it is a two-input and two-output control system. To simplify the stability analysis of the MSFW rotation modes, variable reconstruction is employed in this paper.From the coordinate definition of the MSFW shown in Figure1, α leads β 90 degrees. Accordingly, we can define φ=α+jβ, where j is the imaginary number unit and j2=-1. Note that the dissymmetry of (4) and multiplying the second equation by j and then adding the result to the fourth equation yield the differential equations:(5)Jrφ¨t-jHφ˙t-2khlm2φt+Pφt-τ+Dφ˙t-τ+jkcHφ˙t-τ=0,where P=2lmlskikskp,D=2lmlskikskd.According to the control system theory, the necessary and sufficient conditions of the rotation modes stability are as follows:(6)limt→∞αt=0,limt→∞βt=0.Accordingly,(7)limt→∞αt+jβt=0.As for complex-coefficient system (5), the necessary and sufficient condition of the stable control system is (8)limt→∞ϕt=0.Note thatϕ(t)=α(t)+jβ(t); that is, limt→∞α(t)=0 and limt→∞β(t)=0.Therefore, the stability of system (4) is equivalent to that of (5). Then, a two-input two-output system is converted into a single-input single-output control system with complex coefficient. That is, we can get the stability of system (4) by analyzing the stability of system (5).


## 2.1. Modeling of the MSFW System

Here, we review the model of a rigid rotor supported by magnetic bearings, which has been established in [7]. Figure 1 shows the coordinate system definition of a MSFW, where O is the geometric center of the magnetic bearing stator; lm refers to the distance between the central point of radial magnetic bearing and point O; Ω is the rotor speed, X-, Y-, and Z-axis form the generalized coordinate system of the rotor position; px and py are the torques in the X and Y directions, respectively; α and β are the rotor angular displacements about the X- and Y-axis. fax, fay, fbx, fby, faz, and fbz are the magnetic forces along the magnetic bearing coordinate system AX-, AY-, BX-, BY-, AZ-, and BZ-axis.Figure 1
Force and coordinate system definition of a MSFW.Then, the dynamic model of a rigid rotor supported by PMHMBs can be described as(1)mx¨=fax+fbx,Jyβ¨-JzΩα˙=py=lmfax-fbx,my¨=fay+fby,Jxα¨+JzΩβ˙=px=lmfby-fay,mz¨=faz+fbz,where m is the mass of the flywheel rotor; Jx, Jy, and Jz are the moments of inertia of the rotor about the X-, Y-, and Z-axis, respectively, and Jx=Jy=Jr; and x, y, and z are the linear displacements of the mass center of the rotor from point O in the X-, Y-, and Z-axis, respectively.Figure2 shows the schematic of the cross feedback control system, where gb(s), gcr(s), and gf(s) are the transfer functions of the basic controller (such as a PID controller), the cross feedback controller, and the antialias filter, respectively.Figure 2
Schematic diagram of the decentralized PD plus filter cross feedback controller.As for an axisymmetric rotor supported by PMHMBs, the magnetic forcefλ can be linearized as [7](2)fλ=kiλiλ+khλhλ,λ=ax,ay,bx,by,az,bz,where kiλ and khλ denote, respectively, the current stiffness and displacement stiffness of the λ channel. Note that the four radial channels of the MSFW are often symmetric; for simplicity, ki and kh denote the current and displacement stiffness of the radial magnetic bearings, respectively.For PMHMBs, there is no time delay associated with the bias field since it is produced by a permanent magnet [24], and the delay only exists in the electromagnetic field. Accordingly, the closed-loop model of the control system with the cross feedback controller can be given by [33, 35](3)mx¨t-2khxt=-2kiksgbxt-τ,Jyβ¨t-JzΩα˙t-2khlm2βt=-2lmlskiksgbβt-τ+gcrαt-τ,my¨t-2khyt=-2kiksgbyt-τ,Jxα¨t+JzΩβ˙t-2khlm2αt=-2lmlskiksgbαt-τ-gcrβt-τ,mz¨t-2khzzt=-2kizksgbzzt-τ,where ga denotes the power amplifier of the HMBs, ks is the proportional gain of the displacement sensor, ls is the distance from the geometric center of the rotor to the sensor, and τ is the control system delay.It can be seen from (3) that the cross feedback control only influences the system rotation modes stability and does not endanger the system translation modes stability. As for the translation modes of the MSFW, their stability analysis can be equivalent to that of the single-DOF magnetic bearing control system, which has been clearly analyzed. Accordingly, this paper focuses on the rotation modes stability analysis of the MSFW system with cross feedback control system.The rotation modes model can be described as(4)Jyβ¨t-JzΩα˙t-2khlm2βt=-2lmlskiksgbβt-τ+gcrαt-τ,Jxα¨t+JzΩβ˙t-2khlm2αt=-2lmlskiksgbαt-τ-gcrβt-τ.


## 2.2. Variable Reconstruction and Its Stability Equivalence Analysis

As for the rotation modes model (4), it is a two-input and two-output control system. To simplify the stability analysis of the MSFW rotation modes, variable reconstruction is employed in this paper.From the coordinate definition of the MSFW shown in Figure1, α leads β 90 degrees. Accordingly, we can define φ=α+jβ, where j is the imaginary number unit and j2=-1. Note that the dissymmetry of (4) and multiplying the second equation by j and then adding the result to the fourth equation yield the differential equations:(5)Jrφ¨t-jHφ˙t-2khlm2φt+Pφt-τ+Dφ˙t-τ+jkcHφ˙t-τ=0,where P=2lmlskikskp,D=2lmlskikskd.According to the control system theory, the necessary and sufficient conditions of the rotation modes stability are as follows:(6)limt→∞αt=0,limt→∞βt=0.Accordingly,(7)limt→∞αt+jβt=0.As for complex-coefficient system (5), the necessary and sufficient condition of the stable control system is (8)limt→∞ϕt=0.Note thatϕ(t)=α(t)+jβ(t); that is, limt→∞α(t)=0 and limt→∞β(t)=0.Therefore, the stability of system (4) is equivalent to that of (5). Then, a two-input two-output system is converted into a single-input single-output control system with complex coefficient. That is, we can get the stability of system (4) by analyzing the stability of system (5).


## 3. Stability Analysis of Rotation Modes of the MSFW

### 3.1. Stability Analysis of the Rotation Modes

According to the classical small signal stability analysis method, assuming an eigenvalue solution of the formφ(t)=Aeλt, then (9)φt-τ=Aeλt-τ.Substituting (9) into (5) yields the characteristic equation:(10)Jrλ2-jHλ-2khlm2+Pe-λτ+Dλe-λτ+jkcHλe-λτ=0.Define the nondimensional parameters as follows:(11)μ=λτ,τ02=Pτ2Jr,2ξ1τ0=DτJr,2ξ2τ0=kcHτJr,τn2=2khlm2τ2Jr,τη=τHJr.Then, the eigenvalue equation can be simplified as follows:(12)μ2-jτημ-τn2eu+2τ0ξ1+jξ2μ+τ02=0.The eigenvalue equation has the similar solution to that of the translation mode under the condition the rotor speed of the MSFW is zero. If the eigenvalue is complex, then the imaginary part of the eigenvalue is the circular frequency of the vibration that multiplies the time delay; that is(13)μ=α+jβ=α+jωτ.Substituting (13) into (12) yields(14)Im⁡μ2-jτημ-τn2eu+2τ0ξ1β+2τ0ξ2α=0,Reμ2-jτημ-τn2eu+2τ0ξ1α-2τ0ξ2β+τ02=0.Solving the real and imaginary portions of the eigenvalue equation, the nondimension damping obtains(15)2τ0ξ1=-Im⁡μ2-jτημ-τn2euβ-2τ0ξ2αβ.And the nondimensional stiffness achieves(16)τ02=-Reμ2-jτημ-τn2eu+2τ0ε2β.If the real part of the eigenvalueα is positive, the vibrations grow in time and the system is unstable. Contrarily, if α is negative, the vibrations decay in time and the system is stable. The case that α equals zero defines the stability boundary as follows:(17)2τ0ξ1=β2-τηβ+τn2sin⁡ββ,τ02=β2-τηβ+τn2cos⁡β+ξ2ξ1sin⁡β.Note that this paper mainly focuses on the nutation modes stability analysis and phase compensation for MSFW with cross feedback controller and time delay, the possible hidden oscillations [36], Bifurcation and global stability analysis being the next topics in the near future.


### 3.2. Relationship between Cross Feedback Control System Stability and Time Delay

To further analyze how the control system delay influences the nutation stability, the following analysis has been developed, and an effective phase compensation strategy is proposed.The eigenvalue equation of (5) can be resolved as(18)Jrλ2-jHλ-2khlm2+Pe-λτ+Dλe-λτ+jHλe-λτ=0.Suppose its eigenvalue isλ=σ+jω; then, substituting this into (18) yields(19)Jrσ2-ω2+Hω-2khlm2+Pe-στ1cos⁡ωτ1+Dσe-στ1cos⁡ωτ1+Mωe-στ1sin⁡ωτ1-Hωcos⁡ωτ2+Hσsin⁡ωτ2=0,2Jrσω-Hσ-Pe-στ1sin⁡ωτ1+De-στ1ωcos⁡ωτ1-σsin⁡ωτ1+Hσcos⁡ωτ2+Hωsin⁡ωτ2=0.(1)
Ifσ=ω=0, then we have L=2khlm2. The system has a simple steady state bifurcation.(2)
Ifσ=0,ω≠0, the system is critically stable; (19) can be simplified as(20)Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2=Jrω2-Hω+2khlm2-Psin⁡ωτ1+Dωcos⁡ωτ1+Hωsin⁡ωτ2=0.Moving the terms of trigonometric functions into the left-hand side and then squaring and adding these two equations yield a quadratic equation forω:(21)P2+D2ω2+H2ω2-2PHωcos⁡τ1-τ2-2DHω2sin⁡ωτ1-ωτ2=Jrω2-Hω+2khlm22.Whenτ1=τ2, the above equation can be further simplified as(22)P2+D2ω2+H2ω2-2PHω=Jrω2-Hω+2khlm22.Sumptuously, multiplyingsin⁡ωτ2 and cos⁡ωτ2 to the two equations of (19) and then adding them together yield(23)Psin⁡ωτ1-ωτ2-Dωcos⁡ωτ1-ωτ2+Jrω2-Hω+2khlm2sin⁡ωτ2=0.Note thatτ1=τ2 as, for a MSFW control system without phase compensation, (23) can be rewritten as (24)Dω=Jrω2-Hω+2khlm2sin⁡ωτ2.Accordingly, the second equation of (20) achieves(25)Hω-Psin⁡ωτ+Dωcos⁡ωτ=0.Note that, for the MSFW with high rotor speed,Hω≫P. As for positive ω, we have(26)cos⁡ωτ<0,sin⁡ωτ>0.That is, ωτ∈2kπ+π/2,2kπ+π,k=0,1,2,….Combining (24) and (26) yields(27)Jrω2-Hω+2khlm2>0.Combing (25) and (26) obtains(28)τ2=π-arctan⁡Dω/Hω-Pω.Deriving the above equation with ω yields(29)dτ2dω=DPω/d2ω2+Hω-P2+arctan⁡Dω/Hω-P-πω2<DPω/2DωHω-P+arctan⁡Dω/Hω-P-πω2.Note thatHω≫P and H≫D in practice; therefore, dτ2/dω<0. That is, the larger the delay time is, the smaller the critical nutation frequency is.Theorem 1.
If the rotor speed is limited to zero,limφ→0-∂φ/∂τ2/∂φ/∂τ1≪1, and if the rotor speed is high enough, limφ→0-∂φ/∂τ2/∂φ/∂τ1≫1.Proof.
Multiplying the first and the second item of (25) by sin⁡(ωτ2)and cos⁡(ωτ2), respectively, and then adding them together yields(30)φ=arctan⁡Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2+nπ.Then,(31)∂φ∂τ1=11+ξ2-Dω2sin⁡ωτ1-Pωcos⁡ωτ1γ--Psin⁡ωτ1+Dω2cos⁡ωτ1ηγ2,where(32)γ=Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2,η=Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2,ξ=Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2.
Therefore,(33)limφ→0-ηγ=0.
Hence,(34)limφ→0-∂φ/∂τ2∂φ/∂τ1=Hω2cos⁡ωτ-Dω2sin⁡ωτ-Pωcos⁡ωτ.
Whenφ=0, we have(35)Dωcos⁡ωτ-Psin⁡ωτ+Hωsin⁡ωτ=0.Noting that cos⁡(ωτ)≠0 yields(36)tan⁡ωτ=DωP-Hω.
Substituting (36) into (34) yields(37)limφ→0-∂φ/∂τ2∂φ/∂τ1=Hω2-Dω2tan⁡ωτ-Pω=Hω2-D2ω3/P-Hω-Pω=HωP-HωPH-P2-D2ω2.
Obviously, if the rotor speed is limited to zero,limφ→0-∂φ/∂τ2/∂φ/∂τ1≪1, compared with the cross-channel, the decentralized channel is the main factor which endangers the system nutation stability. If the rotor speed is high enough, limφ→0-∂φ/∂τ2/∂φ/∂τ1≫1, the delay time of the cross-channel is the main factor which influences the system nutation stability. This completes the proof.According to Theorem1, we can further draw the following corollary.Corollary 2.
As for the control system (23), if τ2=0, then the stability of the rotation modes of the MSFW has nothing to do with the rotor speed.Proof.
Sinceτ2=0, (24) can be simplified as (38)Jrφ¨t-2khlm2φt+Pφt-τ1+Dφ˙t-τ1=0.The stability analysis is similar to that of the translation modes of the MSFW. This completes the proof.Accordingly, as for the MSFW system with high rotor speed, it is the delay time of the cross-channel that influences seriously the system nutation stability instead of the decentralized channel. As a result, cross-channel phase compensation is more effective than that of the decentralized channel at high speeds and the nutation stability of the MSFW can be further improved.Similarly, it can be proven that the cross-channel compensation method is more effective than that of the combinational compensation method. Accordingly, at high rotor speeds, cross-channel compensation is the most effective phase compensation method of the three compensation methods.


## 3.1. Stability Analysis of the Rotation Modes

According to the classical small signal stability analysis method, assuming an eigenvalue solution of the formφ(t)=Aeλt, then (9)φt-τ=Aeλt-τ.Substituting (9) into (5) yields the characteristic equation:(10)Jrλ2-jHλ-2khlm2+Pe-λτ+Dλe-λτ+jkcHλe-λτ=0.Define the nondimensional parameters as follows:(11)μ=λτ,τ02=Pτ2Jr,2ξ1τ0=DτJr,2ξ2τ0=kcHτJr,τn2=2khlm2τ2Jr,τη=τHJr.Then, the eigenvalue equation can be simplified as follows:(12)μ2-jτημ-τn2eu+2τ0ξ1+jξ2μ+τ02=0.The eigenvalue equation has the similar solution to that of the translation mode under the condition the rotor speed of the MSFW is zero. If the eigenvalue is complex, then the imaginary part of the eigenvalue is the circular frequency of the vibration that multiplies the time delay; that is(13)μ=α+jβ=α+jωτ.Substituting (13) into (12) yields(14)Im⁡μ2-jτημ-τn2eu+2τ0ξ1β+2τ0ξ2α=0,Reμ2-jτημ-τn2eu+2τ0ξ1α-2τ0ξ2β+τ02=0.Solving the real and imaginary portions of the eigenvalue equation, the nondimension damping obtains(15)2τ0ξ1=-Im⁡μ2-jτημ-τn2euβ-2τ0ξ2αβ.And the nondimensional stiffness achieves(16)τ02=-Reμ2-jτημ-τn2eu+2τ0ε2β.If the real part of the eigenvalueα is positive, the vibrations grow in time and the system is unstable. Contrarily, if α is negative, the vibrations decay in time and the system is stable. The case that α equals zero defines the stability boundary as follows:(17)2τ0ξ1=β2-τηβ+τn2sin⁡ββ,τ02=β2-τηβ+τn2cos⁡β+ξ2ξ1sin⁡β.Note that this paper mainly focuses on the nutation modes stability analysis and phase compensation for MSFW with cross feedback controller and time delay, the possible hidden oscillations [36], Bifurcation and global stability analysis being the next topics in the near future.


## 3.2. Relationship between Cross Feedback Control System Stability and Time Delay

To further analyze how the control system delay influences the nutation stability, the following analysis has been developed, and an effective phase compensation strategy is proposed.The eigenvalue equation of (5) can be resolved as(18)Jrλ2-jHλ-2khlm2+Pe-λτ+Dλe-λτ+jHλe-λτ=0.Suppose its eigenvalue isλ=σ+jω; then, substituting this into (18) yields(19)Jrσ2-ω2+Hω-2khlm2+Pe-στ1cos⁡ωτ1+Dσe-στ1cos⁡ωτ1+Mωe-στ1sin⁡ωτ1-Hωcos⁡ωτ2+Hσsin⁡ωτ2=0,2Jrσω-Hσ-Pe-στ1sin⁡ωτ1+De-στ1ωcos⁡ωτ1-σsin⁡ωτ1+Hσcos⁡ωτ2+Hωsin⁡ωτ2=0.(1)
Ifσ=ω=0, then we have L=2khlm2. The system has a simple steady state bifurcation.(2)
Ifσ=0,ω≠0, the system is critically stable; (19) can be simplified as(20)Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2=Jrω2-Hω+2khlm2-Psin⁡ωτ1+Dωcos⁡ωτ1+Hωsin⁡ωτ2=0.Moving the terms of trigonometric functions into the left-hand side and then squaring and adding these two equations yield a quadratic equation forω:(21)P2+D2ω2+H2ω2-2PHωcos⁡τ1-τ2-2DHω2sin⁡ωτ1-ωτ2=Jrω2-Hω+2khlm22.Whenτ1=τ2, the above equation can be further simplified as(22)P2+D2ω2+H2ω2-2PHω=Jrω2-Hω+2khlm22.Sumptuously, multiplyingsin⁡ωτ2 and cos⁡ωτ2 to the two equations of (19) and then adding them together yield(23)Psin⁡ωτ1-ωτ2-Dωcos⁡ωτ1-ωτ2+Jrω2-Hω+2khlm2sin⁡ωτ2=0.Note thatτ1=τ2 as, for a MSFW control system without phase compensation, (23) can be rewritten as (24)Dω=Jrω2-Hω+2khlm2sin⁡ωτ2.Accordingly, the second equation of (20) achieves(25)Hω-Psin⁡ωτ+Dωcos⁡ωτ=0.Note that, for the MSFW with high rotor speed,Hω≫P. As for positive ω, we have(26)cos⁡ωτ<0,sin⁡ωτ>0.That is, ωτ∈2kπ+π/2,2kπ+π,k=0,1,2,….Combining (24) and (26) yields(27)Jrω2-Hω+2khlm2>0.Combing (25) and (26) obtains(28)τ2=π-arctan⁡Dω/Hω-Pω.Deriving the above equation with ω yields(29)dτ2dω=DPω/d2ω2+Hω-P2+arctan⁡Dω/Hω-P-πω2<DPω/2DωHω-P+arctan⁡Dω/Hω-P-πω2.Note thatHω≫P and H≫D in practice; therefore, dτ2/dω<0. That is, the larger the delay time is, the smaller the critical nutation frequency is.Theorem 1.
If the rotor speed is limited to zero,limφ→0-∂φ/∂τ2/∂φ/∂τ1≪1, and if the rotor speed is high enough, limφ→0-∂φ/∂τ2/∂φ/∂τ1≫1.Proof.
Multiplying the first and the second item of (25) by sin⁡(ωτ2)and cos⁡(ωτ2), respectively, and then adding them together yields(30)φ=arctan⁡Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2+nπ.Then,(31)∂φ∂τ1=11+ξ2-Dω2sin⁡ωτ1-Pωcos⁡ωτ1γ--Psin⁡ωτ1+Dω2cos⁡ωτ1ηγ2,where(32)γ=Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2,η=Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2,ξ=Dωcos⁡ωτ1-Psin⁡ωτ1+Hωsin⁡ωτ2Pcos⁡ωτ1+Dωsin⁡ωτ1-Hωcos⁡ωτ2.
Therefore,(33)limφ→0-ηγ=0.
Hence,(34)limφ→0-∂φ/∂τ2∂φ/∂τ1=Hω2cos⁡ωτ-Dω2sin⁡ωτ-Pωcos⁡ωτ.
Whenφ=0, we have(35)Dωcos⁡ωτ-Psin⁡ωτ+Hωsin⁡ωτ=0.Noting that cos⁡(ωτ)≠0 yields(36)tan⁡ωτ=DωP-Hω.
Substituting (36) into (34) yields(37)limφ→0-∂φ/∂τ2∂φ/∂τ1=Hω2-Dω2tan⁡ωτ-Pω=Hω2-D2ω3/P-Hω-Pω=HωP-HωPH-P2-D2ω2.
Obviously, if the rotor speed is limited to zero,limφ→0-∂φ/∂τ2/∂φ/∂τ1≪1, compared with the cross-channel, the decentralized channel is the main factor which endangers the system nutation stability. If the rotor speed is high enough, limφ→0-∂φ/∂τ2/∂φ/∂τ1≫1, the delay time of the cross-channel is the main factor which influences the system nutation stability. This completes the proof.According to Theorem1, we can further draw the following corollary.Corollary 2.
As for the control system (23), if τ2=0, then the stability of the rotation modes of the MSFW has nothing to do with the rotor speed.Proof.
Sinceτ2=0, (24) can be simplified as (38)Jrφ¨t-2khlm2φt+Pφt-τ1+Dφ˙t-τ1=0.The stability analysis is similar to that of the translation modes of the MSFW. This completes the proof.Accordingly, as for the MSFW system with high rotor speed, it is the delay time of the cross-channel that influences seriously the system nutation stability instead of the decentralized channel. As a result, cross-channel phase compensation is more effective than that of the decentralized channel at high speeds and the nutation stability of the MSFW can be further improved.Similarly, it can be proven that the cross-channel compensation method is more effective than that of the combinational compensation method. Accordingly, at high rotor speeds, cross-channel compensation is the most effective phase compensation method of the three compensation methods.


## 4. Simulation and Experimental Study

In order to demonstrate the effectiveness and validity of the proposed stability analysis and phase compensation strategy, and to reveal how closely the theory represents the physical system, both simulation and experiments have been developed.


### 4.1. Simulation and Experimental Setup

The large magnetically suspended control moment gyroscope (MSCMG), as shown in Figure3, is employed as the simulation and experimental setup since its MSFW has strong gyroscopic effects.Figure 3
Photograph of the fabricated laboratory setup.The well-known decentralized PD plus filter cross feedback controller is employed. The main parameters of the MSCMG and the controllers used in simulation and experiments have been listed in Tables1 and 2. In Table 1, x0 is the air gap of the radial magnetic bearings; ki0 and kh0 are the current stiffness and displacement stiffness of the radial magnetic bearings under the conditions iγ=0 and hmγ=0(γ=ax,ay,bx,by); and R and L are the coil resistance and inductance of the radial magnetic bearings. In Table 2, kam and ic are the proportional and feedback coefficients of the current controller for the MSFW system, which are the same as those of the proposed controller; kp and kd are the proportional, integral, and differential coefficients of the PD plus cross feedback controller for the MSFW system; khc and klc are the cross coefficients of the high-pass and low-pass filters, respectively; kc is the total cross coefficient; and fh and fl are the cut-off frequencies of the two-order high-pass and low-pass filters.Table 1
System parameters of the MSCMG.
Parameter
Value
m
56 kg
J
x
0.6366 kg⋅m2
J
y
0.6366 kg⋅m2
J
z
0.7958 kg⋅m2
k
i
0
1140 N/A
x
0
250μm
L
60.5 mH
R
2.5Ω
R
r
12000 r/min
l
s
0.177 m
l
m
0.113 m
k
h
0
2.83 N/μmTable 2
Parameters of the controller.
Parameter
Value
k
a
m
2.0
k
p
2.5
k
h
c
0.0045
k
c
1.0
f
l
30 Hz
k
r
i
0.04
k
c
p
30
i
c
10
k
d
0.008
k
l
c
0.00015
f
h
300 Hz
k
r
p
0.8
k
r
d
0.3Considering the bandwidths of interest of the current loop and the position loop of the MSFW and the simplicity of the realization, second-order antialias filters are introduced into both the displacement and current AD sampling of the MSFW system, whose transfer function is given as(39)gfs=13.3×10-5s+12.As far as the realization and implementation in hardware are concerned, a TMS320C32 digital signal processor (DSP) is employed. Both the sampling time and servo time are set to 150μs.


### 4.2. Simulation Results

Figure4 shows a plot of the nondimensional damping versus nondimensional stiffness for different growth rate. Gyroscopic effects result in either nutation (forward whirl) or precession (backward whirl). It also shows the growth rate, α (real part of the eigenvalue), and the frequency of oscillation, β (imaginary part of the eigenvalue). The nondimensional stiffness is given by (16) and the nondimensional damping is given by (15). If the real part of the eigenvalue α is defined as positive, the vibrations grow in time and the system is unstable. If α is defined as negative, the vibrations decay in time and the system is stable. The condition that α equals zero defines the stability boundary. The imaginary part of the eigenvalue, β, is a parameter which varies from 0 to π/2. The nondimensional speed of the flywheel is 0.1 for nutation, 0 for zero rotor speed, and −0.1 for precession. Nondimensional negative stiffness of 0.2 is assumed. From Figure 4, it can be drawn that a stable region is similar to that of the single-DOF magnetic bearing except it is now limited by nutation.Figure 4
Stable region for the MSFW using a centralized controller with time delay.Figure5 shows the dynamic stability map for a MSFW with time delay for various values of rotor speed, Ω. Note that the stiffness must be greater than the negative stiffness for the system to be stable. The area from the trajectory to the left line (τ02=0.2) represents the stability limits resulting from the forward whirl (the real part of the eigenvalue, σ, equals zero). The stable region becomes smaller for higher nondimensional rotor speeds. As for the traditional decentralized PD control, there is no stable region for nondimensional rotor speeds higher than 0.9. Contrarily, as for the decentralized PD plus cross feedback control, the nondimensional rotor speeds can reach about 1.7, and there is larger stable region than that of the case without cross feedback control at every rotor speed. This demonstrates that the cross feedback control algorithm can effectively improve the system stability, which is in accordance with the analysis made above.Figure 5
Dynamic stability map for a MSFW with time delay for various values of speed with different controller. (a) Decentralized PD controller. (b) Decentralized PD control plus cross feedback control.
(a)
(b)To further prove the correctness of Theorem1, that is, the phase compensation method, the positive and negative frequency characteristics analysis method [33] is employed. At Ω=40Hz, three different phase compensation methods have been compared, that is, the decentralized channels compensation, the cross-channel compensation, and the combinational compensation. To be fair, the three methods employ the same compensation filters: (40)gcs=2.2s2+3080s+1.1×1050.2s2+2000s+1.13×105.The control channel transfer function of the MSFW system without phase compensation can be described as(41)G0s=gasgfsgbs+jgcrse-τs.According to the positive and negative frequency characteristics stability criterion [33], the control channel transfer function with the decentralized channel can be described as(42)G1s=gasgfsgbsgcs+jgcrse-τs.Similarly, the control channel transfer function with the cross-channel compensation can be described as(43)G2s=gasgfsgbs+jgcrsgcse-τs.And the control channel transfer function with the combinational phase compensation can be described as(44)G2s=gasgfsgbs+jgcrsgcse-τs.Figure5 shows the comparative simulation results among different compensation methods, where “without compensation” demonstrates that the control system does not introduce any phase compensation. According to the nutation and precession stability criterion of the magnetically suspended with strong gyroscopic effects [33], at rotor speed Ω, the rotor speed stability margin Ωr can be described as(45)Ωr=JrωnJz+2khlm2-2lmlskiksinJzωn-Ω,where in is the gain of the control channel at the frequency ωn.From Figure6, the critical whirling frequencies and their corresponding gain can be obtained: ωn0=916rad/s, ωn1=1080rad/s, ωn2=1140rad/s, ωn3=1250rad/s, in0=0.0931, in1=0.162, i2=0.0812, and i3=0.161. According to the nutation stability criterion [33], it can be resolved that Ωr0=22.0Hz, Ωr1=7.0Hz, Ωr2=68.5Hz, and Ωr3=41.6Hz. It is obvious that, among the three phase compensation methods, the cross-channel phase compensation improves the system nutation stability most greatly. It is necessary to note that the traditional combinational compensation method does have some effect in improving the system nutation stability. However, it is not the most effective way compared with the cross-channel compensation. Additionally, the decentralized channel compensation method endangers the system nutation stability margin instead of improving it. All in all, among the three compensation methods, the cross-channel compensation is the most effective way, which is in accordance with the analysis made above.Figure 6
Comparative positive frequency characteristics curves among different compensation methods.


### 4.3. Experimental Results

To further prove the effectiveness and superiority of the phase compensation method, a comparison of experiments between the cross-channel compensation and the decentralized channels compensation has been performed. AtΩ=40Hz, these two different compensation methods are performed. Figures 7 and 8 show the relationship between the nutation amplitude and the rotor speed with the two different compensation methods, respectively.Figure 7
Nutation amplitude with rotor speed by the decentralized channel compensation.Figure 8
Nutation amplitude with temperature by the cross-channel compensation.As shown in Figure7, the nutation amplitude increases heavily with the rotor speed with the decentralized compensation method. And it increases to about −49 dB at 45 Hz, which demonstrates that the system will lose its nutation stability if the rotor speed continues to increase. On the contrary, the nutation amplitude keeps at about −54 dB in the rotor speed of 40 to 70 Hz with the cross-channel compensation method as shown in Figure 8. This shows that the system nutation stability has been improved greatly by the cross-channel compensation instead of the decentralized channel compensation.


## 4.1. Simulation and Experimental Setup

The large magnetically suspended control moment gyroscope (MSCMG), as shown in Figure3, is employed as the simulation and experimental setup since its MSFW has strong gyroscopic effects.Figure 3
Photograph of the fabricated laboratory setup.The well-known decentralized PD plus filter cross feedback controller is employed. The main parameters of the MSCMG and the controllers used in simulation and experiments have been listed in Tables1 and 2. In Table 1, x0 is the air gap of the radial magnetic bearings; ki0 and kh0 are the current stiffness and displacement stiffness of the radial magnetic bearings under the conditions iγ=0 and hmγ=0(γ=ax,ay,bx,by); and R and L are the coil resistance and inductance of the radial magnetic bearings. In Table 2, kam and ic are the proportional and feedback coefficients of the current controller for the MSFW system, which are the same as those of the proposed controller; kp and kd are the proportional, integral, and differential coefficients of the PD plus cross feedback controller for the MSFW system; khc and klc are the cross coefficients of the high-pass and low-pass filters, respectively; kc is the total cross coefficient; and fh and fl are the cut-off frequencies of the two-order high-pass and low-pass filters.Table 1
System parameters of the MSCMG.
Parameter
Value
m
56 kg
J
x
0.6366 kg⋅m2
J
y
0.6366 kg⋅m2
J
z
0.7958 kg⋅m2
k
i
0
1140 N/A
x
0
250μm
L
60.5 mH
R
2.5Ω
R
r
12000 r/min
l
s
0.177 m
l
m
0.113 m
k
h
0
2.83 N/μmTable 2
Parameters of the controller.
Parameter
Value
k
a
m
2.0
k
p
2.5
k
h
c
0.0045
k
c
1.0
f
l
30 Hz
k
r
i
0.04
k
c
p
30
i
c
10
k
d
0.008
k
l
c
0.00015
f
h
300 Hz
k
r
p
0.8
k
r
d
0.3Considering the bandwidths of interest of the current loop and the position loop of the MSFW and the simplicity of the realization, second-order antialias filters are introduced into both the displacement and current AD sampling of the MSFW system, whose transfer function is given as(39)gfs=13.3×10-5s+12.As far as the realization and implementation in hardware are concerned, a TMS320C32 digital signal processor (DSP) is employed. Both the sampling time and servo time are set to 150μs.


## 4.2. Simulation Results

Figure4 shows a plot of the nondimensional damping versus nondimensional stiffness for different growth rate. Gyroscopic effects result in either nutation (forward whirl) or precession (backward whirl). It also shows the growth rate, α (real part of the eigenvalue), and the frequency of oscillation, β (imaginary part of the eigenvalue). The nondimensional stiffness is given by (16) and the nondimensional damping is given by (15). If the real part of the eigenvalue α is defined as positive, the vibrations grow in time and the system is unstable. If α is defined as negative, the vibrations decay in time and the system is stable. The condition that α equals zero defines the stability boundary. The imaginary part of the eigenvalue, β, is a parameter which varies from 0 to π/2. The nondimensional speed of the flywheel is 0.1 for nutation, 0 for zero rotor speed, and −0.1 for precession. Nondimensional negative stiffness of 0.2 is assumed. From Figure 4, it can be drawn that a stable region is similar to that of the single-DOF magnetic bearing except it is now limited by nutation.Figure 4
Stable region for the MSFW using a centralized controller with time delay.Figure5 shows the dynamic stability map for a MSFW with time delay for various values of rotor speed, Ω. Note that the stiffness must be greater than the negative stiffness for the system to be stable. The area from the trajectory to the left line (τ02=0.2) represents the stability limits resulting from the forward whirl (the real part of the eigenvalue, σ, equals zero). The stable region becomes smaller for higher nondimensional rotor speeds. As for the traditional decentralized PD control, there is no stable region for nondimensional rotor speeds higher than 0.9. Contrarily, as for the decentralized PD plus cross feedback control, the nondimensional rotor speeds can reach about 1.7, and there is larger stable region than that of the case without cross feedback control at every rotor speed. This demonstrates that the cross feedback control algorithm can effectively improve the system stability, which is in accordance with the analysis made above.Figure 5
Dynamic stability map for a MSFW with time delay for various values of speed with different controller. (a) Decentralized PD controller. (b) Decentralized PD control plus cross feedback control.
(a)
(b)To further prove the correctness of Theorem1, that is, the phase compensation method, the positive and negative frequency characteristics analysis method [33] is employed. At Ω=40Hz, three different phase compensation methods have been compared, that is, the decentralized channels compensation, the cross-channel compensation, and the combinational compensation. To be fair, the three methods employ the same compensation filters: (40)gcs=2.2s2+3080s+1.1×1050.2s2+2000s+1.13×105.The control channel transfer function of the MSFW system without phase compensation can be described as(41)G0s=gasgfsgbs+jgcrse-τs.According to the positive and negative frequency characteristics stability criterion [33], the control channel transfer function with the decentralized channel can be described as(42)G1s=gasgfsgbsgcs+jgcrse-τs.Similarly, the control channel transfer function with the cross-channel compensation can be described as(43)G2s=gasgfsgbs+jgcrsgcse-τs.And the control channel transfer function with the combinational phase compensation can be described as(44)G2s=gasgfsgbs+jgcrsgcse-τs.Figure5 shows the comparative simulation results among different compensation methods, where “without compensation” demonstrates that the control system does not introduce any phase compensation. According to the nutation and precession stability criterion of the magnetically suspended with strong gyroscopic effects [33], at rotor speed Ω, the rotor speed stability margin Ωr can be described as(45)Ωr=JrωnJz+2khlm2-2lmlskiksinJzωn-Ω,where in is the gain of the control channel at the frequency ωn.From Figure6, the critical whirling frequencies and their corresponding gain can be obtained: ωn0=916rad/s, ωn1=1080rad/s, ωn2=1140rad/s, ωn3=1250rad/s, in0=0.0931, in1=0.162, i2=0.0812, and i3=0.161. According to the nutation stability criterion [33], it can be resolved that Ωr0=22.0Hz, Ωr1=7.0Hz, Ωr2=68.5Hz, and Ωr3=41.6Hz. It is obvious that, among the three phase compensation methods, the cross-channel phase compensation improves the system nutation stability most greatly. It is necessary to note that the traditional combinational compensation method does have some effect in improving the system nutation stability. However, it is not the most effective way compared with the cross-channel compensation. Additionally, the decentralized channel compensation method endangers the system nutation stability margin instead of improving it. All in all, among the three compensation methods, the cross-channel compensation is the most effective way, which is in accordance with the analysis made above.Figure 6
Comparative positive frequency characteristics curves among different compensation methods.


## 4.3. Experimental Results

To further prove the effectiveness and superiority of the phase compensation method, a comparison of experiments between the cross-channel compensation and the decentralized channels compensation has been performed. AtΩ=40Hz, these two different compensation methods are performed. Figures 7 and 8 show the relationship between the nutation amplitude and the rotor speed with the two different compensation methods, respectively.Figure 7
Nutation amplitude with rotor speed by the decentralized channel compensation.Figure 8
Nutation amplitude with temperature by the cross-channel compensation.As shown in Figure7, the nutation amplitude increases heavily with the rotor speed with the decentralized compensation method. And it increases to about −49 dB at 45 Hz, which demonstrates that the system will lose its nutation stability if the rotor speed continues to increase. On the contrary, the nutation amplitude keeps at about −54 dB in the rotor speed of 40 to 70 Hz with the cross-channel compensation method as shown in Figure 8. This shows that the system nutation stability has been improved greatly by the cross-channel compensation instead of the decentralized channel compensation.


## 5. Conclusion

This paper studied the effects of time delay on the stability of the rotation modes and presented an effective phase compensation method based on the cross-channels for the MSFW with strong gyroscopic effects. Analysis, simulation, and experimental results yield the following conclusion:(1)
The larger the rotor speed is, the smaller the stability region is.(2)
The forward whirl limited the region of stability.(3)
As for high speed MSFW, cross feedback control can effectively reject the nutation and precession instability.(4)
The delay of the cross-channel is the main factor influencing the nutation stability of the high-speed MSFW system, and the cross-channel compensation is the most effective way to improve the nutation stability.To sum up, the presented local stability analysis method can effectively demonstrate the system stability region of the MSFW with strong gyroscopic effects and time delay, and the presented phase compensation method can greatly improve the system nutation stability. The global stability and semistable oscillations are the next topics in the near future.



---

*Source: 3783740-2016-03-10.xml*