Application of active magnetic bearings for in situ flexible rotor residual balancing using a novel generalized influence coefficient method

ABSTRACT

An identification algorithm for the estimation of residual unbalances is developed for a flexible rotor-bearing system integrated with an active magnetic bearing (AMB). The main advantage of the present methodology is that it can be applied in situ and also with safe traversal of resonances during run-up. The flexible rotor-AMB system is modelled by the finite element method, considering randomly distributed residual unbalances along with disc discrete residual unbalances. In order to obtain equivalent discrete unbalances at limited disc locations for the balancing of the flexible rotor system, a generalized influence coefficient method (GICM) is proposed. Equivalent unbalances obtained at discrete balancing planes are used to check the reduction in unbalance responses. The proposed algorithm is compared with the conventional influence coefficient method (CICM) without AMB to check the accuracy and consistency of the estimated values of residual unbalances. The values estimated by both methods are found to be consistent, while the GICM has an advantage of using AMB for the safe operation of the unbalanced rotor at high speeds during measurements. In order to have a suitable discrete position of AMB (for effective estimation of residual unbalances), various rotor system models are considered with different positions of the AMB.

KEYWORDS:

• Flexible rotor
• balancing
• active magnetic bearing
• generalized influence coefficient method
• residual unbalance

MSC SUBJECT CLASSIFICATIONS:

• 74A99
• 70E99
• 70E55
• 70J35

1. Introduction

Nowadays, in industries, the use of high-speed rotating machines has increased. So, there is a need for improving the characteristics and performance of such rotating systems. The problems that arise in these machines are the unbalances, misalignments, fatigue cracks, and wear and tear. The major problems are the residual unbalances, which cause excessive vibration in machines and lead to fatigue failures. The residual unbalances (inherent unbalances) in machineries occur due to manufacturing errors, improper commissioning, wear and tear, and so on. There is a need to actively control the rotor system vibration and it can be done effectively by the use of active magnetic bearing (AMB), which provides contactless actuation and helps in reduction of unnecessary vibrations during traversing of critical speeds or unstable regions [1–3]. However, the design of AMB is complex as it involves the design of controllers, actuators, amplifiers, filters and sensors [4–6]. AMB reduces excessive vibration by applying magnetic forces in opposite directions to unbalance forces. Nevertheless, instead of applying this corrective force all the time, if these corrective forces from AMBs can be used for finding the residual unbalance states of the rotor that will save unnecessary power consumption of AMBs, especially when rotors are still mounted on conventional bearings.

To reduce excessive vibration of rotors mounted on conventional bearings, many balancing methods have been developed by researchers, notably the conventional influence coefficient method (CICM) and modal balancing method [1]. Kang et al. [7] derived influence coefficient matrices from equations of motion for unsymmetrical rotors by using a complex coordinate representation and the finite element method (FEM). They predicted that the modified influence coefficients are found to be correlated to unbalance forces. Ehyaei et al. [8] investigated a system of unbalanced flexible shaft equipped with ‘n’ automatic ball balancers. They concluded that if the auto ball balancers are closer to the unbalanced masses distributed in the shaft, a better vibration reduction is achieved. Lee et al. [9] studied an electromagnetic type active balancing device and a balancing programme using the influence coefficient method for both simple and reliable structures. Sudhakar and Sekhar [10] identified faults in the rotor system by minimizing the difference between equivalent loads estimated in the system. Pennacchi et al. [11] balanced the rotor system using high breakdown point and bounded influence estimators.

Several of the research works have been done to determine the magnetic forces and stiffness parameters of AMBs. Autila et al. [12] predicted the performance and parameters of a linearized model of radial magnetic bearings by using a nonlinear two-dimensional FEM. Tiwari and Chakravarthy [13] estimated bearing dynamic parameters and residual unbalances simultaneously using impulse response measurements for a multi-degree-of-freedom flexible rotor system. Imlach et al. [14] estimated force and stiffness of radial AMBs using the linear magnetic circuit theory. Zhou et al. [15] obtained the stiffness and damping coefficients of an AMB by using the magnetic bearing as a sinusoidal exciter. Tiwari and Chougale [16] developed an identification algorithm for estimating speed-dependent parameters of AMB using flexible rotors fully levitated on AMBs. They further validated their work with the experimental data.

Lots of researches have been carried out in the design of the AMB controller for the rotor system. Various control principles have been applied to the magnetic bearing from the traditional proportional-derivative controller to the contemporary advanced adaptive ones [17–19]. Binder et al. [20] were involved in the modelling, simulation, and control of AMB system by using the state-space formulation and PID controllers. Behzad et al. [21] designed and implemented both linear and nonlinear control methods for a three-pole AMB. For a high-speed rotating system, the AMB can also be used to determine the balancing masses to physically balance the system as it reduces the controller effort and power consumed by the controlling system. In literature, the AMB has been used as the bearing support (for the static load of a rotor), controller (for dynamic internal disturbances) and exciter (for the system response independent dynamic external excitation), individually or in combinations. The literatures on the use of AMB for condition monitoring are rare. Conventionally, the AMB is used as an exciter to give a known excitation and based on the input and output information the system fault parameters are estimated. AMBs have also been utilized for fully levitating rotors for very high-speed applications but still, several applications require rotors supported mainly on conventional bearings. Especially for flexible rotor applications and balancing them in situ is still a great challenge since it requires measurements up to operating speed before balancing the rotors and herein the AMB as a controller can help in such dynamic balancing of flexible rotors.

In the present work, the AMB is used as a controller in a flexible rotor-bearing system and the controlling force and associated suppressed rotor responses are used to find the unbalance state of the rotor. The advantage of it is to suppress the vibration during measurement even at high speeds and then this suppression force along with responses are used to get the estimate of residual unbalances. Section 2 describes the rotor system and finite element formulation for the rotor-bearing-AMB system. Section 3 illustrates the formulation of CICM and proposed generalized influence coefficient method (GICM) for the balancing of flexible rotors. Section 4 describes the numerical simulation involved in the prediction of unbalances and balancing of the rotor system. Section 5 shows the estimation of unbalances with different rotor systems and selection of the appropriate position of AMB in the rotor system. Section 4 concludes the present work.

2. Rotor-bearing-AMB model descriptions

A flexible rotor-bearing system integrated with AMBs, as shown in Figure 1, has been considered for the study. It consists of rigid discs mounted on flexible shafts that are connected by couplings. The AMBs are used as controlling devices to prevent excessive vibrations at high speeds. Flexible shafts are supported on conventional bearings with cross-coupled stiffness and damping terms. AMBs have linearized cross-coupled displacement and current stiffnesses. In this rotor system, discs are assumed to be balancing planes and residual unbalances are estimated in these planes. A proportional–integral–derivative (PID) controller is used for AMBs.

