Lumped parameters models for the stability analysis of rotors supported on gas bearings

In this paper different lumped parameters models are presented for the stability analysis of rotors supported on gas bearings. The analysis is carried out taking into account the damping and stiffness coefficients of journal bearings, calculated with the perturbation method. Lumped parameters models of different complexity are discussed for the description of the main rigid modes of the spindle. In case of symmetric systems, it is shown that the complexity of the model can be reduced by halving the number of the degrees of freedom. The analogy with the single mass approach is demonstrated. The considerations discussed in the paper are preliminary for the stability analysis of more complex systems, such as the ones with non-fixed bushes.


Introduction
Gas bearings systems suffer from the so-called half speed whirl, which is a self-excited phenomenon which occurs at high rotational speeds [1].It involves the increase of the radius of the rotor orbits until a possible crash of the rotor with the bearings.The onset of the half-speed whirl occurs when the rotational speed exceeds the stability threshold.This phenomenon can be studied with two alternative approaches: the socalled orbit method and the linearized coefficients method.The first method solves the time dependent Reynolds equation (RE) coupled with the equations of motion of the rotor [2].It takes into account the non-linear effects in bearings, but it is time expensive and it is unsuitable for creating stability maps.Explicit or implicit schemes can be employed to solve the time marching problem.A noteworthy method is the Alternating Direction Implicit (ADI) scheme [3], which involves the solution of a tridiagonal system and is less costly than the implicit scheme which involves a penta-diagonal system.
The second method is based on the linearization of the bearings forces around a given operational point [4] and on the solution of the eigenproblem resulting from the rotor equations of motion [5], [6].The linearized stiffness and damping coefficients are computed with the so-called perturbation method [7].The improved version of the method, presented by the same author in [8], is useful for rotor-dynamics calculation of unbalance response and stability.The importance of this method was confirmed by other authors, which implemented it to carry out the stability analysis of gas bearings systems.In [9] it was employed to analyze the stability of a foil bearing and obtain stability maps, which typically show the stability threshold curve on a plot where a mass parameter is represented vs the bearing number.The classical perturbation method and an extended perturbation procedure incorporating the foil compliance explicitly were developed in [10].The linearization method was used in [11] to compute the gas film coefficients for a high speed cryogenic turboexpander and in [12] for a heat pump turbocompressor.An interesting analysis of tilting pads is given in [13].Many other works on the topic can be found in literature, such as [14]- [17].The linearized coefficients method is more convenient from the computational point of view than the time advancing methods as the time needed to compute the linearized coefficients and evaluate the stability is much smaller than the time needed to solve the time dependent problem.Moreover, it is capable of providing a better knowledge about the underlying sources of instability, if compared with the orbit method.Stability charts can be easily plotted to show quantitatively how much the system approaches the stability threshold [6], [7], [18]- [20].
The stability evaluation is often carried out considering the so-called single-mass approach [21]- [23], which simplifies the complexity of the full rotor-bearings model and is often employed to have a first approximation of the stability threshold.Anyway, it is not evident in which cases this approach is more convenient and in which cases the full rotor model is to be preferred.To clarify this aspect, the linearized coefficients method to evaluate the stability of rotors supported on gas bearings is discussed in this paper.The full rotor model with four degrees-of-freedom (4 DOFs) system is presented, showing in which cases it is possible to reduce the number of DOFs and consider separately the cylindrical and the conical whirl motions.The 2 DOFs models are also compared with the single point-mass approach and the analogy between the two is demonstrated.Such considerations are preliminary for the stability analysis of more complex systems, such as the ones with non-fixed bushes.

The lumped parameters models
The dynamic equations of motion of the spindle considered rigid are considered in case of fixed bushings.The bearings forces are linearized around the central shaft positions and the coefficients of stiffness and damping are collected in matrixes [K] and [C].The system dynamics is described by eq. ( 1) where [M] and [G] are the stiffness and the gyroscopic matrixes respectively, q is the state vector of the generalized displacements and Funb is the vector of the actions due to the residual unbalance of the shaft.
Matrixes [K] and [C] depend on the linearized coefficients of the radial air films in journal bearings, which are frequency dependent.

Rigid rotor with 4 DOFs
The rigid rotor model with 4 DOFs involves all the radial degrees of freedom of the shaft.The axial degree of freedom is not considered as the axial mode and the radial modes can be considered decoupled in a first approximation; the thrust bearing is supposed to have a little influence on the whirl stability of the electrospindle due to its small radius compared with the axial distance of the journal bearings.It is possible to consider two alternative versions of the generalized displacements vector: vector q involves both the radial displacement of the center of mass of the shaft and the rotations of the shaft axis around axes x and y: Vector q' involves the radial displacements of the shaft measured in the middle radial plane of the journal bearings: The following kinematic relationship is valid where l1 and l2 are the axial distances between the shaft center of mass and the middle plane of journal bearings, see Figure 1.Notice that in general the system is not symmetric and distances l1 and l2 do not coincide.
Figure 1: sketch of the system's geometry In case q is used, the system is described by eq. ( 1), while in case q' is used, by eq (3) The expressions of the matrixes and of the unbalance vector are reported in appendix A1.It is worth noticing that in Equation ( 1) the mass matrix is diagonal, while stiffness and damping matrixes are full.This means that in general the four equations are coupled.Similarly for system (2), in which also the mass matrix is not diagonal.

