Dynamic analysis of the variable stiffness support rotor system with elastic rings

Elastic rings are common rotor supporting structure, which have been widely used in aeroengine rotor support system. However, large inertia force and gyroscopic moment may occur during the operation of aeroengine, which may lead to contact between elastic ring and bearing pedestal, and then introduce variable stiffness into the rotor support system. In this paper, a mathematical model of variable stiffness of elastic ring is proposed and this model is subsequently verified by comparison with simulation analysis and experimental results. Based on this model, a variable stiffness model of an elastic ring-supported rotor is developed by coupling the kinetic equations of the rotor with the deformation of the combined support. Then, the spectrum cascades are used to analyze the dynamic characteristics of the rotor system. In addition, the influences of the variable stiffness of elastic ring on the critical speed of the system are also examined. Finally, some simulation results are verified by experiments on a combined test bench of an elastic ring-supported rotor.

Then, the spectrum cascades are used to analyze the dynamic characteristics of the rotor system. In addition, the influences of the variable stiffness of elastic ring on the critical speed of the system are also examined. Finally, some simulation results are verified by experiments on a combined test bench of an elastic ring-supported rotor. Translating inertial matrix of disk element, rotating inertial matrix of disk element, gyroscopic moment matrix of disk element M r , G r Mass matrix of rotor substructure, gyroscopic matrix of rotor substructure K r , C Stiffness matrix of rotor substructure, damping matrix of rotor substructure F e , F b , F g Unbalance vector, bearing force vector, gravity vector f n1 , f n2 First and second natural frequency of the rotor system n 1 , n 2 First and second modal damping ratios k B , k sq Contact stiffness of bearing, stiffness of squirrel cage r 0 , N b , h j Bearing clearance, number of ball elements, angle location x r , x o , y r , y o Translation of shaft and outer rings of the bearing along the x-and y-axes w c Angular velocity of the cage Radius of outer race, radius of inner race n sq , b sq , h sq , L sq Number of cage strips, section width of cage strips, section height of cage strips, length of cage strips f erx , f ery Elastic ring force in the x and y directions f sqx , f sqy Squirrel cage force in the x and y directions F sq , F er , F b Squirrel cage force vector, elastic ring force vector, bearing force vector Supports with elastic rings (ER) have been widely used in aeroengine rotor system [1], which are applied to adjust the rotor critical speed from the operating mode in order to reduce the vibration. To limit the deformation of elastic supports, limiting amplitude structures are set up [2]. However, the aeroengine rotor system with ER is usually exposed to great evolutionary overloading [3]. If the limited structure parameters are selected too small or the circumferential distribution is uneven, in extreme cases, due to the action of large inertia forces, the elastic supports will contact with the limiting structures, and the supporting stiffness will change, so that the system will become a rotor support system with variable stiffness [4]. Hence, in order to study the variable stiffness characteristics effectively, it is necessary to develop the combined support rotor system model with ER and analyze the dynamic characteristics of the variable stiffness rotor system.
Among the existing ER solution methods, IS [3], Artemov, Hronin and finite element model (FE) are the most common methods to solve the stiffness [5]. However, the above methods only consider the deformation of ER within the allowable range. Thus, these methods are not suitable for the rotor with great unbalanced forces and huge weight. To limit the deformation of elastic support, Ma et al. [2] presented a novel structure with the amplitude limiting bulges. However, in extreme cases, the ER will contact the base, resulting in nonlinear stiffness of the ER. Luo et al. [6] proposed the finite element model of ER and the static stiffness test device of ER was designed based on the Guo's test [7]. Then, the characteristic of piecewise linear stiffness of ER was revealed by comparison between experiment and FE simulation. As ER exists contact stiffness, it is necessary to study contact mode. Pereira et al. [8] and Machado et al. [9] established a non-coordinated contact model of the shaft-hole with clearance. Persson et al. [10] and Ciavarella et al. [11,12] ignored the influence of hole's thickness on the contact characteristics of the shaft-hole and established the contact force model of the shaft-hole. Li et al. [13] considered the hole's thickness and established the shaft-hole coordinated contact model. This model could provide reference for the study of contact and collision between the shaft and the hole.
From the above literature, emphasis is put on the stiffness of ER and the dynamic analysis of the rotor system is not considered. Li et al. [14,15] adopted numerical and experimental methods to analyze the dynamic properties and stability of a rotor system with bolted joint. Sun et al. [16] investigated the responses of a dual-rotor system with rub-impact and analyzed the stability of system. Considering the influences of clearance and unbalance on rotor system, Villa et al. [17] established the rolling bearing model with the internal clearance and Hertz nonlinearity and presented the dynamic analysis of a rotor system with ball bearings. Hou et al. [18] analyzed the primary resonance of a dual-rotor system with ball bearings. In order to analyze the complex rotor system with combined support, Luo et al. [19] established a combined support rotor model and analyzed the dynamic response at different speeds. Chen et al. [20] proposed a novel coupling model and investigated the vibration performances of the aeroengine. In recent years, the elastic ring squeeze film damper (ERSFD), which is installed in many modern aeroengine, has been grown by many scholars. Zhang et al. [21] established a coupled oil film Reynolds equation and dynamic equation of a rotor system with ERSFD and investigated the dynamic parametric characteristics of the ERSFD. Considering the contact between the journal and the ER, Wang et al. [22,23] proposed an analytical method to investigate the dynamic properties of the ERSFD. Chen et al. [24] established the spiral bevel gear drive model with ERSFD and studied the dynamic characteristics of this system. Han et al. [25,26] proposed a semi-analytic method to analyze the oil film force of ERSFD and established a Jeffcott rotor model with ERSFD. Subsequently, this model was solved by Runge-Kutta method and the dynamic response of the ERSFD-rotor was analyzed. But there are relatively few studies on the dynamic analysis of the ER-rotor with variable stiffness [6].
Aiming at the existing issues, a variable stiffness model of ER is presented and the Newmark-b method is applied to reveal the dynamic characteristics of the variable stiffness ER-rotor system. Meanwhile, the correctness of the model is evaluated by the variable stiffness experiment of ER and some dynamic properties of ER-rotor are substantiated by numerical and experimental researches. Moreover, the research results of ER-rotor system can be a stepping stone in other rotor systems containing elastic rings.
The structure of this paper is as follows. The theoretical framework and basic idea of the variable stiffness model of ER are introduced and the kinetic equations of the ER-rotor with variable stiffness are established in Sect. 2. Then in Sect. 3, the nonlinear characteristics of rotor system caused by the ER are analyzed in detail. In Sect. 4, two experimental studies are illustrated to verify the ER model and some of the simulation results. Finally, the conclusions are summarized in Sect. 5.