Figure 1. A schematic of a flexible rotor-bearing system integrated with active magnetic bearings.

2.1. Finite element formulation of the rotor-AMB system

In the FEM, the rotor system is discretized into smaller elements (Figure 2) to simplify the formulation. Equations of motion of each element are considered separately and reassembled to get the global equations of motion of the complete assembly. Then boundary conditions are applied to the global equations to get the desired governing equations. In this model, unbalance forces (both shaft and disc unbalances), inertia forces, gyroscopic effect (discs only), AMB forces, and the stiffness and damping characteristics of conventional bearings are considered. Three orthogonal axes have been selected, i.e. z-axis along the shaft axis, x-axis and y-axes are perpendicular to the shaft axis. Shaft elements are having four degrees-of-freedom at each node (two transverse translational and two transverse rotational).

Figure 2. An element of shaft showing different segments having unbalance forces.

Equations of motion for the rotor system are given as (1) $\mathbf{\text{M}}\ddot{\mathit{\eta }}\left(t\right)+\left(\mathbf{\text{C}}-\omega \mathbf{\text{G}}\right)\dot{\mathit{\eta }}\left(t\right)+\mathbf{\text{K}}\mathit{\eta }\left(t\right)=\mathbf{\text{f}}\left(t\right)$(1) with (2) $\mathbf{\text{f}}\left(t\right)\mathbf{\text{=}}{\mathbf{\text{f}}}_{\mathbf{\text{unb}}}\left(t\right)-{\mathbf{\text{f}}}_{\mathbf{\text{mb}}}\left(t\right),$(2) where ω is the spin speed of the rotor, M, C, G and K represent the mass, damping, gyroscopic and stiffness matrices, respectively [2]. Details of the force vector, f, and the displacement vector, η, have been provided in the Appendix. The above equation is rearranged as follows, to obtain the displacement vector from the SIMULINK model. (3) $\ddot{\mathit{\eta }}={\mathbf{\text{M}}}^{-1}\left(\mathbf{\text{f}}-\mathbf{\text{C}}\dot{\mathit{\eta }}+\omega \mathbf{\text{G}}\dot{\mathit{\eta }}-\mathbf{\text{K}}\mathit{\eta }\right).$(3) The displacement and magnetic force vectors obtained from the SIMULINK software have been used in GICM that is explained next to obtain equivalent residual balance masses of the flexible rotor-bearing system.

3. Formulation of the conventional and GICMs

Two different types of technique are developed for the purpose of balancing of the flexible rotor system. The first one is modal balancing technique developed by Bishop and Gladwell [23] and the second is influence coefficient method developed by Drechler [24]. The CICM works well if the rotor is mounted on only conventional bearings. In the system having AMB, conventional methods cannot be used as the responses produced by trial masses may not directly reflect because of the AMB. In the present work, an AMB is used as controlling device, so in order to obtain the balance masses at disc locations for the balancing, the flexible rotor system a GICM has been proposed. The details of CICM and GICM have been discussed in next subsections.

3.1. Conventional influence coefficient method

The CICM is widely used approach for balancing of rotor systems in industrial applications. The balancing procedure is easy to carry out with the help of computer codes and requires very little knowledge of the rotating machineries and physical-mathematical modelling. The CICM is used to predict residual unbalances in a rigid as well as flexible rotor system with the use of displacement values at different measuring planes. Displacements values are obtained for different trial unbalances at balancing planes for estimation of influence coefficients. Then, these are used in obtaining residual unbalances at balancing planes.

In the CICM, n numbers of balancing planes are selected for trial unbalances and p numbers of measurement planes are considered for the measurement of displacements of the rotor system. Apart from this, a common reference signal is considered to have phases of all measured signals. The displacement vector, d(ω), and unbalance vector, w, are related with influence coefficient matrix, B(ω), for a particular rotor speed, as (4) $\mathbf{d}\left(\omega \right)\mathbf{=}\mathbf{B}\left(\omega \right)\mathbf{w}$(4) With $\begin{array}{rl}\mathbf{w}& ={\left\{\begin{array}{c}{w}_1\\ w_2\\ ⋮\\ w_n\end{array}\right\}}_{n\times 1};\mathbf{d}\left(\omega \right)={\left\{\begin{array}{c}{x}_1\left(\omega \right)\\ y_1\left(\omega \right)\\ x_2\left(\omega \right)\\ y_2\left(\omega \right)\\ ⋮\\ x_p\left(\omega \right)\\ y_p\left(\omega \right)\end{array}\right\}}_{\left(2p\right)\times 1};\\ \mathbf{B}\left(\omega \right)& ={\left[\begin{array}{cccc}{\beta }_{x_1,1}\left(\omega \right)& {\beta }_{x_1,2}\left(\omega \right)& \dots & {\beta }_{x_1,n}\left(\omega \right)\\ {\beta }_{y_1,1}\left(\omega \right)& {\beta }_{y_1,2}\left(\omega \right)& \dots & {\beta }_{y_1,n}\left(\omega \right)\\ {\beta }_{x_2,1}\left(\omega \right)& {\beta }_{x_2,2}\left(\omega \right)& \dots & {\beta }_{x_2,n}\left(\omega \right)\\ {\beta }_{y_2,1}\left(\omega \right)& {\beta }_{y_2,2}\left(\omega \right)& \dots & {\beta }_{y_2,n}\left(\omega \right)\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}\left(\omega \right)& {\beta }_{x_p,2}\left(\omega \right)& \dots & {\beta }_{x_p,n}\left(\omega \right)\\ {\beta }_{y_p,1}\left(\omega \right)& {\beta }_{y_p,2}\left(\omega \right)& \dots & {\beta }_{y_p,n}\left(\omega \right)\end{array}\right]}_{\left(2p\right)\times n}\end{array}$ $x_p\left(\omega \right)=R_{xp}{e}^{\mathrm{j}{\phi }_{xp}}$ and $y_p\left(\omega \right)=R_{yp}{e}^{\mathrm{j}{\phi }_{yp}}$