Rigid rotor with 2 DOFs, cylindrical mode
This model involves just the radial translation of the shaft, neglecting the tilting motions around axes x and y.The state vector is The equations of motion in matrix form are where e is the radial distance between the shaft center of mass and the center of journals (residual static unbalance), m is the mass of the shaft and kij and cij the linearized coefficients of the air films.

Rigid rotor with 2 DOFs, conical mode
In this case only the rotation of the shaft axis around x and y fixed axes are taken into account, neglecting the radial motion of the center of mass of the shaft.The state vector is The equations of motion in matrix form are where I and Ip are the transversal and the polar moments of inertia of the shaft and  is the angle between the principal axis of inertia and the rotational axis of the shaft (residual dynamic unbalance).

Computation of gas film coefficients
The stiffness and damping coefficients of the air films are computed with the so-called perturbation method.Starting from the static pressure distribution p0, obtained solving the discretized Reynolds equation with finite difference technique, the dynamic pressure distribution p is obtained as a result of the sinusoidal perturbation x imposed to the shaft.As sketched in Figure 2, the coefficients are both function of the rotational speed  and the perturbation frequency .The computation of such coefficients involves the solution of a linear system, resulting after discretization of the Reynolds equation in perturbed form ( 6): Eq. ( 6) represents a balance of flow in the perturbed condition; also the input flow through the supply orifices is considered in the balance, taking into account its partial derivatives with respect to pressure p (  ′ ) and x (  ′ ).

The eigenvalue problem and the iterative procedure
The dynamic problem is written in the state space form, resulting in a linear system of size 2n, where n is the number of DOFs: The stability problem is faced calculating the eigenvalues of the state space matrix A. Complex eigenvalues appear in pairs, with conjugate imaginary parts: The real part represents damping, while the imaginary part the damped frequency.The natural frequencies and the damping factors of the modes are computed with: As the coefficients of the air films are frequency dependent, an iterative procedure is necessary to find the convergence between the damped frequency   and the perturbation frequency .The following algorithm is thus implemented for each mode: An example of the output of such algorithm is given in figure 3; the damped frequency and the damping factors of the four modes are plotted vs the rotational speed considering the 4 DOFs model with fixed bushings.The negative damping factor of mode 4 indicates at which rotational speed the whirl becomes unstable (about th=130 krpm).It is sufficient that one mode becomes unstable to have the unstable whirl.

Discussion
In  The stability threshold obtained from the models are compared in table 2.       The cylindrical and the conical motions are uncoupled in case the Coriolis terms are neglected and the system's geometry satisfies condition (10) or condition (11).In the first case there is perfect symmetry: distances l1 and l2 must coincide and the journal bearings must have the same geometry and the same radial film thickness (same air film coefficients): Condition ( 11) is more general, as the center of mass of the shaft could be not in the middle of bearings.
This is evident from the stiffness and damping matrixes of the 4dofs model (see appendix A1) which reduce to band diagonal matrixes.For instance, in condition (10) [K] is reduced into this form: It is thus demonstrated in which cases the full rotor model can be reduced: in case of complete symmetry (condition 10) or more in general according to condition (11).Of course, in case the system does not exactly satisfy such conditions but the asymmetry is small or condition (11) is almost satisfied, the two modes can be still studied separately, being the out-of-diagonal elements negligible.In the other cases, the full model is to be preferred as it better estimates the stability threshold.

