# Theoretical and Experimental Study on the Transient Time-Frequency Characteristics of the Bending-Torsional Coupling Motions of a Rub-Impact Dual-Rotor System

Authors: Tian Gao; Shuqian Cao; Tiancheng Zhang
Journal: Shock and Vibration
 (2021)

---

## Abstract

This paper focuses on the fault characteristics of the bending and torsional motions of a rub-impact dual-rotor system caused by aircraft flight maneuvers. The equations of the bending-torsional coupling motion of a dual-rotor system are established considering a low-pressure rotor rub-impact fault and the transient barrel roll flight of an aircraft. The 4th Runge-Kutta method with varied steps is used to obtain the bending and torsional responses. Then, the influences of the system parameters, including the rub-impact stiffness, friction coefficient, and rotating speed, on the bending and torsional motions of the dual-rotor system are investigated in detail. At last, a rotor rubbing experiment is carried out, verifying the validity of the simulation results. The results show that the rub-impact stiffness affects bending vibration significantly and the torsional motion is sensitive to the friction coefficient. Correspondingly, the torsional responses show apparent fractional fault frequencies and rotating fault frequencies within the whole region of the rub-impact stiffness. The bending responses can only display fault frequencies at certain rub-impact stiffness. As for the rotating speed, the torsional responses are also more effective than the bending responses for the rub-impact fault detection at the low- and high-speed regions. The results will contribute to a comprehensive basis for the rub-impact fault detection.


---

## Body

## 1. Introduction

The rub-impact of rotor-to-stator, which is one of the most common faults in aeroengine rotor systems, could lead to catastrophic failure sometimes, such as blade fracture and engine flameout [1, 2]. Therefore, the early fault detection of rub-impact is crucial for avoiding the disasters. Many researchers have made a profound study on the periodic motions [3, 4], chaotic-responses [5, 6], and harmonic frequency characteristics [7–9] of the bending responses of the rub-impact rotor system. However, under the condition of the rub-impact, the torsional responses also displayed abundant fault dynamic responses because of the friction torque caused by the rub-impact of the rotor-to-stator [10]. Thus, it is important to investigate the fault characteristics of the bending and torsional motions of the rub-impact rotor system, so as to provide a comprehensive theoretical basis for the state detection and fault diagnosis.The mass imbalance of the rotor system is one of the major factors leading to the bending-torsional coupling motions. The analyses of the bending-torsional coupling responses, caused by the mass imbalance, were mainly aimed at the slender flexible rotor, such as turbogenerator unit.. He et al. [11] built the differential equations of the bending-torsional coupling vibration for a Jeffcott rotor model with unbalance mass. They found that the combined resonance of the coupled bending-torsional motion happened when the rotating frequency of the rotor system was equal to the sum (or the difference) of the natural frequencies of the bending and torsional vibration. Li and Chen [12] studied the subsynchronous resonance and combined resonance of the bending and torsional motions of a low-pressure cylinder-generator rotor system through the singularity theory for the two-state variable system. The transition sets and bifurcation diagrams were obtained which could be used for rotor system optimization. Recently, much attention has been paid to the bending-torsional coupling motions of the fault rotor systems. The effect of the parallel misalignment and angular misalignment on the bending and torsional responses of two rotating shafts was examined through the theoretical and numerical analyses [13, 14]. It was concluded that the parallel misalignment could lead to the bending-torsional coupling motions and the angular misalignment increased the volume of the instability region. Darpe et al. [15, 16] investigated the coupling between the transversal, torsional, and longitudinal vibration of a transverse cracked rotor, in which the excitation in one mode resulted in an interaction between all the modes. Patel and Darpe [17] studied the coupled bending-torsional vibrations of a rotor system with rub and crack. The rub excitations and breathing influences of the crack affected the rotor response in both bending and torsional modes. Lin and Chu [18] carried out the dynamic analysis of a rotor system with a slant crack on the shaft. They found that there were coupling stiffness of the bending-torsion, bending-tension, and torsion-tension for a slant-cracked shaft. Wang and Jiang [19] derived the motion equations of a dual-rotor system considering unbalance-misalignment coupling faults. The complicated vibration responses affected by rotational speed, mass eccentricity, misalignment angle, and parallel misalignment were investigated.Regarding the bending-torsional motions of the rub-impact rotor system, Lu and Chu [20] found they consisted of the similar frequency components under the conditions of no rub, full annual rub, and partial rub. Yang et al. [21] studied the bending-torsional coupled vibration of an aeroengine model subjected to blade-casing rub in presence of nonuniform initial gap. The results showed that the nonuniform initial gap led to rich vibration phenomena in the dynamic responses, especially in the torsional responses. Zhang and Ding [22, 23] investigated the influence of rub-impact on the bending and torsional vibrations of a two-disk rotor-stator system. The results indicated that the spectrum analyses of the torsional vibrations were more suitable to determine the axial contact-rubs. Xu et al. [24] studied the vibration frequency spectrum characteristics of a rub-impact dual-rotor system based on a finite element model. Zhang et al. [25] analyzed the coupled bending and torsional motions of a rub-impact dual-rotor system through numerical simulation and experiment. The bending and torsional vibrations both contained multiple frequency, fractional frequency, and combined frequency of the rotating frequency and multiple/fractional frequency.The maneuver load during flight maneuvers of an aircraft could lead to a rub-impact phenomenon of the aeroengine rotor system, accompanied with complex nonlinear behavior [4, 26–28]. However, the rub-impact fault characteristics of the bending vibration do not behave evidently under some certain system parameters. In this study, the equations of the bending-torsional coupling motion of a rub-impact dual-rotor system with maneuver load are established by using Lagrange equation firstly. The transition of the rotor system from nonrubbing to rub-impact, caused by the transient maneuvering flight, is obtained through the numerical integration method. Then, the time-frequency features of the bending and torsional vibrations of the rub-impact dual-rotor system are investigated during the flight of an aircraft. A variety of comparisons for the rub-impact fault features of the bending and torsional motions are observed under varied rub-impact stiffness, friction coefficients, and rotating speeds. At last, a rotor rubbing experiment, under maneuvering flight, is carried out, verifying the validity of the simulation results. The research of this study provides a deeper and more comprehensive theoretical basis for the state detection and fault diagnosis of the aeroengine rotor system during flight maneuvers.


## 2. Mathematical Modeling

### 2.1. Geometrical Model of Dual-Rotor System

AL-31F is one of the most common aero-engines. Based on the basic structure of the AL-31F, the aeroengine rotor system is simplified as a dual-rotor model, which consists of a low-pressure rotor and a high-pressure rotor, as shown in Figure1. The low-pressure rotor adopts the three-point supporting structures in which the mounting positions are at the 1#, 2#, and 5#. The high-pressure rotor absorbs the two-point supporting schemes in which the mounting positions are at the 3# and 4#. The low-pressure rotor and the high-pressure rotor are coupled through the intermediate ball bearing at 4#. The components of the low-pressure rotor system include the lumped mass of the disks m2 and m5 and the lumped mass of the shaft m1, m3, m4, and m6. The high-pressure rotor is regarded as a rigid rotor which only consists of the lumped mass m7. All the supporting structures, except the intermediate ball bearing 4#, are considered as elastic supporting systems with elastic supports and rolling bearings in series. They are simplified to spring-damping systems, wherein the stiffness and damping of the elastic supports are defined as k1, k2, k3, k5, c1, c2, c3, and c5, respectively. All the disks of the dual-rotor system have an unbalance magnitude and the eccentricity of the disks are described by ei (i = 2, 5, and 7) separately. In addition, the rub-impact fault occurs at the disk 2.Figure 1
Schematic diagram of the dual-rotor system.


### 2.2. Equation of Motion of Dual-Rotor System with Rub-Impact Fault in Maneuvering Flight

Lagrange’s equation is used to derive the equations of motion of the dual-rotor system. It is defined as [29](1)ddt∂L∂q˙j−∂L∂qj+∂D∂q˙j=Qj,j=1,2,…,k,where L = T-V, T and V are the kinetic energy and potential energy functions of the dual-rotor system, D is the dissipative energy function, qj is the generalized coordinate in the jth degree of freedom, and Qj is the generalized force corresponding to qj. The process of the derivation is as follows. First, the kinetic energy, potential energy and dissipative energy of the dual-rotor system are obtained. Then, the generalized forces of the ball bearing, rub-impact fault, gravity, and maneuver load are given. At last, the equations of the bending-torsional coupling motion are calculated by inserting the kinetic energy, potential energy, dissipative energy, and generalized forces into Lagrange’s equation given in equation (1).


#### 2.2.1. Energy Functions of Dual-Rotor System