where w₁, w₂, … , w_n are residual unbalances at the n number of balancing planes and these are independent of the rotor speed. R_xpand R_yp are amplitudes of measured displacements in the x and y-axis directions, respectively; φ_xpand φ_yp are phase angles of the displacements in the x and y-axis directions, respectively, with respect to a fixed point of the shaft; and $j=\sqrt{-1}$. The first column of influence coefficient matrix: ${\beta }_{x_1,1}\left(\omega \right)$, ${\beta }_{y_1,1}\left(\omega \right)$, … , ${\beta }_{x_p,1}\left(\omega \right)$ and ${\beta }_{y_p,1}\left(\omega \right)$ are defined as the displacements at different measuring planes in x and y directions, respectively, per unit residual unbalance (w₁). Similarly, other columns of B(ω) can be defined corresponding to different residual unbalances. The influence coefficients matrix, B(ω) can be obtained considering different trial unbalance in each balancing plane at a time and estimating displacements at various measuring planes as (5) $\left\{\begin{array}{c}{x}_{11}\\ y_{11}\\ x_{21}\\ y_{21}\\ ⋮\\ x_{p1}\\ y_{p1}\end{array}\right\}=\left[\begin{array}{cccc}{\beta }_{x_1,1}& {\beta }_{x_1,2}& \dots & {\beta }_{x_1,n}\\ {\beta }_{y_1,1}& {\beta }_{y_1,2}& \dots & {\beta }_{y_1,n}\\ {\beta }_{x_2,1}& {\beta }_{x_2,2}& \dots & {\beta }_{x_2,n}\\ {\beta }_{y_2,1}& {\beta }_{y_2,2}& \dots & {\beta }_{y_2,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}& {\beta }_{x_p,2}& \dots & {\beta }_{x_p,n}\\ {\beta }_{y_p,1}& {\beta }_{y_p,2}& \dots & {\beta }_{y_p,n}\end{array}\right]\left\{\begin{array}{c}{w}_1+T_R\\ w_2\\ ⋮\\ w_n\end{array}\right\}$(5) with $T_R=m_1{e}_1{e}^{j{\phi }_1},$ where T_R is the trial unbalance added at one of the measuring plane, and m₁ is the trial mass, φ₁ is corresponding phase and e₁ is the eccentricity. The second subscript in displacements (x₁₁, y₁₁, … , x_p1, y_p1) represents corresponding to balancing plane in which trial unbalance is kept. It is assumed that the influence coefficient does not change by adding small trial unbalance. On subtracting Equation (4) from Equation (5), we get (6) $\left\{\begin{array}{c}{x}_{11}-x_1\\ y_{11}-y_1\\ x_{21}-x_2\\ y_{21}-y_2\\ ⋮\\ x_{p1}-x_p\\ y_{p1}-y_p\end{array}\right\}=\left[\begin{array}{cccc}{\beta }_{x_1,1}& {\beta }_{x_1,2}& \dots & {\beta }_{x_1,n}\\ {\beta }_{y_1,1}& {\beta }_{y_1,2}& \dots & {\beta }_{y_1,n}\\ {\beta }_{x_2,1}& {\beta }_{x_2,2}& \dots & {\beta }_{x_2,n}\\ {\beta }_{y_2,1}& {\beta }_{y_2,2}& \dots & {\beta }_{y_2,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}& {\beta }_{x_p,2}& \dots & {\beta }_{x_p,n}\\ {\beta }_{y_p,1}& {\beta }_{y_p,2}& \dots & {\beta }_{y_p,n}\end{array}\right]\left\{\begin{array}{c}{T}_R\\ 0\\ ⋮\\ 0\end{array}\right\}$(6) which gives the first column of the influence coefficients as (7) ${\beta }_{x_1,1}=\frac{x_{11}-x_1}{T_R},{\beta }_{y_1,1}=\frac{y_{11}-y_1}{T_R},\dots ,{\beta }_{x_p,1}=\frac{x_{p1}-x_p}{T_R}\text{and}{\beta }_{y_p,1}=\frac{y_{p1}-y_p}{T_R}.$(7)

Similarly, the second column of the influence coefficient matrix can be obtained by adding a trial unbalance on second balancing plane, and similarly other column influence coefficients can be obtained. For flexible rotors influence coefficients change with rotor speed. For different speeds in a particular run-up, above equations will take the form, as (8) $\mathbf{d}\mathbf{=}\mathbf{B}\mathbf{w}$(8)

with $\mathbf{d}\mathbf{=}{\left\{\begin{array}{c}\mathbf{d}\left({\omega }_1\right)\\ \mathbf{d}\left({\omega }_2\right)\\ ⋮\\ \mathbf{d}\left({\omega }_f\right)\end{array}\right\}}_{1\times \left(2pf\right)}^{};\mathbf{B}\mathbf{=}{\left[\begin{array}{c}\mathbf{B}\left({\omega }_1\right)\\ \mathbf{B}\left({\omega }_2\right)\\ ⋮\\ \mathbf{B}\left({\omega }_f\right)\end{array}\right]}_{n\times \left(2pf\right)}^{}\text{and}\mathbf{w}={\left\{\begin{array}{c}{w}_1\\ w_2\\ ⋮\\ w_n\end{array}\right\}}_{n\times 1}$ where the superscript in ω represents a particular speed of the rotor. The vectors d(ω₁), d(ω₂), … , d(ω_f) are the complex displacements corresponding to different speeds obtained from FFT of the time domain data (the procedure in discussed subsequently). The influence coefficients can be obtained for different speed by taking different trial unbalances as shown in Equations (4)–(6) for a single speed. Residual unbalances in the rotor system are obtained as (9) $\mathbf{w}=({\mathbf{B}}^T\mathbf{B}{)}^{-1}{\mathbf{B}}^T\mathbf{d}\mathbf{.}$(9)

Now the extension of the CICM for the rotor-bearing-AMB system will be presented.

3.2. Generalized influence coefficient method

The GICM has been proposed to predict residual unbalances in a rotor-AMB system with the use of displacement values at various measuring planes and magnetic force values of AMBs, while rotors are still supported on conventional bearings. Displacements and magnetic forces are obtained for different trial unbalances at disc locations (considered as balancing planes) for the estimation of generalized influence coefficients. Then, these are used in obtaining residual balance masses at balancing planes. Equations involved in the formulation of GICM are described now.