Mathematical model and equations of motion
To conduct the dynamic investigation of a rotor system with variable stiffness elastic rings, the kinetic equations of the ER-rotor system shown in Fig. 3 are established in this section. Moreover, the variable stiffness model of elastic ring is divided into two parts according to the displacement of the elastic ring: the linear stiffness model (Fig. 1) and the contact stiffness model (Fig. 2). Finally, according to the literature [19,20,[27][28][29], the differential equations of the combined support rotor system with ER are developed.

The variable stiffness model of elastic ring
The displacement of the inner bulge of the elastic ring is related to the displacement of the rotor system. When the displacement of the inner bulge is less than the outer bulge height as shown in Eq. (1), the elastic ring is linear stiffness and the force analysis is shown in Fig. 1. One inner bulge is selected for analysis, and the other inner bulges are similar. F 1 represents the radial force of the inner bulge, and its components in the x-and y-axes are F x and F y , respectively. o 1 and o 2 are the centroid of ER (bearing pedestal) and journal, respectively.
where e represents the radial displacement of the journal; x and y indicate the journal's displacement in the x and y directions; D is the height of the outer bulge.
The stiffness of the ER bulges is much greater than the rest of the ER, so it is assumed that the bulges of the ER are rigid body [24]. Meanwhile, the inner and outer bulges of the elastic ring are in directly contact with the inner and outer bushing. Therefore, the displacement of the inner bulges under the operating mode can be expressed as follows: where U i represents the displacement of the inner bulges; i is the ith inner bulge; and n in is the number of the inner bulges; a indicates the angle between inner and outer bulge. The flexibility of the ring section is assumed as the flexibility at the central section of the fixed beam at both ends.
where k represents the flexibility at the central section of the ring section; L is the length of the ring section; E is the elastic modulus of ring material; I indicates the moment of inertia of ring section; b is the width of the ring; h is the thickness of the ring. The sum of the components of the force of each ring segment in the x and y directions are F x and F y , and the magnitude of the components can be expressed as: Finally, the stiffness of the linear segment of the elastic ring can be written as: When the displacement of the inner bulge is greater than the outer bulge height, the elastic ring is in contact with the bearing pedestal and its stiffness becomes contact stiffness. Before establishing the contact model between the elastic ring and the bearing pedestal, the following assumptions are put forward: (1) Since the height of the bulge is much smaller than the elastic ring, the bulge height is ignored. (2) The journal, inner bushing and elastic ring are approximated as shaft model. (3) The outer bushing and bearing pedestal are approximated as hole model. Then the shaft-hole contact model is introduced into the elastic ring contact model and the force analysis of the ER contact model is shown in Fig. 2. F represents the radial force, and its components in the x and y axes are F x and F y , respectively. o 1 , o 2 and o 3 are the centroid of the bearing pedestal, the journal without contact deformation and the journal under contact deformation, respectively.
According to the literature [13], the nonlinear contact stiffness of the ER can be expressed as follows: where K 2 represents the contact stiffness of the ER; R 1 is the ER radius and R 2 is the bearing pedestal radius; E is the elastic modulus of ring material; d indicates the depth of the contact; DR indicates the radius difference between bearing pedestal and ER; h min is the minimum wall thickness of bearing pedestal. Then, when the displacement is greater than the height of the outer bulge, the overall force of the elastic ring can be expressed as: The overall force's components in the x and y directions are F x and F y , and the magnitude of the components can be written as: Finally, the variable stiffness of elastic ring can be expressed as:

Motion equations of the ER-rotor system
In order to analyze the rotor system with variable stiffness combined support, the rotor system is divided into rotor substructure and support substructure as shown in Fig. 3. Then, the finite element method and the lumped mass method [29] are used to establish rotor substructure model and support substructure model, respectively, where k b , k sq and k er denote the stiffness of bearing, squirrel cage and ER, respectively. m i is the mass of bearing inner race and m o represents the mass of bearing outer race. The gyroscopic moment and rotational inertia of the shaft is considered to investigate the bending vibration of the shaft depicted in Fig. 3. Meanwhile, the rotor substructure is divided into several finite element units and each unit node has four DOFs (x, y, h x , h y ). The equation for the elastic shaft element is given as: where X,u s and Q s represent the rotating speed, the displacement vector and the excitation vector of the shaft; M s T ,M s R ,G s ,K s B are the translating inertial matrix, the rotating inertial matrix, the gyroscopic moment matrix and the stiffness matrix of the beam element, which were described in references [30][31][32].
Because the disk is assumed as a rigid disk, the strain energy is ignored and only consider the kinetic energy. Then, the differential equation of the disk which is obtained by Lagrange's equations is given as: where u d and Q d represent the displacement vector and the excitation vector of the disk; M d T ,M d R ,G d are the translating inertial matrix, the rotating inertial matrix and the gyroscopic moment matrix of the disk element, which were given in reference [19].
Finally, the dynamic equation of the rotor substructure is obtained by assembling the motion equation of the disk, the motion equation of the shaft and the supporting force.
where the displacement vector u r can be written as: and M r ,G r ,K r and C are mass matrix, gyroscopic matrix, stiffness matrix and damping matrix of the rotor substructure; F e is the unbalanced force vector; F b is the bearing force; F g is the gravity vector and n is the number of the nodes of the rotor substructure.
The matrix of Rayleigh damping (C) can be expressed as: where where f n1 , f n2 represent the first and second natural frequency of the rotor system and n 1 , n 2 are the first and second modal damping ratios. The combined support substructure consists of the squirrel cage, ER and the rolling bearing. Meanwhile, the supporting substructure and the rotor substructure are connected by the displacement and the interaction force of the bearing inner race and bearing outer race. Then, based on the theory of Hertz contact [33,34], the bearing force can be expressed as: where k B is the contact stiffness of bearing; r 0 is the clearance of bearing; N b is the number of ball elements; H is the Heaviside function; x r , x o , y r , y o are the translation of shaft and outer rings of the bearing along the x-and y-axes; h j is the angle of the jth ball, which can be expressed as [35,36]: where w c is the angular velocity of the cage; w c = w 9 r b /(R b ? r b ); R b and r b are the radius of bearing inner and outer race, respectively. The nonlinear force of the ER is computed based on the variable stiffness model. It can be expressed as: The squirrel cage force can be expressed as [37][38][39]: where n sq is the number of cage strips; E is the elastic modulus of material; b sq and h sq are the width and height of cage strips, respectively; L sq is the length of cage strips. According to the coupling relationship among the squirrel cage, ER, the bearing and rotor substructure as shown in Fig. 3, the kinetic equation of the combined support rotor system can be written as: where where F sq , F er and F b are the squirrel cage force vector, the ER force vector and the bearing force vector of the support substructure and F sq , F er and F b could be acquired according from Eqs. (19) to (23). The solution flowchart of elastic ring rotor system is shown in Fig. 4. By solving the elastic ring stiffness model and the combined support rotor system model, the dynamic characteristics of the variable stiffness combined support rotor system are analyzed.

Numerical study and discussion
In this section, the linear stiffness model and the variable stiffness model of the ER are introduced into the rotor system, and the Newmark-b method is used to solve the linear and nonlinear ER-rotor system, respectively. Then, by analyzing the characteristics of linear and nonlinear ER-rotor systems, the influence of variable stiffness on dynamic characteristics of ERrotor system is illustrated. In Sect. 3.1, the spectrum cascades under different motor speeds of the two system are presented to illustrate the influence of variable stiffness on frequency domain characteristics. The influence of unbalance force on critical rotational speed of the variable stiffness rotor system is illustrated by the amplitude frequency responses of the two systems under different unbalance forces (Sect. 3.2). Finally, the influence of clearance on critical rotational speed of ER-rotor system is analyzed in Sect. 3.3.

Spectrum cascades analysis
As illustrated in Fig. 3, the combined support ER-rotor system is adopted in this section, which is modeled by using the FE method [40][41][42][43]. The ER-rotor system is divided into 6 sections, 13 beam elements of Timoshenko and 2 support elements. The Timoshenko beam element has 2 translational and 2 rotational degrees of freedom, as illustrated in Fig. 3. The dimensional parameter values of the ER-rotor system as illustrated in Fig. 3 are given in Table 1. The disk is located at node 7 and the bearings are considered as nonlinear, which are, respectively, located at node 4 and node 12. The configurations of the ball bearings are shown in Table 2. The squirrel cages are simplified as spring elements which are computed from Eqs. (23) and the configurations of the squirrel cages are shown in Table 3. The model used for ER is shown in Sect. 2.1 and the dimensional parameter values of the ER are given in Table 7.
The ER-rotor system models with linear stiffness (the ER is simplified as a linear spring model.) and variable stiffness of ER (the ER is considered as variable stiffness model.) are established to compare the frequency domain characteristics of these two models through the spectrum cascades. In addition, the dimensional parameter values of these two models are the same. The spectrum cascades for the ER-rotor system with rotational speed, varying from 500 rev/ min to 2800 rev/min, are acquired as shown in Fig. 5.
The results from the spectrum cascades show that only the fundamental frequencies are prominent in the spectrum cascades of the model with linear stiffness (see Fig. 5a). As a contrast, for the model with variable stiffness, the spectrum cascades not only have the fundamental frequencies, but also emerge other  frequency components (1.1f n , 2f n and 3f n ). Meanwhile, f ba is also prominent in the spectral cascades, which is the natural frequency of the ER-rotor system. Thus, the complex frequency components will have the significant influence on the ER-rotor system.

Influence of the unbalanced mass
In order to study the influence of the unbalanced mass on the critical speed of the ER-rotor system, in this section, the configurations of the ER-rotor system are the same except for the parameters of the unbalanced mass. The different parameter combinations of the unbalanced mass are given in Table 4. The amplitude frequency responses under four different unbalanced masses are acquired, which is shown in Fig. 6. Figure 6a depicts the amplitude frequency responses (the ER model with linear stiffness) corresponding to different unbalanced mass. It can be found that the first-order critical speed does not change and the amplitude increases with the increase of unbalance mass. As illustrated in Fig. 6b, the amplitude frequency responses (the ER model with variable stiffness) are different from Fig. 6a. The critical speed remains unchanged in the relatively low unbalanced mass (U m 3 and U m 4), but the relatively high unbalanced mass causes the first-order critical speed shift right (U m 1 and U m 2), which is because of the variable stiffness of the ER. Figure 7 shows the stiffness of the ER under different clearances. Meanwhile, the different parameter values of the clearance are shown in Table 5.

Influence of the clearance
Due to the existence of clearance, the stiffness of ER presents the characteristics of piecewise linear stiffness (the ER model with linear stiffness) and variable stiffness (the ER model with variable stiffness). From Fig. 7a, the stiffness of the ER is zero  Length of cage strips L sq (mm) 40 Height of cage strips h sq (mm) 3 Stiffness of squirrel cage (N/m) 2.53 9 10 6 (a) (b) Fig. 5 Spectrum cascades of the ER-rotor system a the ER model with linear stiffness and b the ER model with variable stiffness when the displacement is less than the clearance, while the stiffness is linear when the displacement is greater than the clearance (piecewise linear stiffness).
Different from Fig. 7a, b shows that when the displacement of the ER is less than the sum of the height of the bulge and the clearance, the stiffness  trend of the ER is the same as that in Fig. 7a. When the displacement is greater than the sum of the height of the bulge and the clearance, the stiffness shows a variable stiffness characteristic (variable stiffness). The amplitude frequency responses for the present ER-rotor system with rotational speed (500 rpm * 2800 rpm), are obtained as shown in Fig. 8. The configurations of the ER-rotor system are shown in Tables 1, 2, 3, 4 and 5. The author prefers to assigned the two parameters (M = 0.0118 kg, e = 40 mm) directly to study the effect of different clearance on the critical speed. It can be concluded that with the increase of clearance, the critical speed decreases, while the vibration amplitude increases accordingly (both Fig. 8a, b have the same trend).
The last simulation study concerned on the influence of the unbalanced mass of the critical speed under certain clearance. Define clc = 0.005 mm and e = 40 mm as basic control parameters, and the parameters of the rotor system remain unchanged. The control parameters of the unbalanced mass are given in Table. 6 and the amplitude frequency responses of the unbalanced mass under certain clearance are shown in Fig. 9. In Fig. 9, the gray vertical line represents the critical speed of the rotor system under low unbalanced mass and piecewise stiffness; the red vertical line represents the critical speed of the rotor system under medium unbalanced mass and piecewise stiffness; the black vertical line represents the critical speed of the rotor system under large unbalanced mass and variable stiffness. Obviously, the amplitude frequency response curves exhibit in Fig. 9 show that the clearance which is caused by design and machining will make the critical speed change. From Fig. 9a, the relatively low unbalanced mass makes the critical speed shift right, while the critical speed does not change when the unbalance mass reaches a certain degree. Different from Fig. 9a, b shows that the clearance results in an increase of the critical speed at the relatively low unbalanced mass and the critical speed remains unchanged at the relatively high unbalanced mass. When the  unbalanced mass is large, the critical speed continues to move to the right and this is due to the action of clearance and variable stiffness.

Experimental verification
In this section, two experiments are carried out to verify the variable stiffness model of elastic ring and the partial dynamic characteristics of rotor system with elastic rings. The purpose is to demonstrate the existence and applicability of the variable stiffness model and provide guidance for the application of elastic rings in aeroengine rotor system.

Test verification of variable stiffness of elastic ring
In order to definitively prove the simulation results of the ER model with variable stiffness, the experiment about the displacement of the ER under different load has been completed on an ER stiffness test bench as shown in Fig. 10. The ER stiffness test bench consists of electrohydraulic test machine, ER, inner and outer bushing, rigid shaft, loading frame and clamp. Since the stiffness of the clamp is much larger than that of the ER, the clamp is regarded as the rigid body and its deformation is ignored. Meanwhile, the electrohydraulic test machine can provide the load in the vertical direction and record the displacement of the ER under the corresponding load. Then, the stiffness of ER is obtained by drawing load and displacement curves.
In order to verify the correctness of the variable stiffness model, Table 7 summarizes the parameters of the test device used in this study. Then, the stiffness of elastic ring at different circumferential positions is experimentally studied and three groups of special positions of inner bulges, holes and outer bulges are selected, respectively, for vertical loading, as shown in Fig. 11a, b and c. Figure 11 depicts the stiffness test and simulation comparison diagram of elastic ring at different circumferential positions. It is worth noting that the variable stiffness test of elastic ring has two aspects to emphasize: (1) the correctness of the variable stiffness model is verified by comparing the simulation results with the experimental results. Meanwhile when the displacement is greater than the thickness of the bulges, the elastic ring appears variable stiffness. (2) Through the comparison of test and simulation results at different circumferential positions, it is proved that the variable stiffness model is not only suitable for a single position, but also suitable for all circumferential positions.

Experimental study of ER-rotor system
To verify the correctness of the simulation results for the ER-rotor system with variable stiffness, the experiment about the critical speed of the ER-rotor system under different unbalance mass has been completed on the ER-rotor test bench as shown in Fig. 12a. The ER-rotor test bench consists of a motor, disk, LMS, eddy current sensor, combined support, etc. The combined support consists of bearing, squirrel cage and ER as shown in Fig. 12b. Since the sliding of the ER will affect the analysis results, we designed pin holes at the combined support in Fig. 12b to fix the ER with pins to prevent the sliding of the ER. Meanwhile, the combined support is fixed on the bearing pedestal by bolts. The motor powers the rotor system and the rotating speed varies from 0 to 3000 rpm. The combined supports at both ends of the shaft have the same configuration and the parameter values of the ER-rotor are given in Tables 1, 2, 3 Fig. 13. The critical speed appears near 2240 rpm, which is close to the simulation result 2260 rpm. The correctness of the ER-rotor system model is verified. As concluded in Sect. 3.2, the critical speed will remain unchanged with the relatively low unbalanced mass and the critical speed will shift towards right with the high unbalanced mass, this phenomenon is demonstrated in the experimental study as shown in Fig. 13. It is worth mentioning that the purpose of the simulation and experimental study is to acquire some dynamic characteristics of the variable stiffness support rotor system with elastic rings. Simulation and test results show that the variable stiffness of the ER will cause the critical speed to increase under certain conditions.

Conclusion
Considering that the large inertia force and gyroscopic moment may occur during the operation of aeroengine, which may lead to contact between elastic ring and bearing pedestal. This paper proposes a variable stiffness model of elastic rings and this model is subsequently verified by comparison with simulation analysis and experimental results. Meanwhile, the  dynamic model of the ER-rotor is established, and the nonlinear time-varying stiffness of ER is introduced in the simulated analysis. By comparing the spectrum cascades and the amplitude frequency responses of the two ER-rotor system, the dynamic properties of the variable stiffness ER-rotor system are revealed.
Finally, the correctness of the ER-rotor system model and simulation results are verified by the experimental research on the ER-rotor test bench. The conclusions are as follows: 1. The numerical results show that only the fundamental frequencies are prominent in the ER-rotor  Fig. 12 a The ER-rotor test bench and b the combined support system with linear stiffness, but the ER-rotor system produces complex frequency components (1.1f n , 2f n and 3f n ) due to the variable stiffness of the elastic rings. 2. The variable stiffness of the elastic rings has no effect on the critical speed of the rotating system at lower unbalance mass. However, due to the large unbalance mass, the displacement of the ER-rotor system increases, which leads to contact deformation of the elastic ring, and the critical speed of the ER-rotor system increases accordingly. 3. The clearance of the elastic rings can reduce the critical speed of the system, and the critical speed of the ER-rotor system will be different with the increase of the unbalanced mass. 4. The research in this paper is helpful to understand the influence of variable stiffness of elastic rings on the dynamic characteristics of rotating system, and those characteristics can be used to detect the change of elastic ring stiffness in the ER-rotor system.
Author contribution KS investigated the article, summarized the concepts and methods, wrote the original draft, and verified it with the corresponding software. ZL wrote, reviewed, and edited the article, and provided supervision and analysis. LL carried out software verification and visualization of the article, and managed the data. JL carried out software verification and visualization of the article. FW carried out software verification and visualization of the article.