Based on the modeling method of Reference [25], the following assumption is made for the dual-rotor system. As we know, the low-pressure rotor is a relatively flexible slender shaft. The high-pressure rotor is short and the shaft diameter of which is large. Because the short and thick shaft is harder to display large torsional angle, the torsional vibrations of the high-pressure rotor are neglected. The bending and torsional motions are only considered at the disk mi (i = 2, 5) of the low-pressure rotor. They are described by the generalized coordinates xi,yi,αi (i = 2, 5), where xi and yi are the centroid displacements of the disk and αi is the torsional angle. Meanwhile, the torsional angle satisfies the following function φi=ωlowt+αi+ϕi (i = 2, 5), where φi is the phase angle of the mass centroid of the disk, ωlow is the rotating speed of the low-pressure rotor, and ϕi is the initial phase angle. Regarding the high-pressure rotor mi (i = 7) and the lump mass of the shaft of the low-pressure rotor mi (i = 1, 3, 4, 6), the bending motions are taken into account simply. They are defined as (xi, yi) (i = 1, 3, 4, 6, 7), where xi and yi are the centroid displacements of the high-pressure rotor and the shaft of the low-pressure rotor, respectively. In addition, the rotating speed of the high-pressure rotor is described by ωhigh, which satisfies ωhigh=γωωlow (γω is the rotating speed ratio of the high-pressure rotor to low-pressure rotor). Therefore, there are 16 degrees of freedom of the dual-rotor system which consist of 14 bending motions and 2 torsional motions. They are given by the generalized vector q:(2)q=x1,x2,…,x7,y1,y2,…,y7,α2,α5T,where the bending generalized vector x, y and the torsional generalized vector α are separately expressed as(3)x=x1x2x3x4x5x6x7,y=y1y2y3y4y5y6y7,α=α2α5.Because there are eccentricitiesei (i = 2, 5, 7) of the rotating disks, the mass centroid of the high-pressure rotor and the low-pressure rotor is described by the generalized vectors xc and yc:(4)xc=x1x2+e2cosφ2x3x4x5+e5cosφ5x6x7+e7cosφ7,yc=y1y2+e2sinφ2y3y4y5+e5sinφ5y6y7+e7sinφ7.The kinetic energyT of the dual-rotor system consists of translational kinetic energy Tt and rotational kinetic energy Tr. It is defined as(5)T=Tt+Tr=12x˙cTMcx˙cT+12y˙cTMcy˙cT+12φ˙Jcφ˙,where(6)x˙c=x˙1x˙2−e2φ˙2sinφ2x˙3x˙4x˙5−e5φ˙5sinφ5x˙6x˙7−e7φ˙7sinφ7,y˙c=y˙1y˙2+e2φ˙2cosφ2y˙3y˙4y˙5+e5φ˙5cosφ5y˙6y˙7+e7φ˙7cosφ7,φ˙=ωlow+α˙2ωlow+α˙5ωhigh,Mc=m10000m20000⋱0000m7,Jc=J2+m2e22000J5+m5e52000J7+m7e72.where x˙c and y˙c are the translational velocity vectors of the mass centroid, φ˙ is the rotating velocity vector of the disks, Mc is the mass matrix, and Jc is the moment of inertia relative to the mass centroid of the disks. Moreover, Ji = 1/2miRi2 (i = 2, 5, 7), where Ri are the radii of the disks.The potential energyUe of the dual-rotor system consists of the elastic potential energy of the shaft and the elastic potential energy of the elastic supports. It is defined as(7)Ue=12xTKxx+12yTKyy+12αTKtα,where(8)Kx=Ky=k1+k6−k600000−k6k6+k7−k700000−k7k2+k7+k8−k800000−k8k8+k9−k900000−k9k9+k10−k1000000−k10k5+k100000000k3,Kt=kt1+kt2−kt2−kt2kt2,where Kx and Ky are the bending stiffness matrices of the x- and y-directions, respectively, and Kt is the torsional stiffness matrix. In addition, the bending stiffness and torsional stiffness of the shaft are given by(9)k6=3EIl13,k7=12EIl23,k8=12EIl33,k9=12EIl43,k10=3EIl53,kt1=GIPl1,kt2=GIPl2+l3+l4,where k1, k2, k3, and k5 are the stiffness of the elastic supports, the values of which are obtained from the elastic stiffness of a squirrel cage, k6, k7, k8, k9, and k10 are the bending stiffness corresponding to the shaft segment l1, l2, l3, l4, and l5, respectively, kt1 and kt2 are the torsional stiffness corresponding to the shaft segment l1 and l2 + l3 + l4, separately, E and G are Young’s modulus and the shear modulus of the shaft, and I and IP are the moment of inertia and the polar moment of inertia of the shaft.The dissipative energy of the dual-rotor system can be described by the Rayleigh energy dissipation function:(10)D=12x˙TCxx˙+12y˙TCyy˙+12α˙TCtα˙,where(11)Cx=Cy=c1+c6−c600000−c6c6+c7−c700000−c7c2+c7+c8−c800000−c8c8+c9−c900000−c9c9+c10−c1000000−c10c5+c100000000c3,Ct=ct1+ct2−ct2−ct2ct2,where Cx, Cy, and Ct are the damping matrices corresponding to the bending and torsional motions, respectively. c1, c2, c3, and c5 are the given damping coefficients of the elastic supports, c6, c7, c8, c9, and c10 are the given damping coefficients corresponding to the bending motion of the shaft segments l1, l2, l3, l4, and l5, separately, and ct1 and ct2 are the damping coefficients corresponding to the torsional motions of the shaft segment l1 and l2+ l3+ l4, respectively. In this study, the damping coefficients ci (i = 1, 2, …, 10) are directly given, based on the Rayleigh damping ratio of the rotor system valuing around 0.1. Then, the initial free vibration of the simulation model can be reduced quickly.


#### 2.2.2. Generalized Force of Intermediate Ball Bearing

The model of the intermediate ball bearing is shown in Figure2. The inner race of the ball bearing is directly fixed on the low-pressure rotor and the outer race is directly joined together with the high-pressure rotor. We assume that the balls are arranged equidistantly at the cage of the bearing. Then, the rotating speed of the cage of the bearing can be defined as(12)ωc=ωlowRinner+ωhighRouterRinner+Router,where Rinner and Router are the radii of the inner race and outer race of the ball bearing, respectively, and ωlow and ωhigh are the rotating speeds of the low-pressure rotor and high-pressure rotor, separately.Figure 2
Intermediate ball bearing model.The number of the balls is described byNb. Then the location of the jth ball is given by(13)θj=ωct+2πj−1Nb,j=1,2,…,Nb,The relative displacement between the inner and outer races of the ball bearing is expressed by (x7-x4, y7-y4). The initial bearing clearance is defined as δ0. Then, the contact deformation of the jth ball is expressed as(14)δj=x7−x4cosθj+y7−y4sinθj−δ0.According to the Hertzian contact theory [30], the force of the intermediate ball bearing is given by(15)FbxFby=Cb∑j=1Nbδj3/2Hδjcosθjsinθj,where Cb is the Hertz contact stiffness and H( ) is the Heaviside function.


#### 2.2.3. Force of Rub-Impact Fault

The schematic diagram of the rub-impact fault is shown in Figure3. The relative displacement between the disk and the stator is described by r2=x22+y221/2. The initial clearance between the disk and the stator is defined as δ.Figure 3
Rub-impact fault model.When the rub-impact fault happens, there are normal forceFn and tangential force Fτ applied to the rotor system. According to [31], the rub-impact force can be expressed by(16)FrxFry=−krr2−δr21−μμ1x2y2,r2>δ,where kr and μ are the contact stiffness and the friction coefficient between the disk and stator, respectively. In addition, the friction applied to the disk leads to friction torque which is given by(17)Mr=μkrR2r2−δ,r2>δ,where R2 is the radius of the disk. Moreover, from equation (17), we get that that the friction torque is the function consisting of the bending displacement r2. Meanwhile, we know that the friction torque is the main factor resulting in the torsional motion. Therefore, the friction torque is one of the major factors leading to the bending-torsional coupling motion.


#### 2.2.4. Maneuver Load and Gravity

As shown in Figure4, the influence of the barrel roll flight of an aircraft on the rub-impact dual-rotor system is investigated. The process of the transient barrel roll flight is simplified as three stages: at first (t∈0s,2s), the aircraft keeps the steady-state flight, and it enters into the barrel roll flight at t = 2 s and quits the barrel roll flight at t = 4 s. Then, it returns to the steady flight again (t∈4s,6s).Figure 4
Barrel roll flight of an aircraft.The actual barrel roll flight is very complicated. We have to make some simplified assumptions to study its influence. As shown in Figure5, the barrel roll flight is considered as a circular motion. The rolling speed and rolling radius of the aircraft are defined as ωbz and Rbz, respectively. In addition, because the rotation of the aircraft is faster than the revolution during the barrel roll flight, the rotation speed of the aircraft is defined as 2ωbz. To express the load state of the rotor system conveniently, two Cartesian coordinated references are used. The inertial reference frame O-XY, fixed to the ground is used for describing the flight status of the aircraft. The generalized reference frame o-xy, fixed to the aircraft, is used to describe the vibration of the dual-rotor system. The dual-rotor system is subjected to gravity mig and additional inertial force miωbz2Rbz due to the barrel roll flight of the aircraft. The gravity is always vertically downward. The additional inertial force is directed to the center of the rolling circle. Because the rolling direction of the aircraft is the same as the rotating direction of the dual-rotor system, there is no additional inertia moment applied to rotor system. Then, the generalized forces of the maneuver load and gravity are obtained by decomposing them into the generalized reference frame o-xy:(18)FgxiFgyi=miωbz2Rbzcosθbz−migsin2θbz−miωbz2Rbzsinθbz−migcos2θbz,where θ˙bz=ωbzt∈2s,4s.Figure 5
Force analysis of the rotor system during the barrel roll flight of an aircraft.


#### 2.2.5. Equation of Bending-Torsional Coupling Motion of Dual-Rotor System

The equations of the bending-torsional coupling motion of the rub-impact dual-rotor system can be calculated by substituting the energy functions and generalized forces into equation (1):(19)Mq¨+Cq˙+Kq=Fr+Fb+Fu+Fg,where(20)M=Mc000Mc000Jc′,C=Cx000Cy000Ct,K=Kx000Ky000Kt,Jc′=J2+m2e2200J5+m5e52,Fr=0Frx000000Fry00000Mr0T,Fb=000Fbx00−Fbx000Fby00−Fby00T,Fu=0Fux200Fux50Fux70Fuy200Fuy50Fuy7Fuα1Fuα2T,Fg=Fgx1Fgx2Fgx3Fgx4Fgx5Fgx6Fgx7Fgy1Fgy2Fgy3Fgy4Fgy5Fgy6Fgy700T,Fux2=m2e2φ˙22cosφ2+φ¨2sinφ2,Fux5=m5e5φ˙52cosφ5+φ¨5sinφ5,Fux7=m7e7ωh2cosωht,Fuy2=m2e2φ˙22sinφ2−φ¨2cosφ2,Fuy5=m5e5φ˙52sinφ5−φ¨5cosφ5,Fuy7=m7e7ωh2sinωht,Fuα1=m2e2x¨2sinφ2−y¨2cosφ2,Fuα2=m5e5x¨5sinφ5−y¨5cosφ5,where M, C, and K are the generalized mass matrix, damping matrix, and stiffness matrix, respectively, Fr and Fb are the generalized force of the rub-impact fault and intermediate ball bearing, separately, Fu is the unbalanced force, and Fg is the generalized resultant force of the maneuver load and gravity, φi (i = 2, 5) are the time-varying phase angles of the mass centroid of the disks, which have been defined in the Section 2.2.1. In addition, from the derivation of the Fu, we know that unbalanced force is another factor leading to the bending-torsional coupling motion.


### 2.3. Simulation Parameter and Calculation Process

Because of the strong nonlinearity of the rotor system, the numerical integration method is adopted to solve the (19). The bending-torsional coupling dynamic responses of the rub-impact dual-rotor are calculated via the 4th Runge-Kutta method with varied step. According to Reference [25], the simulation parameters of the dual-rotor system are listed in Table 1.Table 1
Parameters of dual-rotor model.
Physical descriptionParameterValueMassm1, m2, m3, m4, m5, m6, m7 (kg)2, 10, 2, 2, 6, 2, 8Radius of the diskR2, R5 (m)0.075, 0.03Initial phase angle of the diskϕ2,ϕ5 (rad)0, 0Eccentricity of the diske2, e5, e7μm46.35, 46.35, 47.5Stiffness of the elastic supportk1, k2, k3, k5 (N/m)106, 106, 106, 106Length of the shaft segmentl1, l2, l3, l4, l5 (m)0.33, 0.33, 0.6, 0.23, 0.23Material properties of the shaftE, G (Pa)2 × 1011, 7.7 × 1010Physical properties of the shaftI, Ip (m4)1.178 × 10−7, 2.356 × 10−7Damping of the elastic supportc1, c2, c3, c5 (N.s/m)300, 300, 300, 300Damping of the bending motionc6, c7, c8, c9, c10 (Ns/m)150, 150, 150, 150, 150Damping of the torsional motionct1, ct2 (N.m.s/rad)1, 1Parameters of the rub-impact faultkr (N/m), δ(m), μ8 × 106, 1.34 × 10−4, 0.1Parameters of the ball bearingCb (N/m3/2), δ0,Rinner,Router(mm)1.334 × 1010, 0.005, 40.1, 63.9, 8Rotating speeds of the rotor systemωlow,ωhighrad/s,γω1500, 2700, 1.8Parameters of the maneuvering flightωbzrad/s, Rbz (m)π, 3In the process of the dynamic analysis, the relative torsional angleΔα=α2−α5 is chosen to show the torsional characteristics of the dual-rotor system. The bending displacements (x2 and y2) of the disk m2 are selected to reveal the bending characteristics. In addition, the motion of the dual-rotor system consists of some responses with low frequency and large amplitude because of the time-varying influence of the maneuvering flight. They are useless for the analysis of the rub-impact fault features. Therefore, a high-pass filter is applied to the response of the dual-rotor system. The detailed calculating process is shown in Figure 6.Figure 6
Calculating process.