Displacement vectors corresponding to p number of measuring planes and magnetic forces corresponding to q number of AMBs are used to obtain residual unbalances at n number of balancing planes. In fact, magnetic force vectors contain both displacement and current at AMB locations. The complex displacement vector taken at different measuring planes in the x- and y-axis directions is shown in Equation (4). The complex magnetic force vectors at AMB locations are represented as (10) $\mathbf{f}\left(\omega \right)={\left\{\begin{array}{c}{f}_{x_1}\left(\omega \right)\\ f_{y_1}\left(\omega \right)\\ f_{x_2}\left(\omega \right)\\ f_{y_2}\left(\omega \right)\\ ⋮\\ f_{x_q}\left(\omega \right)\\ f_{y_q}\left(\omega \right)\end{array}\right\}}_{\left(2q\right)\times 1}$(10) The complex magnetic force vector, f(ω), the displacement vector, d(ω), and unbalance vector, w, are related with generalized influence coefficient matrix, G(ω), for a particular rotor speed, as (11) $\mathbf{h}\left(\omega \right)\mathbf{=}\mathbf{G}\left(\omega \right)\mathbf{w}$(11) with $\mathbf{w}={\left\{\begin{array}{c}{w}_1\\ w_2\\ ⋮\\ w_n\end{array}\right\}}_{n\times 1};\mathbf{h}\left(\omega \right)={\left\{\begin{array}{c}\mathbf{d}\left(\omega \right)\\ \mathbf{f}\left(\omega \right)\end{array}\right\}}_{2(p+q)\times 1};\mathbf{G}\left(\omega \right)={\left[\begin{array}{cccc}{\beta }_{x_1,1}\left(\omega \right)& {\beta }_{x_1,2}\left(\omega \right)& \dots & {\beta }_{x_1,n}\left(\omega \right)\\ {\beta }_{y_1,1}\left(\omega \right)& {\beta }_{y_1,2}\left(\omega \right)& \dots & {\beta }_{y_1,n}\left(\omega \right)\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}\left(\omega \right)& {\beta }_{x_p,2}\left(\omega \right)& \dots & {\beta }_{x_p,n}\left(\omega \right)\\ {\beta }_{y_p,1}\left(\omega \right)& {\beta }_{y_p,2}\left(\omega \right)& \dots & {\beta }_{y_p,n}\left(\omega \right)\\ {\gamma }_{x_1,1}\left(\omega \right)& {\gamma }_{x_1,2}\left(\omega \right)& \dots & {\gamma }_{x_1,n}\left(\omega \right)\\ {\gamma }_{y_1,1}\left(\omega \right)& {\gamma }_{y_1,2}\left(\omega \right)& \dots & {\gamma }_{y_1,n}\left(\omega \right)\\ ⋮& ⋮& \ddots & ⋮\\ {\gamma }_{x_q,1}\left(\omega \right)& {\gamma }_{x_q,2}\left(\omega \right)& \dots & {\gamma }_{x_q,n}\left(\omega \right)\\ {\gamma }_{y_q,1}\left(\omega \right)& {\gamma }_{y_q,2}\left(\omega \right)& \dots & {\gamma }_{y_q,n}\left(\omega \right)\end{array}\right]}_{2(p+q)\times n}$ where w₁, w₂, … , w_n are residual unbalances at the n number of balancing planes. Herein, it should be noted that now d(ω) is the displacement vector in presence of AMB in the rotor system. The generalized influence coefficients matrix, G(ω) can be obtained considering different trial unbalance in each balancing plane at a time and measuring corresponding displacements at various measuring planes, as (12) $\left\{\begin{array}{c}{x}_{11}\\ y_{11}\\ ⋮\\ x_{p1}\\ y_{p1}\\ f_{x_{11}}\\ f_{y_{11}}\\ ⋮\\ f_{x_{q1}}\\ f_{y_{q1}}\end{array}\right\}=\left[\begin{array}{cccc}{\beta }_{x_1,1}& {\beta }_{x_1,2}& \dots & {\beta }_{x_1,n}\\ {\beta }_{y_1,1}& {\beta }_{y_1,2}& \dots & {\beta }_{y_1,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}& {\beta }_{x_p,2}& \dots & {\beta }_{x_p,n}\\ {\beta }_{y_p,1}& {\beta }_{y_p,2}& \dots & {\beta }_{y_p,n}\\ {\gamma }_{x_1,1}& {\gamma }_{x_1,2}& \dots & {\gamma }_{x_1,n}\\ {\gamma }_{y_1,1}& {\gamma }_{y_1,2}& \dots & {\gamma }_{y_1,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\gamma }_{x_q,1}& {\gamma }_{x_q,2}& \dots & {\gamma }_{x_q,n}\\ {\gamma }_{y_q,1}& {\gamma }_{y_q,2}& \dots & {\gamma }_{y_q,n}\end{array}\right]\left\{\begin{array}{c}{w}_1+T_R\\ w_2\\ ⋮\\ w_n\end{array}\right\}$(12) where T_R is the trial unbalance added at one of the measuring planes as shown in Equation (5). Herein, the second subscript of displacements and magnetic forces represents the plane at which the trial mass was kept. On subtracting Equation (11) from Equation (12), we get (13) $\left\{\begin{array}{c}{x}_{11}-x_1\\ y_{11}-y_1\\ ⋮\\ x_{p1}-x_p\\ y_{p1}-y_p\\ f_{x_{11}}-f_{x_1}\\ f_{y_{11}}-f_{y_1}\\ ⋮\\ f_{x_{q1}}-f_{x_q}\\ f_{y_{q1}}-f_{y_q}\end{array}\right\}=\left[\begin{array}{cccc}{\beta }_{x_1,1}& {\beta }_{x_1,2}& \dots & {\beta }_{x_1,n}\\ {\beta }_{y_1,1}& {\beta }_{y_1,2}& \dots & {\beta }_{y_1,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\beta }_{x_p,1}& {\beta }_{x_p,2}& \dots & {\beta }_{x_p,n}\\ {\beta }_{y_p,1}& {\beta }_{y_p,2}& \dots & {\beta }_{y_p,n}\\ {\gamma }_{x_1,1}& {\gamma }_{x_1,2}& \dots & {\gamma }_{x_1,n}\\ {\gamma }_{y_1,1}& {\gamma }_{y_1,2}& \dots & {\gamma }_{y_1,n}\\ ⋮& ⋮& \ddots & ⋮\\ {\gamma }_{x_q,1}& {\gamma }_{x_q,2}& \dots & {\gamma }_{x_q,n}\\ {\gamma }_{y_q,1}& {\gamma }_{y_q,2}& \dots & {\gamma }_{y_q,n}\end{array}\right]\left\{\begin{array}{c}{T}_R\\ 0\\ ⋮\\ 0\end{array}\right\}$(13) which gives, (14) ${\beta }_{x_1,1}=\frac{x_{11}-x_1}{T_R},{\beta }_{y_1,1}=\frac{y_{11}-y_1}{T_R}\dots ,{\beta }_{x_p,1}=\frac{x_{p1}-x_p}{T_R}\text{and}{\beta }_{y_p,1}=\frac{y_{p1}-y_p}{T_R},$(14) (15) ${\mathit{\gamma }}_{x_1,1}=\frac{f_{x_{11}}-f_{x_1}}{T_R},{\mathit{\gamma }}_{y_1,1}=\frac{f_{y_{11}}-f_{y_1}}{T_R},\dots ,{\mathit{\gamma }}_{x_q,1}=\frac{f_{x_{q1}}-f_{x_q}}{T_R}\text{and}{\gamma }_{y_q,1}=\frac{f_{y_{q1}}-f_{y_q}}{T_R}.$(15) Similarly, the second column of the generalized influence coefficient matrix can be obtained by adding trial unbalance on next balancing plane. For different speeds in a particular run-up, above equations will take the form, as (16) $\mathbf{h}\mathbf{=}\mathbf{G}\mathbf{w}$(16) with $\mathbf{h}={\left\{\begin{array}{c}\mathbf{h}\left({\omega }_1\right)\\ \mathbf{h}\left({\omega }_2\right)\\ ⋮\\ \mathbf{h}\left({\omega }_f\right)\end{array}\right\}}_{2f(p+q)\times 1}$ and $\mathbf{G}\mathbf{=}{\left[\begin{array}{c}\mathbf{G}\left({\omega }_1\right)\\ \mathbf{G}\left({\omega }_2\right)\\ ⋮\\ \mathbf{G}\left({\omega }_f\right)\end{array}\right]}_{2f(p+q)\times n}$

