Abstract
Based on finite element (FE) method, a dynamic model of a single span bladerotorcasing coupling system is established. The bearings, bladed rotor, discs and casing are simulated using springdamping element, beam element, shell element and beam element suitable for curved beam separately (the casing stiffness and damping using springdamping element). The bladecasing rubimpact is simulated using pointpoint contact elements; here the corresponding nodes of the blades and the casing are identified as the contact points. To deal with contact constraint conditions and simulate bladecasing frictional characteristics, the augmented Lagrangian method and the coulomb friction model are adopted, respectively. The responses of the system are analyzed under two loading conditions by using variablestep Newmark$\beta $ integral method combined with NewtonRaphson iteration. The results show that the rubimpact at the second critical speed (case 2) is more serious than that at the first critical speed (case 1). 4× (× denotes the rotating frequency) and its edge frequency components, combination frequency components about 1× can be viewed as distinguishable rubimpact features for cases 1 and 2, respectively. Results can provide theoretical basis for rubimpact diagnosis in bladerotorcasing coupling systems.
1. Introduction
The distance between the rotating blade and the casing inwall is referred to as tip clearance in rotating machinery, which is one of the important parameters influencing the properties of the rotating machinery [1]. To improve the performance and efficiency of the rotating machinery, a key design for such improvements is minimizing the operational clearances between rotating blade and stationary casing. However, the probability of bladecasing rubimpact increases with the decreasing clearance, which may lead to a catastrophic failure. Most often rubimpact occurs at seals. Seldom, but more dangerous, a blade rubs against the stator or a vane [2] because the high tangential velocity, large impact energy, bending deformation easily relative to the rotor, which can lead to the local damage, such as the loss corner faults of the blade, high cycle fatigue at blade root [3], etc.
Rubimpact problems between the shaft/disc and the seal/casing have been studied for a long time and a large number of papers have been published [49]. The present theoretical researches mainly focus on the establishment of rubimpact force model, the dynamic characteristics of the rotor with the blades deformation due to rubimpact and the collision process simulation based on contact dynamics, etc.
For the researches on the rubimpact force model, assuming that the blade fixed at rotor (cantilevered) is flexible and the casing is rigid, Padovan et al. [10] developed a normal blade tip force model based on the mechanical energy balance in a simple configuration blade, and analyzed the nonlinear dynamic characteristics of the system under the condition of single and multiple blade rubimpact. Considering the effects of the bladed disc and rotary centrifugal force on the normal rubimpact force, Jiang et al. [11] developed a new rubimpact force model on the basis of the model of Padovan [10]. Kascak et al. [12] presented two kinds of rubimpact model, smearing and abradable rub models, based on different rubbing forms.
For dynamic characteristics of the rotor with the blades deformation due to rubimpact, assuming the blade as a uniform cantilever beam, Sinha [13] derived the dynamical equations of the rotorbladecasing system by considering the effect of the rotary inertia, gyroscopic moments, internal material damping in the shaft, external damping in the bearing, torque and axial forces and discussed the transient response of the rotor due to rubimpact during the acceleration and deceleration through the resonance. Lesaffre et al. [14] presented a model of fully flexible bladed rotor developed in the rotating frame and the system dynamical equations including gyroscopic effects, spin softening effects and the centrifugal stiffening effects were obtained by an energetic method. Considering the contact problem between the rotor and the stator as a static problem in some rotational speed ranges and assuming the turbo machine casing as an elastic ring, the complete problem of frictionless sliding contact between the blades and the casing was studied. For the particular condition of local rubimpact, some researchers adopted periodic pulse force to simulate the rubimpact force in the collision process of the blade and the casing [1517]. Sinha [15] presented a basic dynamical equation for a rotating radial cantilever Timoshenko beam clamped at the hub in a centrifugal force field, and then adopted a periodic pulse load to simulate the local rubimpact against the outer case, and discussed the transient responses of the beam with the tip deformation due to rubimpact in terms of the frequency shift and nonlinear dynamic response of the rotating beam. Turner et al. [16, 17] developed a simplified modeling approach to study the effect of pulsivetype loads on the stiffness behavior of turbomachinery airfoils and quantified the influence of airfoil flexibility based on the relationship between applied force duration and maximum tip deflection.
For the collision process simulation based on contact dynamics, some researchers adopted lumped mass, beam, shell and solid model to simulate the blade, curved beam, used thin wall cylindrical shell to simulate the casing, and reproduced the rubimpact process of the blade and the casing [3, 1822].
From the analysis of the above literatures, the existing researches seldom consider rotor, blade, casing as a total system, and the research on the effect of the bladecasing rubimpact on the vibration of the rotor system considering the coupling of the bladed rotor and the casing structure is also insufficient. The current paper analyzes the effect of bladecasing rubbing on the nonlinear dynamic characteristics of the flexible bladed rotor. Starting from establishing a dynamic model of a single span bladerotorcasing system based on the finite element method, the bladecasing rubimpact is simulated using pointpoint contact elements, where the corresponding nodes of the blades and the casing are identified as the contact points. To deal with contact constraint conditions and simulate bladecasing frictional characteristics, the augmented Lagrangian method and the coulomb friction model are adopted, respectively. The responses of the system are analyzed under two loading conditions by using variablestep Newmark$\beta $ integral method combined with NewtonRaphson iteration.
The paper is organized as follows: finite element model of the bladerotorcasing system is described in section 2. In section 3, numerical simulation of rubimpact between the blade and the casing is performed under two loading conditions. Finally, the conclusions are given in the last section.
2. Finite element model of the bladerotorcasing system
2.1. Rubimpact model based on contact theory
In order to study conveniently, the finite element model of the bladerotorcasing system is simplified according to the following assumptions:
(a) The blade and the corresponding casing are modeled by constantsection beam and curved beam, respectively;
(b) A relatively small sliding in the process of rubimpact is considered, which is similar to point contact between the blade tip and the casing. The rubimpact can be simulated by a pointpoint contact element with clearance, as is shown in Fig. 1.
For convenience, only part bladerotorcasing structure is depicted in Fig. 1. The master body is set as the blade and the slave one is the casing. $g$ is the gap function of the blade and the casing, ${F}_{N}$ and ${F}_{T}$ the normal and tangential rubimpact forces, respectively. A contact force $\left({F}_{N}\right)$ between two bodies is always compressive and is active only when the gap between the bladetip and the outer casing has closed. Assuming that the crosssection of the disc is perpendicular to the axis of rotation and the rubimpact occurs in the $yoz$ plane. $o\text{'}$ is the center of the shaft, $o$ the origin of the global coordinate system, $\omega $ whirl angular velocity, ${k}_{ys}$ and ${c}_{ys}$ the stiffness and the damping of the casing in $y$ direction, respectively.
In this paper the KuhnTucker conditions for Coulomb friction are used for contact detection. In order to ensure mutual nonpenetration of the contact boundary of the blade and the casing, the augmented Lagrangian method considering the friction law is selected to deal with contact constraint conditions [22, 23]. The method properly defines force normal to the contact surface to make the penetration limited to a specific tolerance. Such force is as follows:
where ${k}_{N}$ is the normal contact stiffness, i.e., the penalty factor, ${\lambda}_{N}^{(i+1)}$ the Lagrange multiplier for (i+1)th iteration, which can be found by:
where $\epsilon $ is the specific penetration tolerance. After specific iteration, if the penetration is still greater than $\epsilon $, the contact stiffness of the contact element will be augmented through the Lagrange multiplier. This procedure is repeated until the penetration is less than $\epsilon $.
Fig. 1Rubimpact schematic of the blade and the casing
2.2. Finite element model of the bladerotorcasing system with rubimpact
Some simplifications about the FE model of the bladerotorcasing system are introduced as follows:
(a) Shaft and discs are rigidly connected;
(b) Two bearings are identical and simulated ideally by linear stiffness and damping, in addition the cross terms are neglected;
(c) Movements in torsional and axial directions are negligible;
(d) The flexible casing is simulated by linear stiffness and damping. It is modeled by beam element.
Considering the bladecasing rubbing and the external force action, the equation of motion can be written as follows:
where $\mathbf{M}$, $\mathbf{G}$, $\mathbf{C}$, $\mathbf{K}$ and $\mathbf{u}$ respectively denote mass, gyroscopic, damping (including bearing, casing and the viscous damping) and stiffness (including bearing, shaft, disc, blade and casing stiffness) matrixes and displacement vector of the global system, $\mathbf{\lambda}$ is a vector about Lagrange multiplier, $\mathbf{B}$ the contact constraint matrix in the normal and tangential directions, ${\mathbf{g}}_{0}$ initial normal gap vector and ${\mathbf{F}}_{u}$ external load vector. In this paper, the Rayleigh damping form is adopted to determine the viscous part $\left({\mathbf{C}}_{s}\right)$ of the total damping matrix $\left(\mathbf{C}\right)$ and its expression can be written as follows:
where $\alpha $, $\beta $ denote the proportion coefficients of the mass matrix and stiffness matrix, respectively:
Here ${\omega}_{n1}$ and ${\omega}_{n2}$ denote the first and second critical speeds (r/min), ${\xi}_{1}$ and ${\xi}_{2}$ the first and second modal damping ratios.
The finite element model of bladerotorcasing system is given in Fig. 2. In the figure:
(a) Nodes 4 and 26 denote left and right journals, nodes 11 and 19 the connection points of the shaft and disc 1, disc 2;
(b) ${k}_{zl},{k}_{yl},{c}_{zl},{c}_{yl},{k}_{zr},{k}_{yr},{c}_{zr}$ and ${c}_{yr}$ denote the stiffness and damping of the bearings, where $k$ and $c$ denote the stiffness and the damping, subscripts $z$, $y$ and $l$, $r$ denote $z$, $y$ directions and left, right bearings, respectively;
(c) The coupling, shaft and blade are simulated by Beam188 element; disc by Shell181 element; the stiffness and the damping of the bearing and the casing by Combi214 element; the casing using Beam189 element which can take transverse shear effect and initial curvature of the beam into account, namely is suitable for simulating the curved beam; rubimpact between the blade and the casing by Conta178 element; the finite element mesh schematic and model parameters of discbladecasing system with rubimpact are shown in Fig. 3 and Table 1.
Fig. 2Finite element model of the bladerotorcasing coupling system
The shaft and the disc are rigidly connected by coupling DOF and the connection of the disc and the blade is using a sharing node. In order to analyze the influence of the bladecasing rubimpact on the vibration of the shaft and the casing, the vibration responses of node 11 (the shaft) and node 1844 (the casing) will be analyzed in the following sections. Variablestep Newmark$\beta $ integral method combined with NewtonRaphson iteration is used to solve the nonlinear differential equations considering contact, namely Eq. (3).
Fig. 3Finite element mesh schematic of the discbladecasing system with rubimpact
Table 1Model parameters of the bladerotorcasing system
Parts  Element type  Element number  Geometric parameters  Material parameters 
Coupling  Beam188  2  The details in literature [25]  $E=$ 210 GPa, $\upsilon =$ 0.3, $\rho =$ 7850 kg/m^{3}, ${\xi}_{1}={\xi}_{2}=$ 0.04. 
Shaft  Beam188  24  The details in literature [25]  
Bearing  Combi214  2  ${k}_{zl}={k}_{zr}=$ 200 MN/m, ${k}_{yl}={k}_{yr}=$ 500 MN/m, ${c}_{zl}={c}_{zr}={c}_{yl}={c}_{yr}$ $=$ 2 kN·s/m  
Disc  Shell181  384  The details in literature [25]  
Blade  Beam188  64  ${l}_{1}=$ 40 mm, ${w}_{1}=$ 15 mm, ${h}_{1}=$2 mm, $N=$ 16  
Casing  Beam189  40  ${w}_{2}=$ 15 mm, ${h}_{2}=$ 5 mm  
The stiffness/damping of the casing  Combi214  1  ${k}_{zs}={k}_{ys}=$ 1 MN/m, ${c}_{zs}={c}_{ys}=$ 1 kN·s/m  
Simulation element of the rubimpact  Conta178  8  ${k}_{N}=$ 80 MN/m 
3. Numerical simulation of the bladecasing rubimpact fault under two loading conditions
Based on the API Standard 617 [26] (Axial and Centrifugal Compressors and Expandercompressors for Petroleum, Chemical and Gas Industry Services, Seventh Edition), the two different unbalance loading conditions are determined according to the modal shape of system, as is shown in Fig. 4. In the figure, ${u}_{1}$ and ${u}_{2}$ denote the amount of unbalance of disc 1 and 2, respectively. For case 1, ${u}_{1}={u}_{2}=$ 1.56×10^{4} kg∙m and the phases of two unbalances are the same; for case 2, ${u}_{1}={u}_{2}=$ 1.56×10^{4} kg∙m and the phases of two unbalances are opposite. Geometric and material parameters of the bladerotorcasing system are shown in Table 1. Without rubimpact, unbalance responses of disc 1 (node 11) in $y$ direction under two loading conditions are shown in Fig. 5. From the figure, it can be seen that the first order resonant response is predominant in case 1, however the second order resonant response is predominant in case 2. Moreover ${\omega}_{n1}=$ 1500 r/min and ${\omega}_{n2}=$ 5700 r/min. It means that the first loading condition excites the first critical speed easily and the second loading condition excites the second critical speed easily. The maximum deformations of the bladerotor system at the first and second critical speeds are shown in Fig. 6. In the next sections, the numerical simulation for the rubimpact between the blade and the casing will be performed at two critical speeds under two loading conditions in consideration of large vibration amplitude of the rotor.
Fig. 4The schematic of two unbalance load conditions: a) case 1: inphase of two unbalances, b) case 2: outofphase of two unbalances
Fig. 5Unbalance responses of the rotordiscblade system
Fig. 6The maximum deformation of the rotordiscbladesystem: a) the maximum deformation of the rotordiscbladesystem at the first critical speed under case 1, b) the maximum deformation of the rotordiscbladesystem at the second critical speed under case 2
a)
b)
3.1. Dynamic characteristics of the rotordiscbladecasing system with rubimpact under case 1
Based on the unbalance response shown in Fig. 5, for case 1, two different initial normal gaps between the blade and the casing (${g}_{0}=$ 2.3 mm and ${g}_{0}=$ 1.9 mm) are selected to simulate light and serious rubimpact, respectively, at the first critical speed.
The system vibration responses are shown in Fig. 7 when the light bladecasing rubimpact occurs at ${g}_{0}=$ 2.3 mm. The bladecasing rubimpact creates an amplitude modulated response periodically, and the peak values are near ±2 mm, as is shown in Fig. 7(a). The vibration of the casing is similar to impulse responses caused by periodic rubimpact force, as is shown in Fig. 7(b). Fig. 7(c) shows the rotor orbit after rubimpact. The red point line denotes the initial clearance. It can be concluded from the whirl orbit that the collision of the multiblade against the casing leads to multiple rebounds of the rotor, which makes the rotor orbit similar to the multiple twining circles. Fig. 7(d) shows logarithmic amplitude spectra of the shaft (node 11) and the blade (node 39). Spectrum structure of the blade is more complicated than that of the shaft due to the appearance of more complicated high frequency components. The common feature of them is a wide continuous spectrum peak near 4×, which includes edge frequency components around 4×, such as the components of (4×–5 Hz) and (4×+5 Hz), however 1× is still predominant. The reason is that 4× is close to the second critical speed.
Fig. 7Vibration responses of the bladerotorcasing system with light rubimpact (g0= 2.3 mm): a) time domain waveform of the shaft (node 11), b) time domain waveform of the casing (node 1844), c) rotor orbit (node 11), d) amplitude spectra of the shaft (node 11) and the blade (node 39), e) normal rubimpact force, f) contact status
The normal rubimpact forces of eight contact elements are shown in Fig. 7(e), where the righthand abscissa is contact element serial number shown in Fig. 3. It can be seen from the figure that the normal rubimpact forces fluctuate regularly and the rubimpact periods are not completely the same for different contact elements, which also reflects the disorder of the rubimpact of multipleblade and the casing. The accurate rubimpact position and the degree of collision can be determined by normal contact force and contact time. The contact statuses of different contact elements with time are shown in Fig. 7(f). For the contact status, “1” denotes noncontact of the blade and the casing, “2” sliding contact and “3” sticking contact (no sliding). It can be concluded from Fig. 7(f) that the sliding and sticking contacts appear alternatively in the rubimpact process, where sticking contact appears at lager normal rubimpact force and sliding contact at smaller normal rubimpact force.
The system vibration responses are shown in Fig. 8 when the serious bladecasing rubimpact occurs at ${g}_{0}=$ 1.9 mm. The vibration displacement decreases due to the limitation of the casing, as is shown in Fig. 8(a). The displacement of the stator increases slightly and the collision times have a bigger increase compared with those under the light rubimpact condition, as is shown in Fig. 8(b). The rotor orbit shown in Fig. 8(c) is similar to that under the light rubimpact condition. 1/2 fractional harmonic components, such as 3×/2, 7×/2, etc., can be observed, and edge frequency components around 4×, such as (4×–5 Hz) and (4×+5 Hz) also appear in amplitude spectra shown in Fig. 8(d). The rubimpact is more serious due to the increase of collision times and normal rubimpact force compared with those under the light rubimpact condition, as is shown in Fig. 8(e). Sticking contact increases and sliding contact decreases when the bladecasing clearance decreases, as is shown in Fig. 8(f).
Fig. 8Vibration responses of the bladerotorcasing system with serious rubimpact (g0= 1.9 mm): a) time domain waveform of the shaft (node 11), b) time domain waveform of the casing (node 1844), c) rotor orbit (node 11), d) amplitude spectra of the shaft (node 11) and the blade (node 39), e) normal rubimpact force, f) contact status
3.2. Dynamic characteristics of the rotordiscblade casing system with rubimpact under case 2
For case 2, the same initial normal gaps between the blade and the casing as case 1 are selected to simulate light and serious rubimpact at the second critical speed. Fig. 9(a) shows that the fluctuation of the vibration amplitude is more violent. The displacement of the stator at case 2 is about 4 times of that at case 1, which shows that the vibration at case 2 is more violent, as is shown in Fig. 9(b). From the rotor orbit shown in Fig. 9(c), it can be seen that the rotor vibration is beyond the clearance of the blade and the casing, which shows that the rotor motion is forward whirl in an annular region and the casing vibrates seriously due to the collision. Combination frequency components about the rotating frequency (1×+95 Hz) appear in the amplitude spectra shown in Fig. 9(d), such as a low frequency component (43 Hz) less than 1× and a high frequency component (147 Hz) greater than 1×, and the sum of low frequency and high frequency is equal to 2×, which shows that the frequency components caused by rubimpact are related to 1×. It can be observed from Figs. 9(e) and 9(f) that the normal rubimpact force (about 1000 N) at case 2 is about five times of that at case 1, the collision is more frequent compared with that at case 1 and the contact status is mainly sticking contact.
Fig. 9Vibration responses of the bladerotorcasing system with light rubimpact (g0= 2.3 mm): a) time domain waveform of the shaft (node 11), b) time domain waveform of the casing (node 1844), c) rotor orbit (node 11), d) amplitude spectra of the shaft (node 11) and the blade (node 39), e) normal rubimpact force, f) contact status
The vibration responses of the system with serious rubimpact are shown in Fig. 10. The fluctuation of the time domain waveforms shown in Fig. 10(a) increases and the stator vibration shown in Fig. 10(b) is more regular compared with that under light rubimpact condition at case 2. For the rotor orbit shown in Fig. 10(c), it is clear that the orbit is similar under light and serious rubimpact conditions at case 2, which also indicates that the elastic casing can absorb impact energy and weaken the rebound energy acting on the blade. Combination frequency components about 1× appear in the amplitude spectra shown in Fig. 10(d). The low frequency component changes to 44 Hz, the high frequency component 146 Hz and the sum of both is still equal to 2×, which also indicates that combination frequency components can change with the different rubimpact intensity. Compared with light rubimpact, the normal rubimpact force shown in Fig. 10(e) increases slightly. Contact status shown in Fig. 10(f) is similar to that under light rubimpact condition and mainly sticking contact under serious rubimpact condition.
Fig. 10Vibration responses of the bladerotorcasing system with serious rubimpact (g0= 1.9 mm): a) time domain waveform of the shaft (node 11), b) time domain waveform of the casing (node 1844), c) rotor orbit (node 11), d) amplitude spectra of the shaft (node 11) and the blade (node 39), e) normal rubimpact force, f) contact status
Table 2Rubimpact feature comparison under two loading conditions
Conditions  Rub impact intensity  Maximum displacement of the casing in y direction (μm) (node 1844)  Typical spectrum features of the shaft (node 11) and the blade (node 39)  Normal minimum/maximum contact force (N) (No. 1 contact element)  Contact status (No. 1 contact element) 
Case 1  Light  51  4× and its edge frequency components  Minimum: 7.0 Maximum: 189.3  Sliding contact: 53 times Sticking contact: 14 times 
Serious  52  3×/2,7×/2, 4× and its edge frequency components  Minimum: 6.9 Maximum: 195.7  Sliding contact: 80 times Sticking contact: 30 times  
Case 2  Light  199  Combination frequency components about 1×, such as 43 Hz, 147 Hz, 2×  Minimum: 8.1 Maximum: 1042.5  Sliding contact: 168 times Sticking contact: 58 times 
Serious  204  Combination frequency components about 1×, such as 44 Hz, 146 Hz, 2×  Minimum: 6.6 Maximum: 1079.6  Sliding contact: 169 times Sticking contact: 51 times 
3.3. Fault features comparison under two loading conditions
The detailed fault features including the maximum vibration displacement of the casing in $y$ direction, the spectrum features of the shaft and the blade, normal minimum/maximum contact forces and contact status are listed in Table 2 under two loading conditions. In the table, times of sliding contact are more than that observed from the contact status figures (see Fig. 7(f), 8(f), 9(f), 10(f)) because that the many points indicating sliding status are closed to overlap, which also shows that a slight bounce of the blade will appear when intensive sliding status appears. The fewer sticking contact status and shorter contact time show that the rebound of the blade aggravates. For case 1, there are some differences between the contact status and spectrum features but fewer differences in other aspects under light and serious rubimpact conditions. For case 2, the fault features are all similar under the light and serious rubimpact conditions. Compared with case 1, it is clear that the rubimpact has a great influence on the vibration of the rotor and the casing due to the larger normal rubimpact force and more collision times at case 2, which can be proven by the vibration of the stator.
4. Conclusions
In this paper, nonlinear vibration characteristics of the rubimpact are investigated by developing the bladerotorcasing finite element model. In the analytical discussions of the transient response and contact status of the light and serious rubbing under two loading conditions, the conclusions can be obtained as follows:
(1) For the first loading condition which mainly arouses the first order bending resonance, 4× and its edge frequency components appear and the contact status of the blade and the casing is mainly sliding contact under light rubimpact condition; under serious rubimpact condition, 3×/2, 7×/2, 4× and its edge frequency components appear and the collision times and sticking contact times have a great increase than those of light rubimpact.
For the second loading condition which mainly arouses the second order bending resonance, the combination frequency components about 1× appear and the combination frequency value changes slightly with increasing rubimpact intensity.
(2) The rubimpact of the blade and the casing has a small influence on the displacement of the casing, rotor whirl orbit and normal rubimpact force under the first loading condition. Compared with the first loading condition, the rubimpact is more dangerous under the second loading condition because of the larger normal rubimpact force, many collision times, the lower frequency components with larger amplitude and excessive stator vibration.
References