## 2.1. Geometrical Model of Dual-Rotor System

AL-31F is one of the most common aero-engines. Based on the basic structure of the AL-31F, the aeroengine rotor system is simplified as a dual-rotor model, which consists of a low-pressure rotor and a high-pressure rotor, as shown in Figure1. The low-pressure rotor adopts the three-point supporting structures in which the mounting positions are at the 1#, 2#, and 5#. The high-pressure rotor absorbs the two-point supporting schemes in which the mounting positions are at the 3# and 4#. The low-pressure rotor and the high-pressure rotor are coupled through the intermediate ball bearing at 4#. The components of the low-pressure rotor system include the lumped mass of the disks m2 and m5 and the lumped mass of the shaft m1, m3, m4, and m6. The high-pressure rotor is regarded as a rigid rotor which only consists of the lumped mass m7. All the supporting structures, except the intermediate ball bearing 4#, are considered as elastic supporting systems with elastic supports and rolling bearings in series. They are simplified to spring-damping systems, wherein the stiffness and damping of the elastic supports are defined as k1, k2, k3, k5, c1, c2, c3, and c5, respectively. All the disks of the dual-rotor system have an unbalance magnitude and the eccentricity of the disks are described by ei (i = 2, 5, and 7) separately. In addition, the rub-impact fault occurs at the disk 2.Figure 1
Schematic diagram of the dual-rotor system.


## 2.2. Equation of Motion of Dual-Rotor System with Rub-Impact Fault in Maneuvering Flight

Lagrange’s equation is used to derive the equations of motion of the dual-rotor system. It is defined as [29](1)ddt∂L∂q˙j−∂L∂qj+∂D∂q˙j=Qj,j=1,2,…,k,where L = T-V, T and V are the kinetic energy and potential energy functions of the dual-rotor system, D is the dissipative energy function, qj is the generalized coordinate in the jth degree of freedom, and Qj is the generalized force corresponding to qj. The process of the derivation is as follows. First, the kinetic energy, potential energy and dissipative energy of the dual-rotor system are obtained. Then, the generalized forces of the ball bearing, rub-impact fault, gravity, and maneuver load are given. At last, the equations of the bending-torsional coupling motion are calculated by inserting the kinetic energy, potential energy, dissipative energy, and generalized forces into Lagrange’s equation given in equation (1).


### 2.2.1. Energy Functions of Dual-Rotor System

Based on the modeling method of Reference [25], the following assumption is made for the dual-rotor system. As we know, the low-pressure rotor is a relatively flexible slender shaft. The high-pressure rotor is short and the shaft diameter of which is large. Because the short and thick shaft is harder to display large torsional angle, the torsional vibrations of the high-pressure rotor are neglected. The bending and torsional motions are only considered at the disk mi (i = 2, 5) of the low-pressure rotor. They are described by the generalized coordinates xi,yi,αi (i = 2, 5), where xi and yi are the centroid displacements of the disk and αi is the torsional angle. Meanwhile, the torsional angle satisfies the following function φi=ωlowt+αi+ϕi (i = 2, 5), where φi is the phase angle of the mass centroid of the disk, ωlow is the rotating speed of the low-pressure rotor, and ϕi is the initial phase angle. Regarding the high-pressure rotor mi (i = 7) and the lump mass of the shaft of the low-pressure rotor mi (i = 1, 3, 4, 6), the bending motions are taken into account simply. They are defined as (xi, yi) (i = 1, 3, 4, 6, 7), where xi and yi are the centroid displacements of the high-pressure rotor and the shaft of the low-pressure rotor, respectively. In addition, the rotating speed of the high-pressure rotor is described by ωhigh, which satisfies ωhigh=γωωlow (γω is the rotating speed ratio of the high-pressure rotor to low-pressure rotor). Therefore, there are 16 degrees of freedom of the dual-rotor system which consist of 14 bending motions and 2 torsional motions. They are given by the generalized vector q:(2)q=x1,x2,…,x7,y1,y2,…,y7,α2,α5T,where the bending generalized vector x, y and the torsional generalized vector α are separately expressed as(3)x=x1x2x3x4x5x6x7,y=y1y2y3y4y5y6y7,α=α2α5.Because there are eccentricitiesei (i = 2, 5, 7) of the rotating disks, the mass centroid of the high-pressure rotor and the low-pressure rotor is described by the generalized vectors xc and yc:(4)xc=x1x2+e2cosφ2x3x4x5+e5cosφ5x6x7+e7cosφ7,yc=y1y2+e2sinφ2y3y4y5+e5sinφ5y6y7+e7sinφ7.The kinetic energyT of the dual-rotor system consists of translational kinetic energy Tt and rotational kinetic energy Tr. It is defined as(5)T=Tt+Tr=12x˙cTMcx˙cT+12y˙cTMcy˙cT+12φ˙Jcφ˙,where(6)x˙c=x˙1x˙2−e2φ˙2sinφ2x˙3x˙4x˙5−e5φ˙5sinφ5x˙6x˙7−e7φ˙7sinφ7,y˙c=y˙1y˙2+e2φ˙2cosφ2y˙3y˙4y˙5+e5φ˙5cosφ5y˙6y˙7+e7φ˙7cosφ7,φ˙=ωlow+α˙2ωlow+α˙5ωhigh,Mc=m10000m20000⋱0000m7,Jc=J2+m2e22000J5+m5e52000J7+m7e72.where x˙c and y˙c are the translational velocity vectors of the mass centroid, φ˙ is the rotating velocity vector of the disks, Mc is the mass matrix, and Jc is the moment of inertia relative to the mass centroid of the disks. Moreover, Ji = 1/2miRi2 (i = 2, 5, 7), where Ri are the radii of the disks.The potential energyUe of the dual-rotor system consists of the elastic potential energy of the shaft and the elastic potential energy of the elastic supports. It is defined as(7)Ue=12xTKxx+12yTKyy+12αTKtα,where(8)Kx=Ky=k1+k6−k600000−k6k6+k7−k700000−k7k2+k7+k8−k800000−k8k8+k9−k900000−k9k9+k10−k1000000−k10k5+k100000000k3,Kt=kt1+kt2−kt2−kt2kt2,where Kx and Ky are the bending stiffness matrices of the x- and y-directions, respectively, and Kt is the torsional stiffness matrix. In addition, the bending stiffness and torsional stiffness of the shaft are given by(9)k6=3EIl13,k7=12EIl23,k8=12EIl33,k9=12EIl43,k10=3EIl53,kt1=GIPl1,kt2=GIPl2+l3+l4,where k1, k2, k3, and k5 are the stiffness of the elastic supports, the values of which are obtained from the elastic stiffness of a squirrel cage, k6, k7, k8, k9, and k10 are the bending stiffness corresponding to the shaft segment l1, l2, l3, l4, and l5, respectively, kt1 and kt2 are the torsional stiffness corresponding to the shaft segment l1 and l2 + l3 + l4, separately, E and G are Young’s modulus and the shear modulus of the shaft, and I and IP are the moment of inertia and the polar moment of inertia of the shaft.The dissipative energy of the dual-rotor system can be described by the Rayleigh energy dissipation function:(10)D=12x˙TCxx˙+12y˙TCyy˙+12α˙TCtα˙,where(11)Cx=Cy=c1+c6−c600000−c6c6+c7−c700000−c7c2+c7+c8−c800000−c8c8+c9−c900000−c9c9+c10−c1000000−c10c5+c100000000c3,Ct=ct1+ct2−ct2−ct2ct2,where Cx, Cy, and Ct are the damping matrices corresponding to the bending and torsional motions, respectively. c1, c2, c3, and c5 are the given damping coefficients of the elastic supports, c6, c7, c8, c9, and c10 are the given damping coefficients corresponding to the bending motion of the shaft segments l1, l2, l3, l4, and l5, separately, and ct1 and ct2 are the damping coefficients corresponding to the torsional motions of the shaft segment l1 and l2+ l3+ l4, respectively. In this study, the damping coefficients ci (i = 1, 2, …, 10) are directly given, based on the Rayleigh damping ratio of the rotor system valuing around 0.1. Then, the initial free vibration of the simulation model can be reduced quickly.


### 2.2.2. Generalized Force of Intermediate Ball Bearing

The model of the intermediate ball bearing is shown in Figure2. The inner race of the ball bearing is directly fixed on the low-pressure rotor and the outer race is directly joined together with the high-pressure rotor. We assume that the balls are arranged equidistantly at the cage of the bearing. Then, the rotating speed of the cage of the bearing can be defined as(12)ωc=ωlowRinner+ωhighRouterRinner+Router,where Rinner and Router are the radii of the inner race and outer race of the ball bearing, respectively, and ωlow and ωhigh are the rotating speeds of the low-pressure rotor and high-pressure rotor, separately.Figure 2
Intermediate ball bearing model.The number of the balls is described byNb. Then the location of the jth ball is given by(13)θj=ωct+2πj−1Nb,j=1,2,…,Nb,The relative displacement between the inner and outer races of the ball bearing is expressed by (x7-x4, y7-y4). The initial bearing clearance is defined as δ0. Then, the contact deformation of the jth ball is expressed as(14)δj=x7−x4cosθj+y7−y4sinθj−δ0.According to the Hertzian contact theory [30], the force of the intermediate ball bearing is given by(15)FbxFby=Cb∑j=1Nbδj3/2Hδjcosθjsinθj,where Cb is the Hertz contact stiffness and H( ) is the Heaviside function.


### 2.2.3. Force of Rub-Impact Fault

The schematic diagram of the rub-impact fault is shown in Figure3. The relative displacement between the disk and the stator is described by r2=x22+y221/2. The initial clearance between the disk and the stator is defined as δ.Figure 3
Rub-impact fault model.When the rub-impact fault happens, there are normal forceFn and tangential force Fτ applied to the rotor system. According to [31], the rub-impact force can be expressed by(16)FrxFry=−krr2−δr21−μμ1x2y2,r2>δ,where kr and μ are the contact stiffness and the friction coefficient between the disk and stator, respectively. In addition, the friction applied to the disk leads to friction torque which is given by(17)Mr=μkrR2r2−δ,r2>δ,where R2 is the radius of the disk. Moreover, from equation (17), we get that that the friction torque is the function consisting of the bending displacement r2. Meanwhile, we know that the friction torque is the main factor resulting in the torsional motion. Therefore, the friction torque is one of the major factors leading to the bending-torsional coupling motion.


### 2.2.4. Maneuver Load and Gravity