where, h(ω₁), h(ω₂), … , h(ω_f) are the complex displacement and magnetic force vectors corresponding to different speeds obtained from FFT of the time domain data and magnetic force vectors. The generalized influence coefficient matrix G is obtained in a similar way as shown in Equations (12)–(15) considering different trial unbalances at the balancing planes. Residual unbalances in the rotor system are obtained as (17) $\mathbf{w}=({\mathbf{G}}^T\mathbf{G}{)}^{-1}{\mathbf{G}}^T\mathbf{h}\mathbf{.}$(17) The unbalances obtained in the previous equation are taken as balance masses with +180° phase and are added at the balancing planes to reduce the vibration of the rotor system.

4. Numerical simulations of balancing methodology

In the present study, the rotor system consists of three rigid discs (balancing planes), a flexible shaft, two conventional bearings (supports) and a magnetic bearing (controlling device) as shown in Figure 3(a). Initially, for the present rotor system (Model 1) the free vibration analysis has been done without considering AMB in the system to analyse the gyroscopic effect of the discs and predict the critical speeds of the rotor system. Equations of motion of the rotor system as given by Equation (1) are considered without any external forces for the free vibration analysis. Physical properties of the shaft, discs and conventional bearings [22] have been taken as mentioned in Tables 1 and 2. Eigen values are obtained for different values of the spin speed of the rotor system and Campbell diagram is plotted as shown in Figure 4. The plot shows that the rotor system crosses its fourth critical speed while rotating from 0 to 700 rad/s. It is also shown that a considerable split of natural whirl frequencies of the rotor system in presence of a gyroscopic effect in the range of operating speed. In each of splitting of natural whirl frequencies, the lower curve corresponds to the backward whirl and the upper curve corresponds to forward whirl in the plot. For the backward whirl, critical speeds obtained from Campbell diagram are 74.53, 197.48, 318.88, 537.19 rad/s and forward whirl are 77.49, 220.56, 369.85, 692.26 rad/s.

Figure 3. The rotor-bearing-AMB system (Model 1).

Figure 4. A Campbell diagram of the rotor-AMB system (Model 1).

Table 1. Physical properties of the rotor (Model 1).

Parameter	Numerical value
Mass density of shaft, ρ	7850 kg/m^3
Young’s modulus of elasticity, E	2.06 × 10^11 N/m^2
Diameter of shaft, d and total length of shaft, L	8 and 615 mm
Length of the elements in the rotor system	96, 139, 165, 126 and 89 mm
Disc 1: m (mass) and Id (diametral mass moment of inertia)	0.78 kg and 4.44 × 10^−4 kg m^2
Disc 2 and Disc 3: m and Id	1.20 kg and 13.1 × 10^−4 kg m^2
Magnetic bearing core: m and Id	0.88 kg and 5.17 × 10^−4 kg m^2
Residual unbalance with angular position at Disc 1	2.34 × 10^−5 kg m, 5.73°
Residual unbalance with angular position at Disc 2	7.20 × 10^−5 kg m, 68.75°
Residual unbalance with angular position at Disc 3	3.60 × 10^−5 kg m, 120.32°
Trial mass with angular position at Disc 1	6.00 × 10^−5 kg m, 63.03°
Trial mass with angular position at Disc 2	13.2 × 10^−5 kg m, 80.21°
Trial mass with angular position at Disc 3	8.40 × 10^−5 kg m, 22.92°

Table 2. Physical properties of conventional bearings.

Parameter	Numerical value
Conventional bearing 1
Stiffness in x and y-direction ($k_{x{x}_1}$,$k_{y{y}_1}$)	8 × 10^5 N/m, 4 × 10^5 N/m
Stiffness in x-y and y-x direction ($k_{x{y}_1}$,$k_{y{x}_1}$)	2 × 10^3 N/m, 2 × 10^3 N/m
Damping coefficient in x and y direction ($c_{x{x}_1}$,$c_{y{y}_1}$)	100 N-s/m, 200 N-s/m
Damping coefficient in x-y and y-x direction ($c_{x{y}_1}$,$c_{y{x}_1}$)	8 N-s/m, 9 N-s/m
Conventional bearing 2
Stiffness in x and y-direction ($k_{x{x}_2}$,$k_{y{y}_2}$)	7.5 × 10^5 N/m, 4.5 × 10^5 N/m
Stiffness in x-y and y-x direction ($k_{x{y}_2}$,$k_{y{x}_2}$)	2.1 × 10^3 N/m, 2.1 × 10^3 N/m
Damping coefficient in x and y direction ($c_{x{x}_2}$,$c_{y{y}_2}$)	110 N-s/m, 190 N-s/m
Damping coefficient in x-y and y-x direction ($c_{x{y}_2}$,$c_{y{x}_2}$)	7.5 N-s/m, 8 N-s/m

The proposed identification algorithm is tested with a flexible rotor system using a numerical simulation and detailed steps of it is illustrated with the help of the flowchart in Figure 5. The numerical rotor model is generated with the help of SIMULINK shown in Figure 6. The PID controller has following properties: K_P= 4200 A/m, K_D= 10 A-s/m and K_I= 8000 A/m-s [16]. In order to generate responses using the numerical simulation, properties of the conventional and magnetic bearings for all speeds are assumed and mentioned in Table 2 and Table 3. The speed is varied from 5 to 700 rad/s with an angular acceleration of 139 rad/s². The mass, stiffness and damping matrices, corresponding to Equation (1), are defined in MATLAB, which are then accessed by the SIMULINK model to generate unbalance forces (both of shaft and discs). Displacements and magnetic forces are then estimated in both orthogonal directions (x and y) at two conventional bearing locations and the AMB, respectively. In displacements and magnetic forces, 1% of noise (random values) has been added, to check the robustness of the algorithm against noise.

Figure 5. Flowchart of the flexible rotor balancing using generalized influence coefficients for the rotor-AMB system.