The single mass approach
In literature the so-called single mass approach is known, which allows a quick estimation of the stability of rotor-bearings systems.As explained in [23], [24], the approach is quite common to perform the stability analysis of a rotor-journal bearings systems.According to this approach, the shaft mass is distributed in correspondence of the journal bearings considering a cylindrical motion or a conical motion.The masses ms1 and ms2 change in case the cylindrical or the conical whirl are considered, see figure 9.Such expressions are compared with the masses that appear in [M'] in case a pure cylindrical or a pure conical motion are considered, see appendix A2 for the detailed passages.It is thus demonstrated the analogy of the single mass approach with the simplified 2 DOFs models, in case of proper symmetry of the system.

Conclusions
This paper illustrates in details the linearized stability analysis of rotor-gas bearings systems.It shows in which cases a 4dofs rotor-bearings model can be simplified in a 2dofs model considering alternatively the cylindrical whirl or the conical whirl.In these cases the stability analysis can be carried out on each separate bearing concentrating on it a fraction of the shaft mass.The analogy of the reduced models with the socalled single mass approach is demonstrated.The models can be further extended to the stability analysis of gas bearings systems with non-fixed bushings.

Appendix A1: 4 DOFs model
The matrixes of the 4 DOFs model of the rigid shaft are here explicitated.In case the generalized displacement vector q is used, the following expressions result: (A1.1)In case the generalized displacement vector q' is used, the following expressions result: (A1.2) Appendix A2: analogy with the single-mass approach Cylindrical whirl Starting from system (3) and imposing the cylindrical whirl it is  1 =  2 =   ;  1 =  2 =   ;  = 0; The system becomes where for brevity the damping matrix is not indicated as it has a structure similar to the stiffness matrix.In case of symmetry ( 1 ≈  2 and  1 ≈  2 , or more generally  1  1 ≈  2  2 ) the last two lines are zeroed and the system becomes where ms1 and ms2 are the masses in correspondence of bearings, see table 8, and e1 and e2 their unbalance.It must be This is the equilibrium of forces considering the two masses together.It is equivalent to the force equilibrium of the masses taken singularly, as no radial force is transmitted by the shaft; in this case, the inertia force of the single mass is supported by the correspondent bearing.The advantage of considering the mass concentrated in correspondence of bearings is that the stability of the single bearing can be evaluated separately, with a reduced number of DOFs.

Conical whirl
Starting from system (3) and imposing the conical whirl it is some elements in the mass matrix are zeroed:  In case of symmetry ( 1 ≈  2 and  1 ≈  2 , or more generally  1  1 ≈  2  2 ) the first two lines are zeroed and the system can be further simplified neglecting the Coriolis terms: where it is ( −   ) = ( 1  1  1 −  2  2  2 ) (A2.8) Equations A2.7 represent to the equilibrium of moments of the masses concentrated in correspondence of bearings.As before, the equilibrium can be written for the single mass taken separately as no radial force is transmitted by the shaft and the inertia force of each mass is balanced by the corresponding bearing.Dividing by the arm li, the moment equilibrium simplifies into a force equilibrium.

Compliance with Ethical Standards
The authors certify that they have NO affiliations with or involvement in any organization or entity with any financial interest, or non-financial interest in the subject matter or materials discussed in this manuscript.

Figure 2 :
Figure 2: layout of the input and output variables involved in the static and in the perturbed problems

•
Step1: kij and cij are evaluated at a first attempt frequency  • Step 2: the eigenvalues of system (3) are calculated, obtaining the damped frequency d • Step 3: kij and cij are evaluated at the new frequency =d • Step 4: in case the new frequency coincides with the damped frequency the algorithm has reached convergence, otherwise come back to step 2.

Figure 3 :
Figure 3: damped frequencies d and damping factors  vs rotational speed; 4 dofs model, case a)

Figure 4 :
Figure 4: damped frequencies d and damping factors  vs rotational speed for case a); 2 dofs model, cylindrical whirl in the left; 2 dofs model, conical whirl in the right

Figure 5 :
Figure 5: damped frequencies d and damping factors  vs rotational speed; 4 dofs model, case b)

Figure 6 : 4
Figure 6: damped frequencies d and damping factors  vs rotational speed for case b); 2 dofs model, cylindrical whirl in the left; 2 DOFs model, conical whirl in the rightThe stability threshold obtained from the models are compared in table4.Table 4: comparison of the stability threshold according to the considered models; case b) Model Stability threshold (krpm) 2 DOFs cylindrical whirl model 16430 rad/s≈157 krpm 2 DOFs conical whirl model 16800 rad/s≈160 krpm 4 DOFs model 15700 rad/s≈150 krpm

Figure 7 :
Figure 7: damped frequencies d and damping factors  vs rotational speed; 4 dofs model, case c)

Figure 8 :
Figure 8: damped frequencies d and damping factors  vs rotational speed for case c); 2 dofs model, cylindrical whirl in the left; 2 DOFs model, conical whirl in the right In this case, due to the complete symmetry of the journal bearings with respect to the center of mass of the shaft, the results of the 2 DOFs models perfectly match with the results of the 4 DOFs model.The curves representing the damped frequencies and the damping factors coincide.Table7: comparison of matrixes [K] and [C] in cases a), b) and c); the coefficients are calculated with the perturbation method supposing =200 krpm and =9000 rad/s.

Table 2 :
comparison of the stability threshold according to the considered models; case a) DOFs models with the 4 DOFs model, the damped frequencies and the damping ratio of the modes do not coincide.The reason is that in this case the cylindrical and the conical whirl motions are coupled and the 2 DOFs models, which take into account only the pure cylindrical or the conical whirl, result to be inaccurate.The complete 4 DOFs model is to be preferred as it treats the coupled cylindrical and conical whirl.

Table 3 :
geometry considered in case b) The 4 DOFs model results are here below depicted.

Table 6 :
comparison of the stability threshold according to the considered models, case c)