As shown in Figure4, the influence of the barrel roll flight of an aircraft on the rub-impact dual-rotor system is investigated. The process of the transient barrel roll flight is simplified as three stages: at first (t∈0s,2s), the aircraft keeps the steady-state flight, and it enters into the barrel roll flight at t = 2 s and quits the barrel roll flight at t = 4 s. Then, it returns to the steady flight again (t∈4s,6s).Figure 4
Barrel roll flight of an aircraft.The actual barrel roll flight is very complicated. We have to make some simplified assumptions to study its influence. As shown in Figure5, the barrel roll flight is considered as a circular motion. The rolling speed and rolling radius of the aircraft are defined as ωbz and Rbz, respectively. In addition, because the rotation of the aircraft is faster than the revolution during the barrel roll flight, the rotation speed of the aircraft is defined as 2ωbz. To express the load state of the rotor system conveniently, two Cartesian coordinated references are used. The inertial reference frame O-XY, fixed to the ground is used for describing the flight status of the aircraft. The generalized reference frame o-xy, fixed to the aircraft, is used to describe the vibration of the dual-rotor system. The dual-rotor system is subjected to gravity mig and additional inertial force miωbz2Rbz due to the barrel roll flight of the aircraft. The gravity is always vertically downward. The additional inertial force is directed to the center of the rolling circle. Because the rolling direction of the aircraft is the same as the rotating direction of the dual-rotor system, there is no additional inertia moment applied to rotor system. Then, the generalized forces of the maneuver load and gravity are obtained by decomposing them into the generalized reference frame o-xy:(18)FgxiFgyi=miωbz2Rbzcosθbz−migsin2θbz−miωbz2Rbzsinθbz−migcos2θbz,where θ˙bz=ωbzt∈2s,4s.Figure 5
Force analysis of the rotor system during the barrel roll flight of an aircraft.


### 2.2.5. Equation of Bending-Torsional Coupling Motion of Dual-Rotor System

The equations of the bending-torsional coupling motion of the rub-impact dual-rotor system can be calculated by substituting the energy functions and generalized forces into equation (1):(19)Mq¨+Cq˙+Kq=Fr+Fb+Fu+Fg,where(20)M=Mc000Mc000Jc′,C=Cx000Cy000Ct,K=Kx000Ky000Kt,Jc′=J2+m2e2200J5+m5e52,Fr=0Frx000000Fry00000Mr0T,Fb=000Fbx00−Fbx000Fby00−Fby00T,Fu=0Fux200Fux50Fux70Fuy200Fuy50Fuy7Fuα1Fuα2T,Fg=Fgx1Fgx2Fgx3Fgx4Fgx5Fgx6Fgx7Fgy1Fgy2Fgy3Fgy4Fgy5Fgy6Fgy700T,Fux2=m2e2φ˙22cosφ2+φ¨2sinφ2,Fux5=m5e5φ˙52cosφ5+φ¨5sinφ5,Fux7=m7e7ωh2cosωht,Fuy2=m2e2φ˙22sinφ2−φ¨2cosφ2,Fuy5=m5e5φ˙52sinφ5−φ¨5cosφ5,Fuy7=m7e7ωh2sinωht,Fuα1=m2e2x¨2sinφ2−y¨2cosφ2,Fuα2=m5e5x¨5sinφ5−y¨5cosφ5,where M, C, and K are the generalized mass matrix, damping matrix, and stiffness matrix, respectively, Fr and Fb are the generalized force of the rub-impact fault and intermediate ball bearing, separately, Fu is the unbalanced force, and Fg is the generalized resultant force of the maneuver load and gravity, φi (i = 2, 5) are the time-varying phase angles of the mass centroid of the disks, which have been defined in the Section 2.2.1. In addition, from the derivation of the Fu, we know that unbalanced force is another factor leading to the bending-torsional coupling motion.


## 2.2.1. Energy Functions of Dual-Rotor System

Based on the modeling method of Reference [25], the following assumption is made for the dual-rotor system. As we know, the low-pressure rotor is a relatively flexible slender shaft. The high-pressure rotor is short and the shaft diameter of which is large. Because the short and thick shaft is harder to display large torsional angle, the torsional vibrations of the high-pressure rotor are neglected. The bending and torsional motions are only considered at the disk mi (i = 2, 5) of the low-pressure rotor. They are described by the generalized coordinates xi,yi,αi (i = 2, 5), where xi and yi are the centroid displacements of the disk and αi is the torsional angle. Meanwhile, the torsional angle satisfies the following function φi=ωlowt+αi+ϕi (i = 2, 5), where φi is the phase angle of the mass centroid of the disk, ωlow is the rotating speed of the low-pressure rotor, and ϕi is the initial phase angle. Regarding the high-pressure rotor mi (i = 7) and the lump mass of the shaft of the low-pressure rotor mi (i = 1, 3, 4, 6), the bending motions are taken into account simply. They are defined as (xi, yi) (i = 1, 3, 4, 6, 7), where xi and yi are the centroid displacements of the high-pressure rotor and the shaft of the low-pressure rotor, respectively. In addition, the rotating speed of the high-pressure rotor is described by ωhigh, which satisfies ωhigh=γωωlow (γω is the rotating speed ratio of the high-pressure rotor to low-pressure rotor). Therefore, there are 16 degrees of freedom of the dual-rotor system which consist of 14 bending motions and 2 torsional motions. They are given by the generalized vector q:(2)q=x1,x2,…,x7,y1,y2,…,y7,α2,α5T,where the bending generalized vector x, y and the torsional generalized vector α are separately expressed as(3)x=x1x2x3x4x5x6x7,y=y1y2y3y4y5y6y7,α=α2α5.Because there are eccentricitiesei (i = 2, 5, 7) of the rotating disks, the mass centroid of the high-pressure rotor and the low-pressure rotor is described by the generalized vectors xc and yc:(4)xc=x1x2+e2cosφ2x3x4x5+e5cosφ5x6x7+e7cosφ7,yc=y1y2+e2sinφ2y3y4y5+e5sinφ5y6y7+e7sinφ7.The kinetic energyT of the dual-rotor system consists of translational kinetic energy Tt and rotational kinetic energy Tr. It is defined as(5)T=Tt+Tr=12x˙cTMcx˙cT+12y˙cTMcy˙cT+12φ˙Jcφ˙,where(6)x˙c=x˙1x˙2−e2φ˙2sinφ2x˙3x˙4x˙5−e5φ˙5sinφ5x˙6x˙7−e7φ˙7sinφ7,y˙c=y˙1y˙2+e2φ˙2cosφ2y˙3y˙4y˙5+e5φ˙5cosφ5y˙6y˙7+e7φ˙7cosφ7,φ˙=ωlow+α˙2ωlow+α˙5ωhigh,Mc=m10000m20000⋱0000m7,Jc=J2+m2e22000J5+m5e52000J7+m7e72.where x˙c and y˙c are the translational velocity vectors of the mass centroid, φ˙ is the rotating velocity vector of the disks, Mc is the mass matrix, and Jc is the moment of inertia relative to the mass centroid of the disks. Moreover, Ji = 1/2miRi2 (i = 2, 5, 7), where Ri are the radii of the disks.The potential energyUe of the dual-rotor system consists of the elastic potential energy of the shaft and the elastic potential energy of the elastic supports. It is defined as(7)Ue=12xTKxx+12yTKyy+12αTKtα,where(8)Kx=Ky=k1+k6−k600000−k6k6+k7−k700000−k7k2+k7+k8−k800000−k8k8+k9−k900000−k9k9+k10−k1000000−k10k5+k100000000k3,Kt=kt1+kt2−kt2−kt2kt2,where Kx and Ky are the bending stiffness matrices of the x- and y-directions, respectively, and Kt is the torsional stiffness matrix. In addition, the bending stiffness and torsional stiffness of the shaft are given by(9)k6=3EIl13,k7=12EIl23,k8=12EIl33,k9=12EIl43,k10=3EIl53,kt1=GIPl1,kt2=GIPl2+l3+l4,where k1, k2, k3, and k5 are the stiffness of the elastic supports, the values of which are obtained from the elastic stiffness of a squirrel cage, k6, k7, k8, k9, and k10 are the bending stiffness corresponding to the shaft segment l1, l2, l3, l4, and l5, respectively, kt1 and kt2 are the torsional stiffness corresponding to the shaft segment l1 and l2 + l3 + l4, separately, E and G are Young’s modulus and the shear modulus of the shaft, and I and IP are the moment of inertia and the polar moment of inertia of the shaft.The dissipative energy of the dual-rotor system can be described by the Rayleigh energy dissipation function:(10)D=12x˙TCxx˙+12y˙TCyy˙+12α˙TCtα˙,where(11)Cx=Cy=c1+c6−c600000−c6c6+c7−c700000−c7c2+c7+c8−c800000−c8c8+c9−c900000−c9c9+c10−c1000000−c10c5+c100000000c3,Ct=ct1+ct2−ct2−ct2ct2,where Cx, Cy, and Ct are the damping matrices corresponding to the bending and torsional motions, respectively. c1, c2, c3, and c5 are the given damping coefficients of the elastic supports, c6, c7, c8, c9, and c10 are the given damping coefficients corresponding to the bending motion of the shaft segments l1, l2, l3, l4, and l5, separately, and ct1 and ct2 are the damping coefficients corresponding to the torsional motions of the shaft segment l1 and l2+ l3+ l4, respectively. In this study, the damping coefficients ci (i = 1, 2, …, 10) are directly given, based on the Rayleigh damping ratio of the rotor system valuing around 0.1. Then, the initial free vibration of the simulation model can be reduced quickly.


## 2.2.2. Generalized Force of Intermediate Ball Bearing

The model of the intermediate ball bearing is shown in Figure2. The inner race of the ball bearing is directly fixed on the low-pressure rotor and the outer race is directly joined together with the high-pressure rotor. We assume that the balls are arranged equidistantly at the cage of the bearing. Then, the rotating speed of the cage of the bearing can be defined as(12)ωc=ωlowRinner+ωhighRouterRinner+Router,where Rinner and Router are the radii of the inner race and outer race of the ball bearing, respectively, and ωlow and ωhigh are the rotating speeds of the low-pressure rotor and high-pressure rotor, separately.Figure 2
Intermediate ball bearing model.The number of the balls is described byNb. Then the location of the jth ball is given by(13)θj=ωct+2πj−1Nb,j=1,2,…,Nb,The relative displacement between the inner and outer races of the ball bearing is expressed by (x7-x4, y7-y4). The initial bearing clearance is defined as δ0. Then, the contact deformation of the jth ball is expressed as(14)δj=x7−x4cosθj+y7−y4sinθj−δ0.According to the Hertzian contact theory [30], the force of the intermediate ball bearing is given by(15)FbxFby=Cb∑j=1Nbδj3/2Hδjcosθjsinθj,where Cb is the Hertz contact stiffness and H( ) is the Heaviside function.


## 2.2.3. Force of Rub-Impact Fault

The schematic diagram of the rub-impact fault is shown in Figure3. The relative displacement between the disk and the stator is described by r2=x22+y221/2. The initial clearance between the disk and the stator is defined as δ.Figure 3
Rub-impact fault model.When the rub-impact fault happens, there are normal forceFn and tangential force Fτ applied to the rotor system. According to [31], the rub-impact force can be expressed by(16)FrxFry=−krr2−δr21−μμ1x2y2,r2>δ,where kr and μ are the contact stiffness and the friction coefficient between the disk and stator, respectively. In addition, the friction applied to the disk leads to friction torque which is given by(17)Mr=μkrR2r2−δ,r2>δ,where R2 is the radius of the disk. Moreover, from equation (17), we get that that the friction torque is the function consisting of the bending displacement r2. Meanwhile, we know that the friction torque is the main factor resulting in the torsional motion. Therefore, the friction torque is one of the major factors leading to the bending-torsional coupling motion.


## 2.2.4. Maneuver Load and Gravity