Figure 6. A block diagram of SIMULINK model for the rotor-AMB system.

Table 3. Physical properties of AMB.

Parameter	Numerical value
Proportionality constant for PID controller, kP	4200 A/m
Derivative constant for PID controller, kD	10 A-s/m
Integral constant for PID controller, kI	8000 A/m-s
Direct force current factor in x and y-direction ($k_{i_{xx}}$,$k_{i_{yy}}$)	42.10 N/A
Cross-coupled force current factor in x-y and y-x direction ($k_{i_{xy}}$,$k_{i_{yx}}$)	4.21 N/A
Direct force displacement factor in x and y-direction ($k_{s_{xx}}$,$k_{s_{yy}}$)	1.052 × 10^5 N/m
Cross-coupled force displacement factor in x-y and y-x direction ($k_{s_{xy}}$,$k_{s_{yx}}$)	0.526 × 10^4 N/m

The time domain run-up displacement (or current) data obtained is made even in size in order to get the data for the half-spectrum. The half of the FFT data from the end has been removed (according to Nyquist criterion) and the rest has been used for obtaining unbalance masses using the GICM or the CICM. But before that, the phases of different run-up displacement (or current) data have to be corrected with respect to a common reference signal. As shown in Figure 7, fast Fourier transform (FFT) of both displacements and reference signal are carried out to obtain the complex displacement and reference signal vectors, respectively, with its amplitude and phase. The phase of the reference signal is subtracted from the phase of displacement (or current) signal to obtain the required phase of the displacement (or current) signal for different run-ups. The same process is carried out with the magnetic force vector obtained at AMB location to obtain complex magnetic force vectors in both x and y directions. The complex displacement and magnetic force vectors are obtained for the rotor system with and without trial unbalances at the balancing planes. With the help of displacements and magnetic forces (which otherwise can be easily measured in practice), residual unbalance mass magnitudes and phases are calculated for the flexible rotor system using GICM (through Equations (10)–(17)). The estimated residual unbalances at the balancing planes are considered as balance masses with 180° phase. Balance masses are then added to the rotor system and displacements at different rotor nodes are obtained in both time domain and frequency domain without considering AMB. Displacements are also obtained without balance masses and AMB to compare the reduction in unbalance responses with the addition of balance masses.

Figure 7. A process chart for the calculation of correct phase in run-up responses using the FFT.

To check the algorithm, first only disc residual (discrete) unbalances are taken and with the help of displacements and magnetic forces, the unbalance masses magnitude and phase are calculated for the rotor system using the GICM (Equations (4)–(9)) as provided in Table 4. As shown in Table 4, the assumed and the estimated values of residual unbalances are found to be the same. Now, to check the reduction in unbalance responses, considering shaft residual (distributed) unbalances along with disc (discrete) unbalances, the magnitude and phase of equivalent unbalances at three balancing planes are obtained by the GICM as given in Table 5. The unbalances obtained have different values of magnitude and phase at the discs than mentioned in Table 4 because of the introduction of the shaft unbalances in the system. Unbalance responses are again generated by adding estimated unbalances obtained as balance masses with 180° phase to the system. The displacement amplitude and phase values of responses have been obtained using algorithm mentioned in Figure 7 in the horizontal direction, i.e. x-axis. The displacement amplitudes consist of only 1× harmonic at all speeds in the run-up condition of the rotor system. The displacement amplitude versus the rotor run-up speed plots have been plotted for the rotor system to check the reduction in unbalances after adding balance masses for rotor Model 1 at disc 1 and conventional bearing 1 positions and without using AMB as shown in Figures 8 and 9. The peaks in Figure 8(b) show the amplitude of vibration while the rotor crosses its different critical speeds. As shown in these plots, unbalance responses are found to be reduced at different critical speeds after balancing in Figures 8(b) and 9(a). The variation of the phase lag (phase of the responses with respect to the reference signal) of the unbalance response with run-up speed is shown in Figures 8(d) and 9(b). The time domain plots before and after balancing for rotor Model 1 at disc 1 is shown in Figure 8(a,c), respectively. As shown in the plots, the amplitude of vibration has been reduced to the greater extent with the addition of the balance masses at the balancing planes.

Figure 8. (a) Displacement versus time curve without balancing, (b) displacement amplitude versus rotor run-up speed with and without balancing, (c) displacement versus time curve with balancing and (d) displacement phase versus rotor run-up speed (for Model 1 at disc 1, without magnetic bearing: magnetic force is considered only in the numerical simulation for obtaining unbalances.

Figure 9. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 1 at conventional bearing 1 without magnetic bearing: the magnetic force is considered only in numerical simulation for obtaining unbalances).

Table 4. Comparison of the assumed and estimated unbalance magnitudes and phases from the numerical simulation with residual unbalances only at discs and without unbalance in shaft for rotor Model 1 (with noise).

	Assumed values		Estimated values
Disc	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)
Disc 1	0.2340×10^−4	0.1000	0.2340×10^−4	0.1000
Disc 2	0.7200×10^−4	1.2000	0.7200×10^−4	1.2000
Disc 3	0.3600×10^−4	2.1000	0.3600×10^−4	2.1000

Table 5. Comparison of estimated unbalance magnitude and phase from GICM with residual unbalances in the shaft as well at discs for all three models shown in Figures 3, 13 and 14.

	Model 1		Model 2		Model 3
Disc	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)
Disc 1	0.3578 × 10^−4	0.2210	0.2564 × 10^−4	0.3593	0.2138 × 10^−4	0.1900
Disc 2	0.4338 × 10^−4	1.9140	0.3837 × 10^−4	1.4734	0.8557 × 10^−4	1.1527
Disc 3	0.4375 × 10^−4	2.2796	0.2982 × 10^−4	1.6659	0.6418 × 10^−4	1.2293

4.1. Comparison of generalized and CICMs

In order to compare the suitability of the GICM and the CICM (while the AMB is not considered in rotor system) in the estimation of unbalances, by the latter method first displacements at two conventional bearing positions are taken and then additional displacements at AMB locations of the GICM (as such the AMB was absent in CICM) is also considered as shown in Table 6. The estimation of unbalances by using only displacements at conventional bearing positions are not appropriate as there is less reduction in unbalance responses after balancing (Figure 10). The variation in estimated unbalances is because of using displacements at less number of locations. This may be due to on considering measurement at intermediate shaft locations, the modal information is captured more accurately. As provided in Table 5 and Table 6, estimated unbalances by two methods (GICM and CICM) are nearly the same. The reduction in unbalance responses for the GICM has been shown in Figures 8, 9 and for the CICM in Figures 11 and 12 in the frequency domain for Model 1 at disc 1 and conventional bearing 1.