Mehalic C. M., Ziemianski J. A. Performance Deterioration of Commercial HighBypass Ratio Turbo Fan Region. ASME Technical Paper Series 80118, 1980.

Muszynska A. Rotortostationary element rubrelated vibration phenomena in rotating machinery – literature survey. Vibration Inst., The Shock and Vibration Digest, Vol. 21, Issue 3, 1989, p. 311.

Liu S. G., Hong J., Chen M. Numerical simulation of the dynamic process of aeroengine blade to case rub impact. Journal of Aerospace Power, Vol. 26, Issue 6, 2011, p. 12821288.

Muszynska A. Rotordynamics. Boca Raton, CRC Press, 2005.

Ahmad S. Rotor casing contact phenomenon in rotor dynamics – literature survey. Journal of Vibration and Control, Vol. 16, Issue 9, 2010, p. 13691377.

Wen B. C., Wu X. H., Ding Q., et al. Theory and Experiment of Nonlinear Dynamics for Rotating Machinery with Faults. Science Press, 2004.

Gao W., Zhang X. J., Zhang Y. Review of nonlinear rotordynamics. Journal of Southeast University (Natural Science Edition), Vol. 32, Issue 3, 2002, p. 19.

Chen Y., Zhang H. Review and prospect on the research of dynamics of complete aeroengine systems. Acta Aeronautica et Astronautica Sinica, Vol. 32, Issue 8, 2011, p. 13711391.