As shown in Figure4, the influence of the barrel roll flight of an aircraft on the rub-impact dual-rotor system is investigated. The process of the transient barrel roll flight is simplified as three stages: at first (t∈0s,2s), the aircraft keeps the steady-state flight, and it enters into the barrel roll flight at t = 2 s and quits the barrel roll flight at t = 4 s. Then, it returns to the steady flight again (t∈4s,6s).Figure 4
Barrel roll flight of an aircraft.The actual barrel roll flight is very complicated. We have to make some simplified assumptions to study its influence. As shown in Figure5, the barrel roll flight is considered as a circular motion. The rolling speed and rolling radius of the aircraft are defined as ωbz and Rbz, respectively. In addition, because the rotation of the aircraft is faster than the revolution during the barrel roll flight, the rotation speed of the aircraft is defined as 2ωbz. To express the load state of the rotor system conveniently, two Cartesian coordinated references are used. The inertial reference frame O-XY, fixed to the ground is used for describing the flight status of the aircraft. The generalized reference frame o-xy, fixed to the aircraft, is used to describe the vibration of the dual-rotor system. The dual-rotor system is subjected to gravity mig and additional inertial force miωbz2Rbz due to the barrel roll flight of the aircraft. The gravity is always vertically downward. The additional inertial force is directed to the center of the rolling circle. Because the rolling direction of the aircraft is the same as the rotating direction of the dual-rotor system, there is no additional inertia moment applied to rotor system. Then, the generalized forces of the maneuver load and gravity are obtained by decomposing them into the generalized reference frame o-xy:(18)FgxiFgyi=miωbz2Rbzcosθbz−migsin2θbz−miωbz2Rbzsinθbz−migcos2θbz,where θ˙bz=ωbzt∈2s,4s.Figure 5
Force analysis of the rotor system during the barrel roll flight of an aircraft.


## 2.2.5. Equation of Bending-Torsional Coupling Motion of Dual-Rotor System

The equations of the bending-torsional coupling motion of the rub-impact dual-rotor system can be calculated by substituting the energy functions and generalized forces into equation (1):(19)Mq¨+Cq˙+Kq=Fr+Fb+Fu+Fg,where(20)M=Mc000Mc000Jc′,C=Cx000Cy000Ct,K=Kx000Ky000Kt,Jc′=J2+m2e2200J5+m5e52,Fr=0Frx000000Fry00000Mr0T,Fb=000Fbx00−Fbx000Fby00−Fby00T,Fu=0Fux200Fux50Fux70Fuy200Fuy50Fuy7Fuα1Fuα2T,Fg=Fgx1Fgx2Fgx3Fgx4Fgx5Fgx6Fgx7Fgy1Fgy2Fgy3Fgy4Fgy5Fgy6Fgy700T,Fux2=m2e2φ˙22cosφ2+φ¨2sinφ2,Fux5=m5e5φ˙52cosφ5+φ¨5sinφ5,Fux7=m7e7ωh2cosωht,Fuy2=m2e2φ˙22sinφ2−φ¨2cosφ2,Fuy5=m5e5φ˙52sinφ5−φ¨5cosφ5,Fuy7=m7e7ωh2sinωht,Fuα1=m2e2x¨2sinφ2−y¨2cosφ2,Fuα2=m5e5x¨5sinφ5−y¨5cosφ5,where M, C, and K are the generalized mass matrix, damping matrix, and stiffness matrix, respectively, Fr and Fb are the generalized force of the rub-impact fault and intermediate ball bearing, separately, Fu is the unbalanced force, and Fg is the generalized resultant force of the maneuver load and gravity, φi (i = 2, 5) are the time-varying phase angles of the mass centroid of the disks, which have been defined in the Section 2.2.1. In addition, from the derivation of the Fu, we know that unbalanced force is another factor leading to the bending-torsional coupling motion.


## 2.3. Simulation Parameter and Calculation Process

Because of the strong nonlinearity of the rotor system, the numerical integration method is adopted to solve the (19). The bending-torsional coupling dynamic responses of the rub-impact dual-rotor are calculated via the 4th Runge-Kutta method with varied step. According to Reference [25], the simulation parameters of the dual-rotor system are listed in Table 1.Table 1
Parameters of dual-rotor model.
Physical descriptionParameterValueMassm1, m2, m3, m4, m5, m6, m7 (kg)2, 10, 2, 2, 6, 2, 8Radius of the diskR2, R5 (m)0.075, 0.03Initial phase angle of the diskϕ2,ϕ5 (rad)0, 0Eccentricity of the diske2, e5, e7μm46.35, 46.35, 47.5Stiffness of the elastic supportk1, k2, k3, k5 (N/m)106, 106, 106, 106Length of the shaft segmentl1, l2, l3, l4, l5 (m)0.33, 0.33, 0.6, 0.23, 0.23Material properties of the shaftE, G (Pa)2 × 1011, 7.7 × 1010Physical properties of the shaftI, Ip (m4)1.178 × 10−7, 2.356 × 10−7Damping of the elastic supportc1, c2, c3, c5 (N.s/m)300, 300, 300, 300Damping of the bending motionc6, c7, c8, c9, c10 (Ns/m)150, 150, 150, 150, 150Damping of the torsional motionct1, ct2 (N.m.s/rad)1, 1Parameters of the rub-impact faultkr (N/m), δ(m), μ8 × 106, 1.34 × 10−4, 0.1Parameters of the ball bearingCb (N/m3/2), δ0,Rinner,Router(mm)1.334 × 1010, 0.005, 40.1, 63.9, 8Rotating speeds of the rotor systemωlow,ωhighrad/s,γω1500, 2700, 1.8Parameters of the maneuvering flightωbzrad/s, Rbz (m)π, 3In the process of the dynamic analysis, the relative torsional angleΔα=α2−α5 is chosen to show the torsional characteristics of the dual-rotor system. The bending displacements (x2 and y2) of the disk m2 are selected to reveal the bending characteristics. In addition, the motion of the dual-rotor system consists of some responses with low frequency and large amplitude because of the time-varying influence of the maneuvering flight. They are useless for the analysis of the rub-impact fault features. Therefore, a high-pass filter is applied to the response of the dual-rotor system. The detailed calculating process is shown in Figure 6.Figure 6
Calculating process.


## 3. Discussion of Numerical Results

### 3.1. Transient Responses without Rub-Impact Fault

The transient responses of the dual-rotor system are studied firstly without considering the rub -impact fault. The generalized resultant force of the maneuver load and gravity applied to the diskm2 is chosen as the example for illustrating the effect of the barrel roll flight. As shown in Figure 7, the rotor system suffers gravity and additional inertial maneuver load during the barrel roll flight (t∈2s,4s). Because the angle between the maneuver load and gravity varies, the generalized resultant force of the maneuver load and gravity is not constant but time-varying. The transient axis orbit of the disk m2 corresponding to this process is given in Figure 8. It can be seen that at first t∈0,2s the rotor system keeps steady response with small amplitude. As the aircraft enters into the barrel roll flight t∈2s,4s, the axis orbit of the rotor system is deviated entirely. The offset of the axis orbit changes over time due to the varying additional generalized force. Therefore, the rub-impact fault of the rotor system can be induced by the barrel roll flight and the levels of the rub-impact are time-varying during the flight. After the barrel roll flight, the rotor system returns to the steady response again.Figure 7
Generalized resultant force of the maneuver load and gravity applied to the disk m2.Figure 8
Transient axis orbit of the disk m2 during the flight of the aircraft.Without considering the rub-impact fault, the transient time-frequency waterfall curves of the bending and torsional motions of the rotor system are shown in Figure9. Form Figure 9(a), we know that the bending responses of the rotor system are mainly based on the rotating frequencies ωlow and ωhigh. In addition, at the start and end of the barrel roll flight (t = 2 s and t = 4 s), the natural bending frequency of the rotor system ωb is aroused because of the sudden disturbance of the flight. As shown in Figure 9(b), the torsional responses of the rotor system display combined rotating frequency ωhigh−ωlow which is caused by the unbalanced force of the rotor system. Moreover, we know that the barrel roll flight has little effect on the time-frequency responses of the torsional motion of the rotor system.Figure 9
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system without the rub-impact fault. (a) Bending motion and (b) torsional motion.
(a)(b)


### 3.2. Transient Responses with Rub-Impact Fault

#### 3.2.1. Effect of Rub-Impact Stiffness

Taking the rub-impact stiffness as the control parameter, the transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 2 × 106, 8 × 106 and 1.4 × 107 N/m are shown in Figures 10–12, respectively. The other parameters are listed in Table 1. A variety of comparisons for the rub-impact fault features of the bending and torsional motions are observed.Figure 10
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 2 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 11
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 8 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 12
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 1.4 × 107 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)From Figure10(a), we can see that for kr = 2 × 106 N/m, the bending motion of the rotor system shows no rub-impact feature. There are only the rotating frequencies ωlow and ωhigh and the natural bending frequency ωb in the waterfall curves which are the same as the responses in Figure 9(a). However, the rub-impact fault does happen. As shown in Figure 10(b), the torsional motion of the rotor system displays apparent rub-impact behavior during the barrel roll flight. Because of the transmission of the bending-torsional coupling motion, the torsional responses show the rotating fault frequencies ωlow and ωhigh during the process of the rub-impact. It demonstrates that at the low rub-impact stiffness, the torsional responses work better than bending responses in rub-impact fault detection.When the rub-impact stiffness increases fromkr = 2 × 106 N/m to kr = 8 × 106 N/m, by comparing Figure 11(a) with Figure 10(a), we can see that bending motion of the rotor system shows the fault frequencies. It consists of the fractional fault frequency 1/2ωlow and the natural bending frequency ωb during the barrel roll flight. The amplitude of the fault frequencies is changed over time due to the variable generalized forces. As shown in Figure 11(b), the torsional motion of the rotor system behaves more complicated fault frequencies. It contains the fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. However, the natural bending frequency ωb does not transmit to the torsional motion. Therefore, at the middle segment of the rub-impact stiffness, the bending responses and torsional responses both can be used for the rub-impact detection.When the rub-impact stiffness continues increasing tokr = 1.4 × 107 N/m, the transient time-frequency waterfall curves of the bending and torsional motions are shown in Figures 12(a) and 12(b), respectively. As shown in Figure 12(a), the bending responses of the rotor system show the fractional fault frequency 1/2ωlow. But, the free vibration component ωb does not arise because the rub-impact stiffness is strengthened. The torsional responses in Figure 12(b) still contain complicated fractional fault frequencies 1/2ωlow and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the high rub-impact stiffness, the bending responses and torsional responses still can be used for the rub-impact detection.The spectrum cascades of the bending and torsional motions with varied rub-impact stiffness are shown in Figure13. For every certain rub-impact stiffness, we chose the transient responses of the rotor system during the maneuver flight t∈2.2s,3.8s to obtain the spectrum. According to the rub-impact stiffness range, the spectrum cascades are divided into three bifurcation regions, in which kr = 5.5 × 106 N/m and 1.15 × 107 N/m are the bifurcation points. When kr < 5.5 × 106 N/m, the bending responses show no rub-impact features, but the torsional responses displays the early rotating fault frequencies. With kr > 5.5 × 106 N/m, the bending responses begin to show the fractional fault frequency and natural bending frequency. In addition, the natural bending frequency of the bending motion disappears at kr > 1.15 × 107 N/m because of the further enhancement of the rub-impact stiffness. Comparing with the bending responses, the torsional responses contain more complicated fault frequencies within the whole region of the rub-impact stiffness. Moreover, the amplitude of the rotating fault frequency ωlow increases linearly with the enhancement of the rub-impact stiffness which can be used for qualitative detection of the rub-impact stiffness. Therefore, the characteristic frequencies of the torsional motion will be more advantageous for the state detection and rub-impact fault diagnosis.Figure 13
Spectrum cascades of the bending and torsional motions with varied rub-impact stiffnesskr∈5×105N/m,2×107N/m. (a) Bending motion and (b) torsional motion.
(a)(b)