Figure 10. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 1 at disc 1 without magnetic bearing: magnetic force is not considered in numerical simulation for obtaining unbalance).

Figure 11. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 1 at disc 1: displacement at the node position of magnetic bearing is taken without considering magnetic bearing in system by CICM).

Figure 12. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 1 at conventional bearing 1: displacement at the node position of magnetic bearing is taken without considering magnetic bearing in system by the CICM method).

Table 6. Comparison of estimated unbalance magnitude and phase using CICM (without magnetic bearing).

With and without unbalances in shaft for Model 1 (without considering displacement at node position of the
magnetic bearing)
	Assumed unbalances in discs		Estimated unbalances without unbalance in shaft		Estimated unbalances with unbalance in shaft also
Disc	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)
Disc 1	0.2340 × 10^−4	0.1000	0.2340 × 10^−4	0.0999	0.3329 × 10^−4	0.1984
Disc 2	0.7200 × 10^−4	1.2000	0.7200 × 10^−4	1.2000	0.4511 × 10^−4	1.9352
Disc 3	0.3600 × 10^−4	2.1000	0.3600 × 10^−4	2.0999	0.4951 × 10^−4	2.3365
With and without unbalances in shaft for Model 1 (with considering displacement at node position of the magnetic bearing)
Disc 1	0.2340 × 10^−4	0.1000	0.2340 × 10^−4	0.0999	0.2488 × 10^−4	0.3320
Disc 2	0.7200 × 10^−4	1.2000	0.7200 × 10^−4	1.2000	0.3966 × 10^−4	1.6565
Disc 3	0.3600 × 10^−4	2.1000	0.3600 × 10^−4	2.0999	0.3222 × 10^−4	1.7812

However, it should be noted that as the balancing is required for the high-speed rotor system, so it would be difficult to run the rotor system supported on conventional bearings up to the operating speed for obtaining displacements. Conventional bearings may get damaged at high speeds while running without balancing for taking measurements. Moreover, measurements may not be reliable for the CICM, which also requires linearity assumption of the system. Herein, the use of AMB is advantageous in obtaining the reliable residual unbalances as it reduces the vibration to have linear system assumption valid. It also prevents the system from getting damaged at high-speed and ensures more reliable data.

Numerical simulation has been done initially for the rotor system with high acceleration (i.e. 139 rad/s²). In the same manner with low acceleration, i.e. at 34.75 rad/s², the rotor system has been tested numerically. In both the cases, the estimated unbalance masses are found to be similar as shown in Table 5 and Table 7.

Table 7. Comparison of estimated unbalance magnitude and phase from GICM with residual unbalances in the shaft as well at discs for Model 1 with low angular acceleration (i.e. 34.75 rad/s²).

	Assumed unbalances in discs		Estimated unbalances without unbalance in shaft		Estimated unbalances with unbalance in shaft also
Disc	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)	Magnitude (kg m)	Phase (rad)
Disc 1	0.2340 × 10^−4	0.1000	0.2340 × 10^−4	0.1000	0.3606 × 10^−4	0.2184
Disc 2	0.7200 × 10^−4	1.2000	0.7200 × 10^−4	1.2000	0.4352 × 10^−4	1.9225
Disc 3	0.3600 × 10^−4	2.1000	0.3600 × 10^−4	2.1000	0.4391 × 10^−4	2.2869

5. Comparison of unbalance estimations with different rotor models

In order to determine the suitable position of AMB, apart from the rotor system in Figure 3, two different models have been selected having different position of the AMB as shown in Figure 13 (Model 2) and Figure 14 (Model 3). The GICM is then used to estimate the unbalances in these two models and estimates of unbalances are compared in Table 5.

Figure 13. The rotor-bearing-AMB system (Model 2).

Figure 14. The rotor-bearing-AMB system (Model 3).

The amplitude of unbalance responses is lower for Model 2 than Model 1 and Model 3 as shown in Figure 15 and Figures 8 and 16, respectively. On the other hand, Model 3 shows the largest amplitude of unbalance response amongst the others. Variation of the displacement amplitude versus the rotor run-up speed have been plotted models at disc 1 and conventional bearing 1 positions for rotor systems to check the reduction in unbalances after adding correction masses for both rotors and without using AMB as shown in Figures 15–18. As shown in these plots, unbalance responses are found to be reduced at different critical speeds after balancing as shown in Figures 15(a)–18(a). The variation of the phase of the unbalance response is shown in Figures 15(b) and 18(b).

Figure 15. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus run-up speed (for Model 2 at disc 1 without magnetic bearing: magnetic force is considered in numerical simulation for estimating unbalances).

Figure 16. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 3 at disc 1 without magnetic bearing: magnetic force is considered in numerical simulation for estimating unbalances).

Figure 17. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 2 at conventional bearing 1 without magnetic bearing: magnetic force is considered in numerical simulation).

Figure 18. (a) Displacement amplitude versus rotor run-up speed with and without balancing and (b) displacement phase versus rotor run-up speed (for Model 3 at conventional bearing 1 without magnetic bearing: magnetic force is considered in numerical simulation).

The percentage deviation in residual unbalances estimation is also determined for 1% and 5% noise in responses as well as for ±1% and ±5% variation in Young’s modulus (E) and density (ρ). The percentage change in the magnitude and phase of residual unbalances is very low with the addition of noise and variation of physical properties of the rotor system as shown in Tables 8 and 9.

Table 8. Percentage deviation in magnitude and phase of unbalances with variation in the noise level.

	1% noise		5% noise
Disc	Magnitude	Phase	Magnitude	Phase
Disc 1	−0.0430	0.0680	−0.2125	0.3400
Disc 2	0	−0.0196	0.0288	−0.1047
Disc 3	0	−0.0100	−0.0253	−0.0652

Table 9. Percentage deviation in magnitude and phase of unbalances with variation in Young’s modulus and density.

	E + 1% of E and ρ + 1% of ρ		E + 5% of E and ρ + 5% of ρ
Disc	Magnitude	Phase	Magnitude	Phase
Disc 1	1.7650	−2.1421	9.3414	−9.9286
Disc 2	−0.3878	1.2759	−1.6880	6.7526
Disc 3	0.5811	0.3162	2.9308	1.7112
	E − 1% of E and ρ − 1% of ρ		E − 5% of E and ρ − 5% of ρ
Disc	Magnitude	Phase	Magnitude	Phase
Disc 1	−1.8080	2.3121	−8.1360	11.7307
Disc 2	0.4334	−1.2825	1.2318	−5.8496
Disc 3	−0.5811	−0.3262	−2.7287	−1.3600

6. Conclusions