Chen G., Li C. G., Wang D. Y. Nonlinear dynamic analysis and experiment verification of rotorball bearingssupportstator coupling system for aeroengine with rubbing coupling faults. Journal of Engineering for Gas Turbines and Power, Vol. 132, Issue 2, 2010, p. 022501.1022501.9.

Padovan J., Choy F. K. Nonlinear dynamics of rotor/blade/casing rub interactions. Journal of Turbomachinery, Vol. 109, Issue 4, 1987, p. 527534.

Jiang J., Ahrens J., Ulbrich H., et al. A contact model of a rotating, rubbing blade. Proceedings of the 5th International Conference on Rotor Dynamics of the IFTOMM, 1998, p. 478489.

Kascak A. F., Tomko J. J. Effects of Different Rub Models on Simulated Rotor Dynamics. National Aeronautics and Space Administration, Cleveland, OH, Lewis Research Center, 1984.

SinhaS. K. Dynamic characteristics of a flexible bladedrotor with Coulomb damping due to tiprub. Journal of Sound and Vibration, Vol. 273, 2004, p. 875919.

Lesaffre N., Sinou J. J., Thouverez F. Contact analysis of a flexible bladedrotor. European Journal of Mechanics – A/Solids, Vol. 26, Issue 3, 2007, p. 541557.