#### 3.2.2. Effect of Rub-Impact Friction Coefficient

Taking the rub-impact friction coefficient as the control parameter, the spectrum cascades of the bending and torsional motions of the rotor system forμ∈0,0.3 are shown in Figure 14(a) and 14(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rub-impact friction coefficient.Figure 14
Spectrum cascades of the bending and torsional motions with varied rub-impact friction coefficientsμ∈0,0.3. (a) Bending motion and (b) torsional motion.
(a)(b)As shown in Figure14(a), the friction coefficient has little effect on the bending responses. The bending responses show the fractional fault frequency 1/2ωlow and the natural bending frequency ωb within the whole region of the friction coefficient. It proves that rub-impact fault characteristics of the bending motion are mainly caused by the rub-impact stiffness but not the friction coefficient. However, the rub-impact friction coefficient affects the torsional motion greatly. From Figure 14(b), we can see that the amplitude of the fault frequencies 1/2ωlow, 3/2ωlow, ωlow, and ωhigh increases with the increase of the friction coefficient. In addition, the rotating fault frequency ωlow arises even at the little friction coefficient. It demonstrates that the torsional motion of the rotor system is sensitive to the rub-impact friction coefficient. The larger the friction coefficient is, the more obvious the fault characteristics of the torsional motion are. Therefore, the friction coefficient plays a decisive role on the transmission of the bending-torsional coupling motion.


#### 3.2.3. Effect of Rotating Speed

To show the rub-impact fault characteristics of the bending and torsional motions with varied rotating speeds, the spectrum cascades forωlow∈500rad/s,3000rad/s are shown in Figure 15(a) and 15(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rotating speed.Figure 15
Spectrum cascades of the bending and torsional motions with varied rotating speedsωlow∈500rad/s,3000rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)According to the frequency features of the bending responses, the rotating speeds are divided into three regions: low-speed region, intermediate-speed region, and high-speed region. The bending responses display fractional rub-impact fault frequency1/2ωlow and natural bending frequency ωb within the intermediate speed region ωlow∈1140rad/s,1740rad/s. However, at the low-speed region and high-speed region, the bending responses show no fault frequency. Comparing with the bending responses, the torsional responses display rub-impact fault characteristics within the whole region of the rotating speed. In addition, there are more complicated fractional fault frequencies at the intermediate speed region.The transient time-frequency waterfall curves of the bending and torsional motions forωlow=1500 rad/s and 1800 rad/s are chosen as the examples for illustrating the effect of the rotating speed.The responses forωlow = 1500 rad/s have been obtained, and the results are shown in Figure 11. It can be seen that at the intermediate speed region, the bending and torsional responses both show the rub-impact fault frequencies. The bending responses display fractional fault frequency 1/2ωlow and the natural bending frequency ωb. The torsional responses behave more complicated fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the intermediate speed region, the bending and torsional responses both can be used for the rub-impact fault detection.However, this conclusion is inapplicable within the low-speed region and high-speed region. The bending and the torsional responses forωlow = 1800 rad/s are shown in Figures 16(a) and 16(b) for explanation. As shown in Figure 16(a), the bending responses show no fault frequency during the maneuvering flight. Hence, at ωlow = 1800 rad/s, the spectrum of the bending responses cannot be used to distinguish the occurrence of the rub-impact fault. From Figure 16(b), we can see that the apparent rotating fault frequency ωlow is transmitted to the torsional motion because of the influence of the rub-impact. Therefore, within the whole region of the rotating speed, the torsional responses will be more effective than the bending responses for the rub-impact fault detection.Figure 16
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forωlow = 1800 rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)


## 3.1. Transient Responses without Rub-Impact Fault

The transient responses of the dual-rotor system are studied firstly without considering the rub -impact fault. The generalized resultant force of the maneuver load and gravity applied to the diskm2 is chosen as the example for illustrating the effect of the barrel roll flight. As shown in Figure 7, the rotor system suffers gravity and additional inertial maneuver load during the barrel roll flight (t∈2s,4s). Because the angle between the maneuver load and gravity varies, the generalized resultant force of the maneuver load and gravity is not constant but time-varying. The transient axis orbit of the disk m2 corresponding to this process is given in Figure 8. It can be seen that at first t∈0,2s the rotor system keeps steady response with small amplitude. As the aircraft enters into the barrel roll flight t∈2s,4s, the axis orbit of the rotor system is deviated entirely. The offset of the axis orbit changes over time due to the varying additional generalized force. Therefore, the rub-impact fault of the rotor system can be induced by the barrel roll flight and the levels of the rub-impact are time-varying during the flight. After the barrel roll flight, the rotor system returns to the steady response again.Figure 7
Generalized resultant force of the maneuver load and gravity applied to the disk m2.Figure 8
Transient axis orbit of the disk m2 during the flight of the aircraft.Without considering the rub-impact fault, the transient time-frequency waterfall curves of the bending and torsional motions of the rotor system are shown in Figure9. Form Figure 9(a), we know that the bending responses of the rotor system are mainly based on the rotating frequencies ωlow and ωhigh. In addition, at the start and end of the barrel roll flight (t = 2 s and t = 4 s), the natural bending frequency of the rotor system ωb is aroused because of the sudden disturbance of the flight. As shown in Figure 9(b), the torsional responses of the rotor system display combined rotating frequency ωhigh−ωlow which is caused by the unbalanced force of the rotor system. Moreover, we know that the barrel roll flight has little effect on the time-frequency responses of the torsional motion of the rotor system.Figure 9
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system without the rub-impact fault. (a) Bending motion and (b) torsional motion.
(a)(b)


## 3.2. Transient Responses with Rub-Impact Fault

### 3.2.1. Effect of Rub-Impact Stiffness

Taking the rub-impact stiffness as the control parameter, the transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 2 × 106, 8 × 106 and 1.4 × 107 N/m are shown in Figures 10–12, respectively. The other parameters are listed in Table 1. A variety of comparisons for the rub-impact fault features of the bending and torsional motions are observed.Figure 10
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 2 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 11
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 8 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 12
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 1.4 × 107 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)From Figure10(a), we can see that for kr = 2 × 106 N/m, the bending motion of the rotor system shows no rub-impact feature. There are only the rotating frequencies ωlow and ωhigh and the natural bending frequency ωb in the waterfall curves which are the same as the responses in Figure 9(a). However, the rub-impact fault does happen. As shown in Figure 10(b), the torsional motion of the rotor system displays apparent rub-impact behavior during the barrel roll flight. Because of the transmission of the bending-torsional coupling motion, the torsional responses show the rotating fault frequencies ωlow and ωhigh during the process of the rub-impact. It demonstrates that at the low rub-impact stiffness, the torsional responses work better than bending responses in rub-impact fault detection.When the rub-impact stiffness increases fromkr = 2 × 106 N/m to kr = 8 × 106 N/m, by comparing Figure 11(a) with Figure 10(a), we can see that bending motion of the rotor system shows the fault frequencies. It consists of the fractional fault frequency 1/2ωlow and the natural bending frequency ωb during the barrel roll flight. The amplitude of the fault frequencies is changed over time due to the variable generalized forces. As shown in Figure 11(b), the torsional motion of the rotor system behaves more complicated fault frequencies. It contains the fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. However, the natural bending frequency ωb does not transmit to the torsional motion. Therefore, at the middle segment of the rub-impact stiffness, the bending responses and torsional responses both can be used for the rub-impact detection.When the rub-impact stiffness continues increasing tokr = 1.4 × 107 N/m, the transient time-frequency waterfall curves of the bending and torsional motions are shown in Figures 12(a) and 12(b), respectively. As shown in Figure 12(a), the bending responses of the rotor system show the fractional fault frequency 1/2ωlow. But, the free vibration component ωb does not arise because the rub-impact stiffness is strengthened. The torsional responses in Figure 12(b) still contain complicated fractional fault frequencies 1/2ωlow and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the high rub-impact stiffness, the bending responses and torsional responses still can be used for the rub-impact detection.The spectrum cascades of the bending and torsional motions with varied rub-impact stiffness are shown in Figure13. For every certain rub-impact stiffness, we chose the transient responses of the rotor system during the maneuver flight t∈2.2s,3.8s to obtain the spectrum. According to the rub-impact stiffness range, the spectrum cascades are divided into three bifurcation regions, in which kr = 5.5 × 106 N/m and 1.15 × 107 N/m are the bifurcation points. When kr < 5.5 × 106 N/m, the bending responses show no rub-impact features, but the torsional responses displays the early rotating fault frequencies. With kr > 5.5 × 106 N/m, the bending responses begin to show the fractional fault frequency and natural bending frequency. In addition, the natural bending frequency of the bending motion disappears at kr > 1.15 × 107 N/m because of the further enhancement of the rub-impact stiffness. Comparing with the bending responses, the torsional responses contain more complicated fault frequencies within the whole region of the rub-impact stiffness. Moreover, the amplitude of the rotating fault frequency ωlow increases linearly with the enhancement of the rub-impact stiffness which can be used for qualitative detection of the rub-impact stiffness. Therefore, the characteristic frequencies of the torsional motion will be more advantageous for the state detection and rub-impact fault diagnosis.Figure 13
Spectrum cascades of the bending and torsional motions with varied rub-impact stiffnesskr∈5×105N/m,2×107N/m. (a) Bending motion and (b) torsional motion.
(a)(b)


### 3.2.2. Effect of Rub-Impact Friction Coefficient

Taking the rub-impact friction coefficient as the control parameter, the spectrum cascades of the bending and torsional motions of the rotor system forμ∈0,0.3 are shown in Figure 14(a) and 14(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rub-impact friction coefficient.Figure 14
Spectrum cascades of the bending and torsional motions with varied rub-impact friction coefficientsμ∈0,0.3. (a) Bending motion and (b) torsional motion.
(a)(b)As shown in Figure14(a), the friction coefficient has little effect on the bending responses. The bending responses show the fractional fault frequency 1/2ωlow and the natural bending frequency ωb within the whole region of the friction coefficient. It proves that rub-impact fault characteristics of the bending motion are mainly caused by the rub-impact stiffness but not the friction coefficient. However, the rub-impact friction coefficient affects the torsional motion greatly. From Figure 14(b), we can see that the amplitude of the fault frequencies 1/2ωlow, 3/2ωlow, ωlow, and ωhigh increases with the increase of the friction coefficient. In addition, the rotating fault frequency ωlow arises even at the little friction coefficient. It demonstrates that the torsional motion of the rotor system is sensitive to the rub-impact friction coefficient. The larger the friction coefficient is, the more obvious the fault characteristics of the torsional motion are. Therefore, the friction coefficient plays a decisive role on the transmission of the bending-torsional coupling motion.


### 3.2.3. Effect of Rotating Speed

To show the rub-impact fault characteristics of the bending and torsional motions with varied rotating speeds, the spectrum cascades forωlow∈500rad/s,3000rad/s are shown in Figure 15(a) and 15(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rotating speed.Figure 15
Spectrum cascades of the bending and torsional motions with varied rotating speedsωlow∈500rad/s,3000rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)According to the frequency features of the bending responses, the rotating speeds are divided into three regions: low-speed region, intermediate-speed region, and high-speed region. The bending responses display fractional rub-impact fault frequency1/2ωlow and natural bending frequency ωb within the intermediate speed region ωlow∈1140rad/s,1740rad/s. However, at the low-speed region and high-speed region, the bending responses show no fault frequency. Comparing with the bending responses, the torsional responses display rub-impact fault characteristics within the whole region of the rotating speed. In addition, there are more complicated fractional fault frequencies at the intermediate speed region.The transient time-frequency waterfall curves of the bending and torsional motions forωlow=1500 rad/s and 1800 rad/s are chosen as the examples for illustrating the effect of the rotating speed.The responses forωlow = 1500 rad/s have been obtained, and the results are shown in Figure 11. It can be seen that at the intermediate speed region, the bending and torsional responses both show the rub-impact fault frequencies. The bending responses display fractional fault frequency 1/2ωlow and the natural bending frequency ωb. The torsional responses behave more complicated fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the intermediate speed region, the bending and torsional responses both can be used for the rub-impact fault detection.However, this conclusion is inapplicable within the low-speed region and high-speed region. The bending and the torsional responses forωlow = 1800 rad/s are shown in Figures 16(a) and 16(b) for explanation. As shown in Figure 16(a), the bending responses show no fault frequency during the maneuvering flight. Hence, at ωlow = 1800 rad/s, the spectrum of the bending responses cannot be used to distinguish the occurrence of the rub-impact fault. From Figure 16(b), we can see that the apparent rotating fault frequency ωlow is transmitted to the torsional motion because of the influence of the rub-impact. Therefore, within the whole region of the rotating speed, the torsional responses will be more effective than the bending responses for the rub-impact fault detection.Figure 16
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forωlow = 1800 rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)