In the present work, an identification algorithm for the estimation of unbalances at selected balancing planes has been developed with the help of AMB as a controller in flexible rotor-bearing systems. Numerical responses obtained show that with the addition of correction masses there is an appreciable reduction in unbalance responses corresponding to different critical speeds. This shows that the present methodology predicts quite accurate values of magnitude and phase of unbalances that is present in the rotor system. The variation in estimated values of unbalances is found to be very low even in the presence of measurement noise and modelling errors. The above results also show that the use of GICM in the estimation of unbalances provides the nearly similar result as the CICM without AMB. However, the GICM also has the advantage of having AMB at high speed during measurement and subsequently while estimation than the CICM in which it is difficult to cross critical speeds without balancing a rotor for taking measurements. So, it is advantageous to use the GICM over conventional method for the rotor system at high speeds. The use of displacements at more number of locations in the estimation of unbalances gives more appropriate results. Among rotor models considered, the Model 2 is found to be more appropriate than other models (since it has much less reduction in unbalance responses after balancing). As we can see in Model 2, the position of AMB is approximately at the middle of the shaft between the conventional bearings, so the placement of the AMB is very much important. Application of the proposed methodology in a test rig or practical rotor would be very interesting. As a future work, physical trial masses can be replaced by ‘virtual unbalance’ through AMB, which would give similar force as of a ‘physical unbalance’. However, it needs multiple AMBs at each balancing planes that would pose practical difficulty.

Disclosure statement

No potential conflict of interest was reported by the authors.

Appendix

Rotor system matrices and vectors

Detailed matrices of equations of motion Equation (1) are given below.

Displacement vector corresponding to each element (A1) ${\mathit{\eta }}_x=\left\{\begin{array}{c}{x}_1\\ {\phi }_{y_1}\\ x_2\\ {\phi }_{y_2}\end{array}\right\}\mathbf{\text{and}}{\mathit{\eta }}_y=\left\{\begin{array}{c}{y}_1\\ {\phi }_{x_1}\\ y_2\\ {\phi }_{x_2}\end{array}\right\},$(A1) where η_x and η_y are the displacement vectors for each element in x- and y-directions, respectively.

Magnetic force vector (A2) ${\mathbf{f}}_{mb}=-{\mathbf{K}}_{\mathit{s}}{\mathit{\eta }}_m\mathbf{+}{\mathbf{K}}_{\mathit{i}}{\mathbf{i}}_c$(A2) with ${\mathbf{K}}_s=\left[\begin{array}{cc}{k}_{s_{xx}}& k_{s_{xy}}\\ k_{s_{yx}}& k_{s_{yy}}\end{array}\right],{\mathbf{K}}_i=\left[\begin{array}{cc}{k}_{i_{xx}}& k_{i_{xy}}\\ k_{i_{yx}}& k_{i_{yy}}\end{array}\right],$ ${\mathit{\eta }}_m=\left\{\begin{array}{c}x\\ y\end{array}\right\},{\mathbf{i}}_c=\left[\begin{array}{cccccc}{k}_P& 0& k_D& 0& k_I& 0\\ 0& k_P& 0& k_D& 0& k_I\end{array}\right]\left\{\begin{array}{c}x\\ y\\ \dot{x}\\ \dot{y}\\ \overline{x}\\ \overline{y}\end{array}\right\},$ where K_s is the displacement stiffness matrix, K_i is the current stiffness matrix, η_m is displacement vector at AMB and i_c is the controlling current vector. Controller parameters k_P represents the proportionality constant, k_D represents the derivative constant and k_I represents the integral constant. The dot and bar on the displacement represents differentiation and integration, respectively.

Unbalance force vector is given as (A3) ${\mathbf{\text{f}}}_{unb}={\mathbf{\text{f}}}_d+{\mathbf{\text{f}}}_s,$(A3) where ${\mathbf{\text{f}}}_d$ is the assembled unbalance force vector for discs and ${\mathbf{\text{f}}}_s$ is the assembled unbalance force vector for the shaft. Force vector for discs in x and y-direction is given as (A4) $f_{d{x}_i}={\omega }^2{m}_i{e}_i\text{cos}(\omega t+{\phi }_i).f_{d{y}_i}={\omega }^2{m}_i{e}_i\text{sin}(\omega t+{\phi }_i)$(A4) where m_i is the mass of unbalance, e_i is the eccentricity, φ_i is the phase, and ω is the spin speed of the rotor.

In order to introduce shaft unbalance in the rotor system, each element of the shaft has been divided into a finite number of segments (n), each segment is considered as discs with some random eccentricity and phase. The unbalance force acting on the element of the shaft is shown in Figure 2. The equation of unbalance force for each segment is written as (A5) $f_i\left(t\right)=\rho A{l}_i{e}_i{\omega }^2\text{cos}(\omega t+{\phi }_i),$(A5) where ρ is the density of shaft, l_i is the length of each shaft segment, e_i is the eccentricity for the segment and A is the area of cross section of the shaft. The choice of unbalances in each shaft segment is arbitrary. The subscript i = (1, 2, … , n). Net unbalance force in each element is given by (A6) ${\mathbf{f}}_{se}={\int }_0^{l_1}\left\{f_1\left(t\right)\delta (z-z_1)+f_2\left(t\right)\delta (z-z_2)+\cdots +f_n\left(t\right)\delta (z-z_n)\right\}\mathit{n}\left(z\right)dz$(A6) with (A7) $\mathit{n}\left(z\right)=\left[\begin{array}{c}1-3\frac{z^2}{l^2}+2\frac{z^3}{l^3}\\ z-2\frac{z^2}{l}+\frac{z^3}{l^2}\\ 3\frac{z^2}{l^2}-2\frac{z^3}{l^3}\\ -\frac{z^2}{l}+\frac{z^3}{l^2}\end{array}\right],$(A7) where f_se is the net unbalance force in each element corresponding to four degrees-of-freedoms (DoFs) in each plane (x-z or y-z), δ is the Dirac delta function which enables the point force (shaft unbalance) as a distributed force description in form of integral, n(z) is the shape function for an 8-DOF shaft element subjected to transverse forces and z is the distance of the location of force from the right end of the element. Unbalance forces due to the shaft for each element are given as (A8) ${\mathit{f}}_{sxe}=Re\left[\left(\rho A\sum _{i=1}^n{e}_i{l}_i\mathit{n}\left(z_i\right){\omega }^2{\text{e}}^{i\left(\omega t+f_i\right)}\right)\right]$(A8) and (A9) ${\mathit{f}}_{sye}=Im\left[\left(\rho A\sum _{i=1}^n{e}_i{l}_i\mathit{n}\left(z_i\right){\omega }^2{\text{e}}^{i\left(\omega t+f_i\right)}\right)\right]$(A9) where vectors f_sxe and f_sye are the unbalance forces due to the shaft in x and y-direction, respectively.