Sinha S. K. Nonlinear dynamic response of a rotating radial Timoshenko beam with periodic pulse loading at the free end. International Journal of NonLinear Mechanics, Vol. 40, Issue 1, 2005, p. 113149.

Turner K., Adams M., Dunn M. Simulation of engine blade tiprub induced vibration. Proceedings of GT2005, RenTahoe, Nevada, USA, 2005.

Turner K., Dunn M., Padova C. Airfoil deflection characteristics during rub events. Journal of Turbomachinery, Vol. 134, Issue 1, 2012, p. 01101810110187.

Legrand M., Pierre C., Peseux B. Structural modal interaction of a four degree of freedom bladed disk and casing model. Journal of Computational and Nonlinear Dynamics, Vol. 5, Issue 4, 2010, p. 1341.

Batailly A., Legrand M., Cartraud P., et al. Assessment of reduced models for the detection of modal interaction through rotor stator contacts. Journal of Sound and Vibration, Vol. 329, Issue 26, 2010, p. 55465562.

Legrand M., Batailly A., Pierre C. Numerical investigation of abradable coating removal in aircraft engines through plastic constitutive law. Journal of Computational and Nonlinear Dynamics, Vol. 7, Issue 11, 2011, p. 011010.1011010.11.

Legrand M., Batailly A., Magnain B., et al. Full threedimensional investigation of structural contact interactions in turbomachines. Journal of Sound and Vibration, Vol. 331, Issue 11, 2012, p. 25782601.

Roques S., Legrand M., Cartraud P., et al. Modeling of a rotor speed transient response with radial rubbing. Journal of Sound and Vibration, Vol. 329, Issue 5, 2010, p. 527546.

Ma H., Shi C. Y., Han Q. K., et al. Fixedpoint rubbing fault characteristic analysis of a rotor system based on contact theory. Mechanical Systems and Signal Processing, Vol. 38, Issue 1, 2013, p. 137153.

Simo J. C., Laursen T. A. An augmented Lagrangian treatment of contact problems involving friction. Computers and Structures, Vol. 42, Issue 1, 1992, p. 97116.

Ma H., Tai X. Y., Sun J., et al. Analysis of dynamic characteristics for a dualdisk rotor system with single rubimpact. Advanced Science Letters, Vol. 4, Issues 810, 2011, p. 27822789.

API 617. Axial and Centrifugal Compressors and Turboexpanders for Petroleum, Chemical and Gas Industry Services. American Petroleum Institute, Washington, D. C., 2002.
About this article
We are grateful to the Fundamental Research Funds for the Central Universities (Grant No. N120203001) and Program for New Century Excellent Talents in University (Grant No. NCET110078) for providing financial support for this work.