## 3.2.1. Effect of Rub-Impact Stiffness

Taking the rub-impact stiffness as the control parameter, the transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 2 × 106, 8 × 106 and 1.4 × 107 N/m are shown in Figures 10–12, respectively. The other parameters are listed in Table 1. A variety of comparisons for the rub-impact fault features of the bending and torsional motions are observed.Figure 10
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 2 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 11
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system for kr = 8 × 106 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)Figure 12
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forkr = 1.4 × 107 N/m. (a) Bending motion and (b) torsional motion.
(a)(b)From Figure10(a), we can see that for kr = 2 × 106 N/m, the bending motion of the rotor system shows no rub-impact feature. There are only the rotating frequencies ωlow and ωhigh and the natural bending frequency ωb in the waterfall curves which are the same as the responses in Figure 9(a). However, the rub-impact fault does happen. As shown in Figure 10(b), the torsional motion of the rotor system displays apparent rub-impact behavior during the barrel roll flight. Because of the transmission of the bending-torsional coupling motion, the torsional responses show the rotating fault frequencies ωlow and ωhigh during the process of the rub-impact. It demonstrates that at the low rub-impact stiffness, the torsional responses work better than bending responses in rub-impact fault detection.When the rub-impact stiffness increases fromkr = 2 × 106 N/m to kr = 8 × 106 N/m, by comparing Figure 11(a) with Figure 10(a), we can see that bending motion of the rotor system shows the fault frequencies. It consists of the fractional fault frequency 1/2ωlow and the natural bending frequency ωb during the barrel roll flight. The amplitude of the fault frequencies is changed over time due to the variable generalized forces. As shown in Figure 11(b), the torsional motion of the rotor system behaves more complicated fault frequencies. It contains the fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. However, the natural bending frequency ωb does not transmit to the torsional motion. Therefore, at the middle segment of the rub-impact stiffness, the bending responses and torsional responses both can be used for the rub-impact detection.When the rub-impact stiffness continues increasing tokr = 1.4 × 107 N/m, the transient time-frequency waterfall curves of the bending and torsional motions are shown in Figures 12(a) and 12(b), respectively. As shown in Figure 12(a), the bending responses of the rotor system show the fractional fault frequency 1/2ωlow. But, the free vibration component ωb does not arise because the rub-impact stiffness is strengthened. The torsional responses in Figure 12(b) still contain complicated fractional fault frequencies 1/2ωlow and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the high rub-impact stiffness, the bending responses and torsional responses still can be used for the rub-impact detection.The spectrum cascades of the bending and torsional motions with varied rub-impact stiffness are shown in Figure13. For every certain rub-impact stiffness, we chose the transient responses of the rotor system during the maneuver flight t∈2.2s,3.8s to obtain the spectrum. According to the rub-impact stiffness range, the spectrum cascades are divided into three bifurcation regions, in which kr = 5.5 × 106 N/m and 1.15 × 107 N/m are the bifurcation points. When kr < 5.5 × 106 N/m, the bending responses show no rub-impact features, but the torsional responses displays the early rotating fault frequencies. With kr > 5.5 × 106 N/m, the bending responses begin to show the fractional fault frequency and natural bending frequency. In addition, the natural bending frequency of the bending motion disappears at kr > 1.15 × 107 N/m because of the further enhancement of the rub-impact stiffness. Comparing with the bending responses, the torsional responses contain more complicated fault frequencies within the whole region of the rub-impact stiffness. Moreover, the amplitude of the rotating fault frequency ωlow increases linearly with the enhancement of the rub-impact stiffness which can be used for qualitative detection of the rub-impact stiffness. Therefore, the characteristic frequencies of the torsional motion will be more advantageous for the state detection and rub-impact fault diagnosis.Figure 13
Spectrum cascades of the bending and torsional motions with varied rub-impact stiffnesskr∈5×105N/m,2×107N/m. (a) Bending motion and (b) torsional motion.
(a)(b)


## 3.2.2. Effect of Rub-Impact Friction Coefficient

Taking the rub-impact friction coefficient as the control parameter, the spectrum cascades of the bending and torsional motions of the rotor system forμ∈0,0.3 are shown in Figure 14(a) and 14(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rub-impact friction coefficient.Figure 14
Spectrum cascades of the bending and torsional motions with varied rub-impact friction coefficientsμ∈0,0.3. (a) Bending motion and (b) torsional motion.
(a)(b)As shown in Figure14(a), the friction coefficient has little effect on the bending responses. The bending responses show the fractional fault frequency 1/2ωlow and the natural bending frequency ωb within the whole region of the friction coefficient. It proves that rub-impact fault characteristics of the bending motion are mainly caused by the rub-impact stiffness but not the friction coefficient. However, the rub-impact friction coefficient affects the torsional motion greatly. From Figure 14(b), we can see that the amplitude of the fault frequencies 1/2ωlow, 3/2ωlow, ωlow, and ωhigh increases with the increase of the friction coefficient. In addition, the rotating fault frequency ωlow arises even at the little friction coefficient. It demonstrates that the torsional motion of the rotor system is sensitive to the rub-impact friction coefficient. The larger the friction coefficient is, the more obvious the fault characteristics of the torsional motion are. Therefore, the friction coefficient plays a decisive role on the transmission of the bending-torsional coupling motion.


## 3.2.3. Effect of Rotating Speed

To show the rub-impact fault characteristics of the bending and torsional motions with varied rotating speeds, the spectrum cascades forωlow∈500rad/s,3000rad/s are shown in Figure 15(a) and 15(b), respectively. The transient rub-impact responses t∈2.2s,3.8s are chosen to obtain the spectrum for the certain rotating speed.Figure 15
Spectrum cascades of the bending and torsional motions with varied rotating speedsωlow∈500rad/s,3000rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)According to the frequency features of the bending responses, the rotating speeds are divided into three regions: low-speed region, intermediate-speed region, and high-speed region. The bending responses display fractional rub-impact fault frequency1/2ωlow and natural bending frequency ωb within the intermediate speed region ωlow∈1140rad/s,1740rad/s. However, at the low-speed region and high-speed region, the bending responses show no fault frequency. Comparing with the bending responses, the torsional responses display rub-impact fault characteristics within the whole region of the rotating speed. In addition, there are more complicated fractional fault frequencies at the intermediate speed region.The transient time-frequency waterfall curves of the bending and torsional motions forωlow=1500 rad/s and 1800 rad/s are chosen as the examples for illustrating the effect of the rotating speed.The responses forωlow = 1500 rad/s have been obtained, and the results are shown in Figure 11. It can be seen that at the intermediate speed region, the bending and torsional responses both show the rub-impact fault frequencies. The bending responses display fractional fault frequency 1/2ωlow and the natural bending frequency ωb. The torsional responses behave more complicated fractional fault frequencies 1/2ωlow, 5/4ωlow, and 3/2ωlow and the rotating fault frequencies ωlow and ωhigh. Therefore, at the intermediate speed region, the bending and torsional responses both can be used for the rub-impact fault detection.However, this conclusion is inapplicable within the low-speed region and high-speed region. The bending and the torsional responses forωlow = 1800 rad/s are shown in Figures 16(a) and 16(b) for explanation. As shown in Figure 16(a), the bending responses show no fault frequency during the maneuvering flight. Hence, at ωlow = 1800 rad/s, the spectrum of the bending responses cannot be used to distinguish the occurrence of the rub-impact fault. From Figure 16(b), we can see that the apparent rotating fault frequency ωlow is transmitted to the torsional motion because of the influence of the rub-impact. Therefore, within the whole region of the rotating speed, the torsional responses will be more effective than the bending responses for the rub-impact fault detection.Figure 16
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor system forωlow = 1800 rad/s. (a) Bending motion and (b) torsional motion.
(a)(b)


## 4. Discussion of Experimental Results

### 4.1. Experimental Device

In this section, an experiment is carried out to verify the results obtained by the numerical simulation. Figure17(a) shows the main experimental devices, which consist of a rotor test rig, a rub-impact device, a base motion platform, two rotary encoders, two displacement sensors, and a set of data acquisition system. Figure 17(b) displays the rub-impact device. The stator of the rub-impact device is a circular truncated cone, whose inner race is slightly larger than that of the disc. The material of the stator is steel. The rub-impact device can be moved from right to left, so the clearance of the stator and disc can be adjusted from 0 to 5 mm. The rub-impact process can be well imitated because of the convenient clearance adjustment.Figure 17
Main experimental devices and rub-impact device. (a) Main experimental devices and (b) rub-impact device.
(a)(b)Figure18 shows the schematic diagram of the main experimental device. The rub-impact rotor system installed on the base motion platform is used for imitating the aeroengine mounted on the aircraft. A single-rotor test rig is adopted to verify the simulation results qualitatively. There are two sets of control system: rotor motor control system and base motion platform control system. Then, with the rotor at the constant speed, the whole rotor system can also rotate around the rolling shaft of the base motion platform, so as to simulate the vibration of the aero-engine during the barrel roll flight. The rolling radius and maximum rolling speed of the base motion platform are 0.4 m and π rad/s, respectively. In addition, to carry out the signal acquisition and the electric power supply of the motor, a slip ring is added to the base motion platform system, through which the above function can be achieved on the premise of the barrel roll flight of the base motion platform.Figure 18
Schematic diagram of the experimental device. 1, rotary encoder, 2, displacement sensors, 3, rub-impact device, 4, rotor system, 5, rotary encoder, 6, motor, 7, rotor control system, 8, data acquisition system, 9, slip ring, 10, base motion platform, 11, base motion platform control system, 12, rolling shaft.As for the data acquisition system, two eddy current displacement sensors are used to measure the bending displacements of the rotating shaft in the horizontal and vertical directions. The sampling frequency of the displacement data is set to 2,000 Hz. The rotary encoder can pick up the real-time rotating speed of the shaft. Therefore, one rotary encoder is installed at the left shaft end of the rotor and the other one is mounted at the right shaft end. The torsional responses of the rotor system are expressed by the relative rotation speed difference of the two shaft endsΔω=ω1−ω5. The sampling frequency of the rotating speeds of the two shaft ends is set to 12,800 Hz. Because of the slowly time-varying influence of the maneuvering flight, a high-pass filter (frequency>2 Hz) is applied to the vibration data of the rotor test rig before the time-frequency processing.


### 4.2. Experimental Results Analysis

The rotating speed of the rotor test rig is set toω = 1200 rpm, and the rolling speed of the base motion platform is chosen as ωb = π rad/s. The experimental process of the transient barrel roll flight of the base motion platform is the same as the simulation process: at first (t∈0s,2s), the base motion platform keeps the steady-state flight, and it enters into the barrel roll flight at t = 2 s and quits the barrel roll flight at t = 4 s. Then, it returns to the steady flight again (t∈4s,6s).In the first place, the model of the barrel roll flight is verified. As shown in Figure19, the rotor test rig behaves steady responses with small amplitude firstly (t∈0s,2s). Then, with the base motion platform entering into the barrel roll flight t∈2s,4s, the axis orbit of the rotor test rig is deviated entirely. After the barrel roll flight t∈4s,6s, the vibration of the rotor test rig returns to steady responses again. The axis orbits of the rotor test rig influenced by the maneuvering flight are basically the same as the simulation results shown in Figure 8, which demonstrates the validity of the mathematical model of the barrel roll flight.Figure 19
Transient axis orbit of the rotor test rig during the flight of the base motion platformt∈0s,6s.In the second place, the bending and torsional responses of the rotor test rig are investigated to validate the simulation results. The transient time histories of the bending and torsional motions are shown in Figures20(a) and 20(b), respectively. From Figure 20(b), we can see that the amplitude of the torsional responses increase evidently during the barrel roll flight which indicates the occurrence of the rub-impact fault.Figure 20
Transient time history of the bending and torsional motions of the rotor test fig during the flight of the base motion platformt∈0s,6s. (a) Bending motion and (b) torsional motion.
(a)(b)To further reveal the rub-impact fault characteristics of the bending and torsional responses, the transient time-frequency waterfall curves of them are shown in Figures21(a) and 21(b), respectively. As shown in Figure 21(a), the bending responses of the rotor test rig show no apparent fault frequency. Thus, the bending responses cannot be used for rub-impact fault detection at this rotating speed. However, the torsional responses of the rotor test rig in Figure 21(b) display rotating fault frequency ω during the barrel roll flight, which indicates the occurrence of the rub-impact fault. According to the simulation analysis in Figure 15, the frequency spectrums of the bending responses do not show apparent rub-impact fault characteristics at the low rotating speed region, but the distinguishable rotating fault frequency ωlow would be transmitted to the torsional motion because of the influence of the rub-impact. Thus, the bending and torsional vibration phenomena of the rotor test rig are coincident with the frequency characteristics of the simulation results, which demonstrates the dynamics model established in the study is effective and suitable.Figure 21
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor test rig during the flight of the base motion platformt∈0s,6s. (a) Bending motion and (b) torsional motion.
(a)(b)The transient torsional responses, att∈1s,1.2s (steady flight) and at t∈2.5s,2.7s (maneuver flight), are separately chosen for the spectrums. As shown in Figure 22, as the rotor system transforms from nonrubbing to rub-impact, the amplitude of the rotating fault frequency ω increases from 1.016 rpm to 8.328 rpm in which the increment reaches 719.7%. Therefore, the torsional responses can be better used for rub-impact fault diagnose at the low speed region. The experimental results verify the validity of the simulation results qualitatively.Figure 22
Spectrums of the transient torsional motions att∈1s,1.2s (steady flight) and at t∈2.5s,2.7s (maneuvering flight). (a) t∈1s,1.2s and (b) t∈2.5s,2.7s.
(a)(b)


## 4.1. Experimental Device

In this section, an experiment is carried out to verify the results obtained by the numerical simulation. Figure17(a) shows the main experimental devices, which consist of a rotor test rig, a rub-impact device, a base motion platform, two rotary encoders, two displacement sensors, and a set of data acquisition system. Figure 17(b) displays the rub-impact device. The stator of the rub-impact device is a circular truncated cone, whose inner race is slightly larger than that of the disc. The material of the stator is steel. The rub-impact device can be moved from right to left, so the clearance of the stator and disc can be adjusted from 0 to 5 mm. The rub-impact process can be well imitated because of the convenient clearance adjustment.Figure 17
Main experimental devices and rub-impact device. (a) Main experimental devices and (b) rub-impact device.
(a)(b)Figure18 shows the schematic diagram of the main experimental device. The rub-impact rotor system installed on the base motion platform is used for imitating the aeroengine mounted on the aircraft. A single-rotor test rig is adopted to verify the simulation results qualitatively. There are two sets of control system: rotor motor control system and base motion platform control system. Then, with the rotor at the constant speed, the whole rotor system can also rotate around the rolling shaft of the base motion platform, so as to simulate the vibration of the aero-engine during the barrel roll flight. The rolling radius and maximum rolling speed of the base motion platform are 0.4 m and π rad/s, respectively. In addition, to carry out the signal acquisition and the electric power supply of the motor, a slip ring is added to the base motion platform system, through which the above function can be achieved on the premise of the barrel roll flight of the base motion platform.Figure 18
Schematic diagram of the experimental device. 1, rotary encoder, 2, displacement sensors, 3, rub-impact device, 4, rotor system, 5, rotary encoder, 6, motor, 7, rotor control system, 8, data acquisition system, 9, slip ring, 10, base motion platform, 11, base motion platform control system, 12, rolling shaft.As for the data acquisition system, two eddy current displacement sensors are used to measure the bending displacements of the rotating shaft in the horizontal and vertical directions. The sampling frequency of the displacement data is set to 2,000 Hz. The rotary encoder can pick up the real-time rotating speed of the shaft. Therefore, one rotary encoder is installed at the left shaft end of the rotor and the other one is mounted at the right shaft end. The torsional responses of the rotor system are expressed by the relative rotation speed difference of the two shaft endsΔω=ω1−ω5. The sampling frequency of the rotating speeds of the two shaft ends is set to 12,800 Hz. Because of the slowly time-varying influence of the maneuvering flight, a high-pass filter (frequency>2 Hz) is applied to the vibration data of the rotor test rig before the time-frequency processing.


## 4.2. Experimental Results Analysis

The rotating speed of the rotor test rig is set toω = 1200 rpm, and the rolling speed of the base motion platform is chosen as ωb = π rad/s. The experimental process of the transient barrel roll flight of the base motion platform is the same as the simulation process: at first (t∈0s,2s), the base motion platform keeps the steady-state flight, and it enters into the barrel roll flight at t = 2 s and quits the barrel roll flight at t = 4 s. Then, it returns to the steady flight again (t∈4s,6s).In the first place, the model of the barrel roll flight is verified. As shown in Figure19, the rotor test rig behaves steady responses with small amplitude firstly (t∈0s,2s). Then, with the base motion platform entering into the barrel roll flight t∈2s,4s, the axis orbit of the rotor test rig is deviated entirely. After the barrel roll flight t∈4s,6s, the vibration of the rotor test rig returns to steady responses again. The axis orbits of the rotor test rig influenced by the maneuvering flight are basically the same as the simulation results shown in Figure 8, which demonstrates the validity of the mathematical model of the barrel roll flight.Figure 19
Transient axis orbit of the rotor test rig during the flight of the base motion platformt∈0s,6s.In the second place, the bending and torsional responses of the rotor test rig are investigated to validate the simulation results. The transient time histories of the bending and torsional motions are shown in Figures20(a) and 20(b), respectively. From Figure 20(b), we can see that the amplitude of the torsional responses increase evidently during the barrel roll flight which indicates the occurrence of the rub-impact fault.Figure 20
Transient time history of the bending and torsional motions of the rotor test fig during the flight of the base motion platformt∈0s,6s. (a) Bending motion and (b) torsional motion.
(a)(b)To further reveal the rub-impact fault characteristics of the bending and torsional responses, the transient time-frequency waterfall curves of them are shown in Figures21(a) and 21(b), respectively. As shown in Figure 21(a), the bending responses of the rotor test rig show no apparent fault frequency. Thus, the bending responses cannot be used for rub-impact fault detection at this rotating speed. However, the torsional responses of the rotor test rig in Figure 21(b) display rotating fault frequency ω during the barrel roll flight, which indicates the occurrence of the rub-impact fault. According to the simulation analysis in Figure 15, the frequency spectrums of the bending responses do not show apparent rub-impact fault characteristics at the low rotating speed region, but the distinguishable rotating fault frequency ωlow would be transmitted to the torsional motion because of the influence of the rub-impact. Thus, the bending and torsional vibration phenomena of the rotor test rig are coincident with the frequency characteristics of the simulation results, which demonstrates the dynamics model established in the study is effective and suitable.Figure 21
Transient time-frequency waterfall curves of the bending and torsional motions of the rotor test rig during the flight of the base motion platformt∈0s,6s. (a) Bending motion and (b) torsional motion.
(a)(b)The transient torsional responses, att∈1s,1.2s (steady flight) and at t∈2.5s,2.7s (maneuver flight), are separately chosen for the spectrums. As shown in Figure 22, as the rotor system transforms from nonrubbing to rub-impact, the amplitude of the rotating fault frequency ω increases from 1.016 rpm to 8.328 rpm in which the increment reaches 719.7%. Therefore, the torsional responses can be better used for rub-impact fault diagnose at the low speed region. The experimental results verify the validity of the simulation results qualitatively.Figure 22
Spectrums of the transient torsional motions att∈1s,1.2s (steady flight) and at t∈2.5s,2.7s (maneuvering flight). (a) t∈1s,1.2s and (b) t∈2.5s,2.7s.
(a)(b)


## 5. Conclusions

In this study, the transient time-frequency characteristics of the bending and torsional responses of a rub-impact dual-rotor system are investigated considering the barrel roll flight of an aircraft. The equations of the bending-torsional coupling motion are established with the consideration of a low-pressure rotor rub-impact fault. The Runge-Kutta method is adopted to detect the bending and torsional responses, which are influenced by the rub-impact stiffness, friction coefficient, and rotating speed. Some conclusions are as follows:(1)
The unbalanced force and the friction torque caused by the rub-impact fault are the factors leading to the bending-torsional coupling motion of the rotor system. Under the influences of the coupled factors, the torsion responses display the fault frequencies similar to the bending responses, including fractional fault frequencies and rotating fault frequencies.(2)
The rub-impact stiffness affects bending vibration significantly and the torsional motion is sensitive to the friction coefficient. Therefore, the torsional responses show fault frequencies within the whole region of the rub-impact stiffness. The bending responses can only display fault characteristic frequencies at certain rub-impact stiffness.(3)
For the rotating speed, the bending vibration of the rub-impact rotor system will go through the bifurcation behavior [5, 25], in which the bending vibrations show no apparent fault frequency at low-speed region and high-speed region. However, the torsional motion displays fault frequencies within the whole region of the rotating speed.(4)
An experiment is carried out, verifying that the torsional responses work better for the rub-impact fault detection at the low speed region.



---

*Source: 8817228-2021-01-07.xml*