# Study on Fatigue Degradation Behavior of a Cracked Rotor Subjected to Lateral Vibration

Authors: Qian Niu; Shi-xi Yang
Journal: Shock and Vibration
 (2018)

---

## Abstract

Fatigue crack in a rotary shaft is a common failure observed in rotor systems. Since vibration of the shaft causes alternating fatigue loads, the crack propagates slowly. Meanwhile, the propagating crack may cause nonlinear or unstable vibration of the rotor system. In fact, growth of the crack and vibration of the shaft are coupled with each other. Hence, it is necessary to study the fatigue degradation behavior of the cracked rotor accounting for this coupling effect. In this paper, a coupling model of rotor vibration and crack growth is established through dynamic theory and fracture mechanics, and a sequential iterative procedure is proposed to solve the coupling model. Then, the competing degradation failure mode of the cracked rotor is analyzed with regard to the rapid crack growth failure and the unstable vibration failure. And degradation measures are proposed based on the competing degradation failure criterion. At last, degradation behaviors with the coupling effect of nonlinear vibration behavior and multiple parameters including rotation speed, unbalance eccentricity and orientation angle, and damping are investigated by numerical simulation. The results indicate that nonlinear vibration behavior and multiple parameters have considerable influence on the degradation behaviors, which present complex regularity. The findings are of significance to guide the safety design of the rotor system for long time operation and help to the further research on prognostics and lifetime prediction.


---

## Body

## 1. Introduction

Rotary machinery is a kind of common mechanical structure, which is widely used in industry equipments. Because of the rotary motion, the rotor suffers from the alternating fatigue loads. The initiation and further propagation of the fatigue crack are caused subsequently. The presence and propagation of the crack lead to nonlinear, complicated, or even unstable vibration behavior and the performance degradation of the rotor, which may cause severe damage to mechanical components or even catastrophic accidents.Many researchers studied the modeling issues of cracked rotors and ulteriorly analyzed the characteristics of vibration behavior. Wauer gave a literature survey on the dynamics of cracked rotors published before 1990 [1]. Dimarogonas reviewed the research on vibration of cracked structures especially cracked rotors [2]. Papadopoulos reviewed the strain energy release approach for modeling cracks in rotors and vibration behavior of cracked rotors [3]. Kumar and Rostagi gave a brief review on various approaches for analyzing dynamic characteristics of cracked rotors [4]. Most studies are based on an assumption that the crack is relatively static; that is, the crack does not propagate or grow during the operation. Since the propagation of a crack is a much slower process than the vibration of the rotor, the above assumption is reasonable when studying the vibration behavior of the rotor during a short period of time. However, when assessing the performance degradation level of the cracked rotor, predicting the fatigue lifetime, conducting the Prognostic Health Management (PHM), and reliability analysis during a long period of the operation, it is required that not only the vibration behavior but also the crack growth behavior is properly addressed.Crack growth causes slow decrease of the stiffness and variation of the dynamic characteristics of a cracked rotor system, which leads to fatigue degradation. During the slow decreasing evolution, the stiffness also changes fast and periodically. Due to the weight and rotation of the rotor, the surface crack in the shaft opens and closes periodically, namely, the breathing effect, which causes the variation of the stiffness in every rotation period. Accounting for the breathing effect, the cracked rotor presents nonlinear vibration behavior. Some researchers focused their studies mainly on unbalance vibration responses of the nonlinear cracked rotor system [5–9] and the nonlinear vibration behavior due to the breathing effect of the crack, including bifurcation, quasi-periodic, and chaotic motions [10–16]. In engineering practice, the crack growth may lead to variation of the vibration behavior; for example, period-one motion may turn to bifurcation or quasi-periodic motion as the crack depth increases. The nonlinear vibration behavior may cause more complicated fatigue loads than the period-one response of constant amplitude. Therefore, when studying the crack growth behavior of a cracked rotor, the practical fatigue loads with respect to the vibration responses should be well addressed. However, when most researchers worked on the issues of the fatigue crack growth in shafts or round bars, the fatigue loads of constant or variable amplitudes were pregiven artificially or hypothetically, rather than the practical service loads in real applications [17–20]. In fact, the slow crack growth affects the fast rotor vibration behavior while the vibration affects the crack growth behavior. Hence, for a cracked rotor system, the slow process of crack growth and the fast process of rotor vibration are interactional and mutually coupled with each other.Apparently, adequate models and approaches are required to account for both the fast rotor dynamics and the slow fatigue degradation of the fast-slow coupled system of the cracked rotor. Sobczyk and Trebicki provided general formulations for solving stochastic response-degradation coupled problem of structures with fatigue-induced stiffness degradation and presented the efficient sequential method for coupled response-degradation analysis [21]. Then, their approach was developed to analyze nonlinear stochastic vibratory systems with stiffness degradation [22] and multidimensional vibratory degrading systems [23]. Luo et al. proposed a model-based prognostic method and used the interacting multiple model to estimate the damage evolution of the coupled system under multiple operation modes [24, 25]. Lei and Zhu analyzed the crack growth behavior of the elastic element in nonlinear systems under random loads, accounting for the stiffness degradation caused by the crack growth [26].However, the coupling models presented by Sobczyk, Luo, and Lei et al. are general ones or only for simple structures. In their study, the essential relationship between degradation measure (e.g., crack depth) and the system stiffness is directly given by experience. Moreover, the general models barely account for the specific characteristics and behavior of the cracked rotor system. Firstly, the dynamic behavior of the cracked rotor is influenced by multiple parameters of the system including not only stiffness but also rotation speed, unbalance eccentricity and orientation angle, and damping. The multiple parameters may have a complicated and coupling effect on the dynamic characteristics and hence the degradation behavior of the rotor system. Secondly, the only consideration of the stiffness degradation caused by crack growth is inadequate because the stiffness also presents periodic characteristic due to the breathing effect of the crack. The nonlinear vibration behavior induced by the breathing crack may have a remarkable impact on the degradation behavior. Thirdly, the crack growth disciplinarian of dynamic rotary shafts is more complicated than that of static simple structures, like bars and plates. Consequently, it is difficult to provide simplified experience formulas illustrating the complex relationship between the crack size and the system stiffness and the relationship between the vibration responses and the fatigue loads. These relationships are required to be revealed by corresponding physical models. Therefore, with all the above considerations, it is of great significance to establish a specific and applicable fast-slow coupling model of the cracked rotor and then to analyze the degradation behavior with the coupling effect of nonlinear vibration behavior and multiple parameters.Moreover, before analyzing the degradation behavior of the cracked rotor, the degradation failure criterion and appropriate degradation measures are necessary to be addressed. Sobczyk and Trebicki defined degradation measure based on the crack depth, and fracture criterion was regarded as the degradation failure criterion [21, 22]. However, it may be inadequate since the degradation failure may occur before the crack grows to the critical depth due to severe or unstable vibration behavior. To insure the stability of the vibration, some researchers worked on the stability issues of the cracked rotor with multiple parameters [27–29]. Nevertheless, unstable vibration failure may not occur even if the crack grows over the critical depth in some cases. Therefore, the competing degradation failure modes caused by both rapid crack growth and unstable vibration should be well considered. Then, degradation measures should be properly proposed accounting for the competing degradation failure criterion. Moreover, the unstable vibration of the cracked rotor presents bifurcations of different types [29], which may have different influences on the degradation behavior.The main objective of this article is twofold. One is to establish the fast-slow coupling model accounting for nonlinear dynamics of the cracked rotor and the fatigue propagation of the transverse surface crack. The other is to analyze the degradation behavior with the coupling effect of the nonlinear vibration behavior and the multiple parameters of the rotor system. Firstly, the coupling model of a Jeffcott rotor with a transverse surface crack is established based on dynamic theory and fracture mechanics. And a sequential iterative procedure is proposed to approximately solve the coupling model. Then, the degradation measure of the cracked rotor is defined based on the competing degradation failure criterion which accounts for both the unstable vibration of the rotor and the rapid growth of the crack. At last, the degradation behaviors with the coupling effect of multiple parameters of the rotor system including rotation speed, unbalance eccentricity and orientation angle, and damping are analyzed through numerical simulation and several new findings are reported.


## 2. Problem Description

We consider a cracked rotor, depicted in Figure1, modeled as a Jeffcott rotor system in which the rotor is simplified to a centered rigid disk on a massless elastic shaft and a transverse surface crack is located at the midspan of the rotor, which propagates as the shaft rotates. The initial of the crack is not considered in this study, and only the fatigue degradation behavior caused by the macroscopic crack growth is analyzed. For simplifying the analysis and calculation, several assumptions are made in this study:(1)
Torsional vibration and axial vibration of the shaft and shearing action are omitted. Only lateral vibration and crack growth under mode I are considered.(2)
The front of the crack is a straight line and maintained during the growth.(3)
The linear elastic fracture mechanics is employed for the analysis due to the applied stress being much less than the yielding stress.(4)
Only the plane strain condition is considered at the crack front due to the geometry constraint.(5)
The thickness of the disk is neglected, and the displacement at the crack location is the same as the disk.Figure 1
Schematic of the cracked rotor.


## 3. Model Establishment

The physical model of a Jeffcott rotor with a propagating surface transverse crack is established in two steps. Firstly, the dynamic model of the cracked rotor is established and then the crack growth model. The two models are coupled with each other due to the growing crack, which will be described in detail in this section. Noting that the main creative study focuses on the coupling effect of rotor dynamics and crack growth on the degradation behavior, the widely-used classical models of the cracked rotor and the crack growth are adopted for simplicity.


### 3.1. Vibration Differential Equations of the Cracked Rotor

The mass of the Jeffcott rotorm is assumed to be concentrated on the disk in the midspan. Let x-y be the fixed coordinates and ξ−η be the rotating coordinates in the crack cross section as shown in Figure 2.Figure 2
Coordinate system.The unbalance eccentricity of the rotor isε, which is oriented at an angle β (unbalance orientation angle) with the ξ coordinate direction. θ is the instantaneous angle of rotation at time t (θ=ωt), and ω is the rotation speed. c is the damping of the rotor system. Due to the existence of the crack, kx and ky are the stiffness of the x and y coordinate directions, respectively, and kxy and kyx are the cross stiffness. The dynamics of the rotor is combined by the translation of the centroid and the rotation around it. The dynamic equations for the cracked Jeffcott rotor can be written in the fixed coordinate system as(1)my¨+cy˙+kyy+kyxx=mεω2⁡cosθ+β−mg,mx¨+cx˙+kxx+kxyy=mεω2⁡sinθ+β.The stiffness coefficients of (1) in the fixed coordinate system can be obtained from the stiffness coefficients kξ, kη, kξη, and kηξ defined in the rotational frame ξ−η, by the transformation matrix T:(2)kykyxkxykx=T−1kξkξηkηξkηT,where the transformation matrix T is(3)T=cos⁡θsin⁡θ−sin⁡θcos⁡θ.Cross section of the cracked rotor at the crack location is showed in Figure3, where a is the depth of the crack, D is the diameter of the rotor, α′ is the depth of the crack at the location where the η coordinate is w, and α′=D2−2w2.The stiffness matrix is transformed from the flexibility matrix which can be derived by the strain energy release approach [3, 9, 16] (Appendix A). The flexibility coefficients of the cracked rotor can be derived as(4)gξ=L348EI+∬128L2α′2α2EπD8Fαα′2dα⁡dw,gη=L348EI+∬512L2w2αEπD8F′αα′2dα⁡dw,gξη=gηξ=∬256L2α′wαEπD8Fαα′F′αα′dα⁡dw.Figure 3
Cross-section at the crack location.It is worth noting that the integration limits in (4) depend on the open part of the crack which varies over time. On the one hand, the crack opens and closes periodically because of the breathing effect. One the other hand, since the crack propagates, the depth of the crack (at) increases with time. Therefore, the integration limits in (4) depend on two aspects: the breathing effect and the propagation of the crack. Temporally, the propagation of the crack is omitted in this section, and only the breathing effect is considered, so the crack depth at becomes a. Several crack breathing models have been proposed, such as the harmonic model proposed by Mayes and Davies [1, 6, 30, 31], the crack closure line (CCL) model proposed by Darpe et al. [9, 16], and some new models with consideration of unbalance force [32]. In most applications, the gravity of the rotor is much larger than the unbalance force. In this situation, where the gravity is dominant, the harmonic model is appropriate and applied here. The stiffness matrix in the rotational frame is as follows:(5)kξ00kη=k000k0−FθΔkξ00Δkη,where(6)Fθ=121−cosθ,(7)Δkξ=k0−k^ξ,Δkη=k0−k^η,in which k0 is the stiffness of the uncracked rotor, and k0=1/g0; k^ξ and k^η are the stiffness in the coordinate directions ξ and η when the crack is totally open. At fully open state of the crack, the flexibility coefficients g^ξ, g^η, and g^ξη can be calculated by (4), in which the integration region is the arched area of the open part of the crack. And the stiffness coefficients are(8)k^ξ=g^ηg^ξg^η−g^ξη2,k^η=g^ξg^ξg^η−g^ξη2.Using (1)–(3) and (5)–(7), the equations of motion for the cracked rotor in a fixed frame x-y can be written as follows:(9)my¨+cy˙+k0−0.5Fθ2k0−k^ξ−k^η+k^η−k^ξcos⁡2θy−0.5Fθk^η−k^ξsin⁡2θx=mεω2cosθ+β−mg,mx¨+cx˙+k0−0.5Fθ2k0−k^ξ−k^η−k^η−k^ξcos⁡2θx−0.5Fθk^η−k^ξsin⁡2θy=mεω2sinθ+β.It is important to note that stiffness coefficients in the dynamic model of the cracked rotor display periodic variation. And the period is equal to the rotation periodT of the rotor, T=2π/ω. Furthermore, stiffness coefficients decrease slowly due to the growth of the crack, which will be illustrated in detail in Section 3.2. Therefore in real applications, stiffness coefficients in the dynamic model of the cracked rotor include two parts: the periodic part of fast variation and the decreasing part of slow variation. Consequently, the dynamic stiffness is a very complex function unable to be expressed in explicit formulations. The complex stiffness coefficients lead to complex vibration and degradation behaviors and bring difficulty to the solution of the dynamic model.


### 3.2. Propagation Equations of the Transverse Surface Crack

During the long time of operation, the transverse crack at the surface of the shaft grows slowly due to the fatigue stresses caused by the vibration of the rotor. The time history of the crack propagation is a very slow process, which can be regarded as the fatigue degradation evolution of the cracked rotor. A large majority of fatigue crack growth experiments on small specimens have shown that the fatigue growth of a crack in the metallic material is mainly controlled by the range of stress intensity factor (ΔK) of the crack front, although many factors, such as mean stress, load ratio and frequency, random loading, multiaxial and complex stress states as well as environments (corrosion and temperature), plastic zone, and microstructure, affect the manner of fatigue crack growth. When ΔK is larger than the threshold of the fatigue crack growth (ΔKth), the crack growth behavior can be described by the famous Paris–Erdogan laws [33, 34]:(10)dadN=CΔKB,where a is the depth of the crack, N is the number of the stress cycles, and C and B are empirical parameters accounting for the material and environmental effects and in practice, can be obtained from small specimens by a regulated procedure, such as that proposed in ASTM E647-15 [35].Δ
K can be calculated by (A.5)–(A.9). The bending force acting on the crack cross section Qξ and Qη can be calculated as follows:(11)Qξ=kξξ+kξηη,Qη=kηη+kηξξ,where the vibration responses ξ and η in the rotation frame can be obtained by a transformation from the rotor responses x and y in the fixed frame:(12)ξη=cos⁡θ−sin⁡θsin⁡θcos⁡θyx.k
ξ, kη, kξη, and kηξ are stiffness of ξ and η coordinate directions and crossing direction, respectively:(13)kξ=gηgξgη−gξη2,kη=gξgξgη−gξη2,kξη=kηξ=−gξηgξgη−gξη2.From (A.6), one can conclude that the SIF at any location of the crack front is related to the local stress, the geometry size of the crack, and the specific position along the crack front. Therefore, ΔK at different locations is unequal, the same as the crack growth rates, which cause the change of the shape of the crack front. When the crack grows, the new shape of the crack front is only determined by the current crack front and has nothing to do with the initial crack front [19, 35]. In [33], it can be found that, under bending fatigue stress, the shape of the crack front tends to be a straight line as the crack grows to the center of the shaft. And many researchers assumed the crack front to be straight for simplification when studying cracked rotors [3, 36–38]. In this paper, for convenience of the analysis, it is assumed the crack front maintains straight when the crack grows.Based on the assumption that the crack front maintains straight, an appropriate equivalent SIF is necessary to be defined to apply the Paris–Erdogan laws. In [39], the researchers derived an average SIF expression as an equivalent SIF from the SIF expressions of the edge and middle points of the crack front through Lagrange polynomial interpolation. In this study, in order to assure the safety of the rotor system, a conservative equivalent SIF is defined, which is the SIF of the point where the crack grows fastest. From [19, 40], regardless of the stress range Δσ, the SIF of the edge point is larger than that of the middle point. Moreover, in the condition of this study, it can be derived from (A.7) that the stress range Δσ of the edge point is larger than that of the middle point. Therefore, ΔK at the edge point of the crack front is the largest, which means the crack at the edge point grows fastest, and the SIF of the edge point is chosen to be the equivalent SIF.The edge point of the crack front actually refers to a point very close to the intersection of the crack front with the free surface of the shaft [19]. When the depth of the crack is a, the η coordinate of the edge point is w=±aD−a2, the crack depth rate α/α′ converges to zero from the positive direction, and the limits of the functions F and F′ are calculated as follows:(14)limα/α′⟶0+Fαα′=limα/α′⟶0+F′αα′=1.1220.From (A.5)–(A.7), the equivalent SIF Ke is written as follows:(15)Ke=1.1220σeπa,where σe is the equivalent stress and expressed as follows, which is the total stress at the edge point:(16)σe=σξ+ση=L4ID2−aQξ+aD−a2Qη.It is worth noting that the depth of the cracka is growing over time actually and should be written as at. From (15) and (16), one can derive that the equivalent SIF Ke is dominant by the crack depth at which varies slowly and the equivalent stress σe which varies fast. Using (10), the crack depth aN after N (N=1,⁡2,⁡3,⁡…) stress cycles (or rotation periods) can be calculated by integration from the crack depth aN−1 as follows:(17)aN=aN−1+∫N−1NCΔKeNB⁡dN,where ΔKeN is the range of the equivalent SIF in the Nth rotation period.


### 3.3. Solution of the Coupling Model

From Sections3.1 and 3.2, the degrading cracked rotor system can be expressed as a hierarchical dynamic system form:(18a)q˙=fq,⁡μΦ,(18b)Φ˙=νgΦ,⁡q,where q=yxT is a fast-time variable; Φ is a slow-time variable and refers to the crack depth a in this study; μ is a function of Φ representing stiffness matrix; ν is a small rate constant which defines a time scale separation between fast and slow-time dynamics; and the overdot denotes differentiation with respect to time t. Over the time scales of Ον, we consider (18a) to be quasi-stationary since drifts in μ are negligible.The rotor vibration equations and the crack propagation equations are coupled with each other, forming the vibration-degradation coupling model. Since it is difficult to describe the functionμ and g by explicit function expressions, the customary approximate analytic procedures cannot be applied to solve the coupling model. Alternatively, a numerical procedure is proposed here, which is based on sequential iteration. In engineering applications, the crack grows very slowly compared to the natural time scale of the rotor vibration process. It means that the crack depth can be regarded as constant in a short time interval. In this study, it is assumed that the crack depth aM−1 is constant in n rotation periods (time interval nT), and after n periods, the crack depth grows to aM. Based on this assumption, the stiffness coefficients of the vibration equations dominated by the crack depth can be calculated in the time interval of n rotation periods, and vibration differential equations (9) can be solved through various existing analytic or numerical procedures. The solutions of (9) are vibration responses (xt,⁡yt) of the cracked rotor and are used to calculate the time history KeMt (tM−1<t<tM) by (11)–(16). Then, using (17), aM is calculated as follows:(19)aM=aM−1+∫M−1MCΔKeMNB⁡dN,where ΔKeMN denotes the SIF range of the Nth rotation period, which is dependent on KeMt. At last, given the initial crack depth a0, the crack growth process described by the crack depth aM (M=1,⁡2,⁡3,⁡…) of each time interval can be obtained by the iteration procedure.It is worth noting thatKeMt is determined by xt⁡and⁡yt in the time interval tM−1<t<tM. When the vibration of the rotor is period-one motion, the amplitude of each period is equal. Therefore, ΔKeMN of every rotation period is also the same, and in this situation, (19) can be written as follows:(20)aM=aM−1+nCΔKeMNB.However, when the vibration behavior of the cracked rotor is complex, such as bifurcation motion, quasi periodic motion, and chaotic motion [10–16], KeMt of each rotation period is not the same anymore. In this situation, a statistical counting method is necessary to be applied to count ΔKeMN from KeMt. In this study, the rain-flow counting method [41, 42] is applied which is widely used in research and real applications. Then, the integration part of (19) can be calculated.


## 3.1. Vibration Differential Equations of the Cracked Rotor

The mass of the Jeffcott rotorm is assumed to be concentrated on the disk in the midspan. Let x-y be the fixed coordinates and ξ−η be the rotating coordinates in the crack cross section as shown in Figure 2.Figure 2
Coordinate system.The unbalance eccentricity of the rotor isε, which is oriented at an angle β (unbalance orientation angle) with the ξ coordinate direction. θ is the instantaneous angle of rotation at time t (θ=ωt), and ω is the rotation speed. c is the damping of the rotor system. Due to the existence of the crack, kx and ky are the stiffness of the x and y coordinate directions, respectively, and kxy and kyx are the cross stiffness. The dynamics of the rotor is combined by the translation of the centroid and the rotation around it. The dynamic equations for the cracked Jeffcott rotor can be written in the fixed coordinate system as(1)my¨+cy˙+kyy+kyxx=mεω2⁡cosθ+β−mg,mx¨+cx˙+kxx+kxyy=mεω2⁡sinθ+β.The stiffness coefficients of (1) in the fixed coordinate system can be obtained from the stiffness coefficients kξ, kη, kξη, and kηξ defined in the rotational frame ξ−η, by the transformation matrix T:(2)kykyxkxykx=T−1kξkξηkηξkηT,where the transformation matrix T is(3)T=cos⁡θsin⁡θ−sin⁡θcos⁡θ.Cross section of the cracked rotor at the crack location is showed in Figure3, where a is the depth of the crack, D is the diameter of the rotor, α′ is the depth of the crack at the location where the η coordinate is w, and α′=D2−2w2.The stiffness matrix is transformed from the flexibility matrix which can be derived by the strain energy release approach [3, 9, 16] (Appendix A). The flexibility coefficients of the cracked rotor can be derived as(4)gξ=L348EI+∬128L2α′2α2EπD8Fαα′2dα⁡dw,gη=L348EI+∬512L2w2αEπD8F′αα′2dα⁡dw,gξη=gηξ=∬256L2α′wαEπD8Fαα′F′αα′dα⁡dw.Figure 3
Cross-section at the crack location.It is worth noting that the integration limits in (4) depend on the open part of the crack which varies over time. On the one hand, the crack opens and closes periodically because of the breathing effect. One the other hand, since the crack propagates, the depth of the crack (at) increases with time. Therefore, the integration limits in (4) depend on two aspects: the breathing effect and the propagation of the crack. Temporally, the propagation of the crack is omitted in this section, and only the breathing effect is considered, so the crack depth at becomes a. Several crack breathing models have been proposed, such as the harmonic model proposed by Mayes and Davies [1, 6, 30, 31], the crack closure line (CCL) model proposed by Darpe et al. [9, 16], and some new models with consideration of unbalance force [32]. In most applications, the gravity of the rotor is much larger than the unbalance force. In this situation, where the gravity is dominant, the harmonic model is appropriate and applied here. The stiffness matrix in the rotational frame is as follows:(5)kξ00kη=k000k0−FθΔkξ00Δkη,where(6)Fθ=121−cosθ,(7)Δkξ=k0−k^ξ,Δkη=k0−k^η,in which k0 is the stiffness of the uncracked rotor, and k0=1/g0; k^ξ and k^η are the stiffness in the coordinate directions ξ and η when the crack is totally open. At fully open state of the crack, the flexibility coefficients g^ξ, g^η, and g^ξη can be calculated by (4), in which the integration region is the arched area of the open part of the crack. And the stiffness coefficients are(8)k^ξ=g^ηg^ξg^η−g^ξη2,k^η=g^ξg^ξg^η−g^ξη2.Using (1)–(3) and (5)–(7), the equations of motion for the cracked rotor in a fixed frame x-y can be written as follows:(9)my¨+cy˙+k0−0.5Fθ2k0−k^ξ−k^η+k^η−k^ξcos⁡2θy−0.5Fθk^η−k^ξsin⁡2θx=mεω2cosθ+β−mg,mx¨+cx˙+k0−0.5Fθ2k0−k^ξ−k^η−k^η−k^ξcos⁡2θx−0.5Fθk^η−k^ξsin⁡2θy=mεω2sinθ+β.It is important to note that stiffness coefficients in the dynamic model of the cracked rotor display periodic variation. And the period is equal to the rotation periodT of the rotor, T=2π/ω. Furthermore, stiffness coefficients decrease slowly due to the growth of the crack, which will be illustrated in detail in Section 3.2. Therefore in real applications, stiffness coefficients in the dynamic model of the cracked rotor include two parts: the periodic part of fast variation and the decreasing part of slow variation. Consequently, the dynamic stiffness is a very complex function unable to be expressed in explicit formulations. The complex stiffness coefficients lead to complex vibration and degradation behaviors and bring difficulty to the solution of the dynamic model.


## 3.2. Propagation Equations of the Transverse Surface Crack

During the long time of operation, the transverse crack at the surface of the shaft grows slowly due to the fatigue stresses caused by the vibration of the rotor. The time history of the crack propagation is a very slow process, which can be regarded as the fatigue degradation evolution of the cracked rotor. A large majority of fatigue crack growth experiments on small specimens have shown that the fatigue growth of a crack in the metallic material is mainly controlled by the range of stress intensity factor (ΔK) of the crack front, although many factors, such as mean stress, load ratio and frequency, random loading, multiaxial and complex stress states as well as environments (corrosion and temperature), plastic zone, and microstructure, affect the manner of fatigue crack growth. When ΔK is larger than the threshold of the fatigue crack growth (ΔKth), the crack growth behavior can be described by the famous Paris–Erdogan laws [33, 34]:(10)dadN=CΔKB,where a is the depth of the crack, N is the number of the stress cycles, and C and B are empirical parameters accounting for the material and environmental effects and in practice, can be obtained from small specimens by a regulated procedure, such as that proposed in ASTM E647-15 [35].Δ
K can be calculated by (A.5)–(A.9). The bending force acting on the crack cross section Qξ and Qη can be calculated as follows:(11)Qξ=kξξ+kξηη,Qη=kηη+kηξξ,where the vibration responses ξ and η in the rotation frame can be obtained by a transformation from the rotor responses x and y in the fixed frame:(12)ξη=cos⁡θ−sin⁡θsin⁡θcos⁡θyx.k
ξ, kη, kξη, and kηξ are stiffness of ξ and η coordinate directions and crossing direction, respectively:(13)kξ=gηgξgη−gξη2,kη=gξgξgη−gξη2,kξη=kηξ=−gξηgξgη−gξη2.From (A.6), one can conclude that the SIF at any location of the crack front is related to the local stress, the geometry size of the crack, and the specific position along the crack front. Therefore, ΔK at different locations is unequal, the same as the crack growth rates, which cause the change of the shape of the crack front. When the crack grows, the new shape of the crack front is only determined by the current crack front and has nothing to do with the initial crack front [19, 35]. In [33], it can be found that, under bending fatigue stress, the shape of the crack front tends to be a straight line as the crack grows to the center of the shaft. And many researchers assumed the crack front to be straight for simplification when studying cracked rotors [3, 36–38]. In this paper, for convenience of the analysis, it is assumed the crack front maintains straight when the crack grows.Based on the assumption that the crack front maintains straight, an appropriate equivalent SIF is necessary to be defined to apply the Paris–Erdogan laws. In [39], the researchers derived an average SIF expression as an equivalent SIF from the SIF expressions of the edge and middle points of the crack front through Lagrange polynomial interpolation. In this study, in order to assure the safety of the rotor system, a conservative equivalent SIF is defined, which is the SIF of the point where the crack grows fastest. From [19, 40], regardless of the stress range Δσ, the SIF of the edge point is larger than that of the middle point. Moreover, in the condition of this study, it can be derived from (A.7) that the stress range Δσ of the edge point is larger than that of the middle point. Therefore, ΔK at the edge point of the crack front is the largest, which means the crack at the edge point grows fastest, and the SIF of the edge point is chosen to be the equivalent SIF.The edge point of the crack front actually refers to a point very close to the intersection of the crack front with the free surface of the shaft [19]. When the depth of the crack is a, the η coordinate of the edge point is w=±aD−a2, the crack depth rate α/α′ converges to zero from the positive direction, and the limits of the functions F and F′ are calculated as follows:(14)limα/α′⟶0+Fαα′=limα/α′⟶0+F′αα′=1.1220.From (A.5)–(A.7), the equivalent SIF Ke is written as follows:(15)Ke=1.1220σeπa,where σe is the equivalent stress and expressed as follows, which is the total stress at the edge point:(16)σe=σξ+ση=L4ID2−aQξ+aD−a2Qη.It is worth noting that the depth of the cracka is growing over time actually and should be written as at. From (15) and (16), one can derive that the equivalent SIF Ke is dominant by the crack depth at which varies slowly and the equivalent stress σe which varies fast. Using (10), the crack depth aN after N (N=1,⁡2,⁡3,⁡…) stress cycles (or rotation periods) can be calculated by integration from the crack depth aN−1 as follows:(17)aN=aN−1+∫N−1NCΔKeNB⁡dN,where ΔKeN is the range of the equivalent SIF in the Nth rotation period.


## 3.3. Solution of the Coupling Model

From Sections3.1 and 3.2, the degrading cracked rotor system can be expressed as a hierarchical dynamic system form:(18a)q˙=fq,⁡μΦ,(18b)Φ˙=νgΦ,⁡q,where q=yxT is a fast-time variable; Φ is a slow-time variable and refers to the crack depth a in this study; μ is a function of Φ representing stiffness matrix; ν is a small rate constant which defines a time scale separation between fast and slow-time dynamics; and the overdot denotes differentiation with respect to time t. Over the time scales of Ον, we consider (18a) to be quasi-stationary since drifts in μ are negligible.The rotor vibration equations and the crack propagation equations are coupled with each other, forming the vibration-degradation coupling model. Since it is difficult to describe the functionμ and g by explicit function expressions, the customary approximate analytic procedures cannot be applied to solve the coupling model. Alternatively, a numerical procedure is proposed here, which is based on sequential iteration. In engineering applications, the crack grows very slowly compared to the natural time scale of the rotor vibration process. It means that the crack depth can be regarded as constant in a short time interval. In this study, it is assumed that the crack depth aM−1 is constant in n rotation periods (time interval nT), and after n periods, the crack depth grows to aM. Based on this assumption, the stiffness coefficients of the vibration equations dominated by the crack depth can be calculated in the time interval of n rotation periods, and vibration differential equations (9) can be solved through various existing analytic or numerical procedures. The solutions of (9) are vibration responses (xt,⁡yt) of the cracked rotor and are used to calculate the time history KeMt (tM−1<t<tM) by (11)–(16). Then, using (17), aM is calculated as follows:(19)aM=aM−1+∫M−1MCΔKeMNB⁡dN,where ΔKeMN denotes the SIF range of the Nth rotation period, which is dependent on KeMt. At last, given the initial crack depth a0, the crack growth process described by the crack depth aM (M=1,⁡2,⁡3,⁡…) of each time interval can be obtained by the iteration procedure.It is worth noting thatKeMt is determined by xt⁡and⁡yt in the time interval tM−1<t<tM. When the vibration of the rotor is period-one motion, the amplitude of each period is equal. Therefore, ΔKeMN of every rotation period is also the same, and in this situation, (19) can be written as follows:(20)aM=aM−1+nCΔKeMNB.However, when the vibration behavior of the cracked rotor is complex, such as bifurcation motion, quasi periodic motion, and chaotic motion [10–16], KeMt of each rotation period is not the same anymore. In this situation, a statistical counting method is necessary to be applied to count ΔKeMN from KeMt. In this study, the rain-flow counting method [41, 42] is applied which is widely used in research and real applications. Then, the integration part of (19) can be calculated.


## 4. Degradation Measures

To evaluate the degradation level and analyze the degradation behavior, quantitative degradation measures are necessary to be proposed. Moreover, the definition of degradation measures is associated with the threshold of the degradation failure or the so-called failure criterion. For the cracked rotor system, a competing degradation failure criterion is described in detail in this section.


### 4.1. Degradation Failure Criterion of the Cracked Rotor

If one wants to define the degradation measure of the cracked rotor, a clear identification of failure criterion is needed. In [39], two failure criteria are considered: (1) geometrical limit: the failure occurs automatically when crack depth is larger than shaft diameter and (2) critical SIF: the fracture occurs when maximum SIF is larger than critical SIF. However in engineering applications, the fracture of the rotor is not allowable in most cases. The rotor should maintain a safe and stable running state. Once the running state becomes risky or unstable, the inspection or maintenance is required. Detecting the dangerous and unstable running states as soon as possible is of great significance, which is an important objective of condition monitoring and prognostics. Obviously, it is unsafe and inappropriate to regard the fracture criterion of the shaft as the degradation failure criterion of the rotor system. That is why the safe factor is necessary in engineering applications.Fatigue experiments illustrate that the crack grows rapidly at increasing growth rate greater than that of the Paris law in a short period time before maximum SIF reaches critical SIF [19]. Consequently, a critical crack depth for rapid crack growth (notated as ac) can be defined through experiments or experience. Once the crack grows deeper than ac, rapid and unstable crack growth occurs. Although fracture does not occur during the time of rapid and unstable crack growth, the running state of the rotor system is risky. Therefore, in this study, a degradation failure criterion based on the growth of crack is defined as a>ac.Moreover, for the vibration behavior of the cracked rotor, unstable vibration occurs when the crack depth reaches some thresholds for some rotation speeds [29], which is undesirable in real applications. The situation possibly exists, in which unstable vibration occurs before the crack depth reaches ac. Thus, it is necessary to address the condition on which unstable vibration occurs. In Section 3.3, it is assumed that the crack depth is constant in the short time interval of nT. In this situation, the stiffness coefficients of the cracked rotor system vary periodically, and Floquet’s theory is applicable for stability analysis [43–45] (Appendix B).Once Floquet multipliers exceed the unit circle on the complex plane, the vibration of the cracked rotor system becomes unstable and bifurcation occurs. As shown in Figure4, different distributions of the Floquet multipliers on the complex plane correspond to different types of bifurcations. There are three types of bifurcations: fold (period doubling), Hopf, and flip (saddle-node), as shown in Figure 4 [29, 37]. When even small damping (e.g., damping ratio ζ = 0.01) exists, the unstable region corresponding to Hopf bifurcation vanishes [37]. Since damping surely exists in engineering applications, only fold and flip are considered in this study.Figure 4
Eigenvalues on the complex plane corresponding to the three types of bifurcation: (a) fold; (b) Hopf; (c) flip.
(a)
(b)
(c)As to period-doubling bifurcation, vibration amplitude does not increase suddenly. In this situation, degradation failure does not occur to the cracked rotor since the vibration amplitude is not overlarge. As to saddle-node bifurcation, vibration amplitude increases suddenly because of the divergence of the vibration responses. In this situation, degradation failure occurs, even if the crack depth has not reachedac. When the vibration behavior reaches the boundary of the unstable region corresponding to the saddle-node bifurcation, the crack depth at this time is regarded as a critical crack depth for unstable vibration (notated as au). Given the necessary parameters in (9) including mass, damping, unbalance eccentricity and orientation angle, and rotation speed, one can address au under this parameter group through Floquet’s theory. Accordingly, another degradation failure criterion based on vibration behavior is defined as a>au.In conclusion, two degradation failure criteria of the cracked rotor are defined based on the crack growth and the rotor vibration, respectively. If one criterion is reached, the degradation failure occurs. Thus, a competing degradation failure criterion is formed.It is worth noting that when the vibration of the cracked rotor is period-one motion, the coupling model can be solved by (20). However, when bifurcation occurs, the rain-flow counting method is acquired to calculate the integration part of (19).


### 4.2. The Proposed Degradation Measure

To define degradation measures, the competing degradation failure criterion accounting for both rapid crack growth and unstable vibration response are considered. According to [21], the degradation measure D (0≤D≤1) is defined based on the evolution of crack growth. When the crack depth is zero, D equals 0, which means no crack exists. When the crack depth reaches ac or au, D equals 1, which means degradation failure occurs. When the crack depth is aM, the degradation measure is written as DM and calculated as follows:(21)DM=aMminacau.Furthermore, the average degradation rate at this time is defined as follows:(22)VM=DM+1−DMn.It is assumed that there is already macroscopic initial crack on the surface of the shaft. The depth of the initial crack is assumed to bea0 and the initial degradation measure is D0 from (21). The fatigue degradation life TD is defined as the number of stress cycles (or rotation periods), during which the cracked rotor degrades from D0 to 1. If DM=1, TD=nM.


## 4.1. Degradation Failure Criterion of the Cracked Rotor

If one wants to define the degradation measure of the cracked rotor, a clear identification of failure criterion is needed. In [39], two failure criteria are considered: (1) geometrical limit: the failure occurs automatically when crack depth is larger than shaft diameter and (2) critical SIF: the fracture occurs when maximum SIF is larger than critical SIF. However in engineering applications, the fracture of the rotor is not allowable in most cases. The rotor should maintain a safe and stable running state. Once the running state becomes risky or unstable, the inspection or maintenance is required. Detecting the dangerous and unstable running states as soon as possible is of great significance, which is an important objective of condition monitoring and prognostics. Obviously, it is unsafe and inappropriate to regard the fracture criterion of the shaft as the degradation failure criterion of the rotor system. That is why the safe factor is necessary in engineering applications.Fatigue experiments illustrate that the crack grows rapidly at increasing growth rate greater than that of the Paris law in a short period time before maximum SIF reaches critical SIF [19]. Consequently, a critical crack depth for rapid crack growth (notated as ac) can be defined through experiments or experience. Once the crack grows deeper than ac, rapid and unstable crack growth occurs. Although fracture does not occur during the time of rapid and unstable crack growth, the running state of the rotor system is risky. Therefore, in this study, a degradation failure criterion based on the growth of crack is defined as a>ac.Moreover, for the vibration behavior of the cracked rotor, unstable vibration occurs when the crack depth reaches some thresholds for some rotation speeds [29], which is undesirable in real applications. The situation possibly exists, in which unstable vibration occurs before the crack depth reaches ac. Thus, it is necessary to address the condition on which unstable vibration occurs. In Section 3.3, it is assumed that the crack depth is constant in the short time interval of nT. In this situation, the stiffness coefficients of the cracked rotor system vary periodically, and Floquet’s theory is applicable for stability analysis [43–45] (Appendix B).Once Floquet multipliers exceed the unit circle on the complex plane, the vibration of the cracked rotor system becomes unstable and bifurcation occurs. As shown in Figure4, different distributions of the Floquet multipliers on the complex plane correspond to different types of bifurcations. There are three types of bifurcations: fold (period doubling), Hopf, and flip (saddle-node), as shown in Figure 4 [29, 37]. When even small damping (e.g., damping ratio ζ = 0.01) exists, the unstable region corresponding to Hopf bifurcation vanishes [37]. Since damping surely exists in engineering applications, only fold and flip are considered in this study.Figure 4
Eigenvalues on the complex plane corresponding to the three types of bifurcation: (a) fold; (b) Hopf; (c) flip.
(a)
(b)
(c)As to period-doubling bifurcation, vibration amplitude does not increase suddenly. In this situation, degradation failure does not occur to the cracked rotor since the vibration amplitude is not overlarge. As to saddle-node bifurcation, vibration amplitude increases suddenly because of the divergence of the vibration responses. In this situation, degradation failure occurs, even if the crack depth has not reachedac. When the vibration behavior reaches the boundary of the unstable region corresponding to the saddle-node bifurcation, the crack depth at this time is regarded as a critical crack depth for unstable vibration (notated as au). Given the necessary parameters in (9) including mass, damping, unbalance eccentricity and orientation angle, and rotation speed, one can address au under this parameter group through Floquet’s theory. Accordingly, another degradation failure criterion based on vibration behavior is defined as a>au.In conclusion, two degradation failure criteria of the cracked rotor are defined based on the crack growth and the rotor vibration, respectively. If one criterion is reached, the degradation failure occurs. Thus, a competing degradation failure criterion is formed.It is worth noting that when the vibration of the cracked rotor is period-one motion, the coupling model can be solved by (20). However, when bifurcation occurs, the rain-flow counting method is acquired to calculate the integration part of (19).


## 4.2. The Proposed Degradation Measure

To define degradation measures, the competing degradation failure criterion accounting for both rapid crack growth and unstable vibration response are considered. According to [21], the degradation measure D (0≤D≤1) is defined based on the evolution of crack growth. When the crack depth is zero, D equals 0, which means no crack exists. When the crack depth reaches ac or au, D equals 1, which means degradation failure occurs. When the crack depth is aM, the degradation measure is written as DM and calculated as follows:(21)DM=aMminacau.Furthermore, the average degradation rate at this time is defined as follows:(22)VM=DM+1−DMn.It is assumed that there is already macroscopic initial crack on the surface of the shaft. The depth of the initial crack is assumed to bea0 and the initial degradation measure is D0 from (21). The fatigue degradation life TD is defined as the number of stress cycles (or rotation periods), during which the cracked rotor degrades from D0 to 1. If DM=1, TD=nM.


## 5. Numerical Simulation and Parametric Sensitivity

### 5.1. Procedure and Parameters of Numerical Simulation

The degradation behavior of the cracked rotor described by the proposed coupling model is addressed by numerical simulation in this study, and the procedure is as follows:(1)
Given the initial crack depth and necessary parameters of the coupling model, the stiffness matrix of the cracked rotor is calculated using (4)–(8). Then, the Floquet multipliers under this parameter group are calculated according to Appendix B.(2)
Based on Floquet’s theory and the distribution of the Floquet multipliers on the complex plane, the vibration state of the cracked rotor is indentified. If the vibration state is stable or period-doubling bifurcation, go on to the next step. If the vibration state is the saddle-node bifurcation, it means degradation failure occurs and the computation stops. The crack depth at this time is the critical depth for unstable vibration (notated asau).(3)
The steady-state responses of the vibration differential equations are derived by the Runge–Kutta fourth order numerical integration method. Then, using (11)–(16), Ket of n rotation periods of time is calculated. The new crack depth after growth is derived by (19). If the vibration state is stable, (20) is used. If not, the rain-flow counting method is applied to get ΔKeMN.(4)
Judge whether the new crack depth reaches the given critical depth for rapid crack growth (notated asac). If so, it means degradation failure occurs and the computation stops. If not, go back to step 1, and continue the computation using the new crack depth until degradation failure occurs.In order to insure the veracity and reliability of the proposed model and simulation, the classical and validated models and theory are adopted, such as the harmonic model for the breathing effect, Paris’s law for the crack growth, Floquet’s theory for the stability analysis, and so on. Moreover, some intermediate results of the proposed simulation are compared to that of the published literatures, which is illustrated in Section6.Physical rotor parameters used for the simulations are as follows: shaft diameter,D = 0.045 m; length, L = 0.7 m; disk mass, m = 25 kg; and Young’s modulus, E = 2.1 × 105 MPa, referring the cracked rotor system in [16]. The empirical parameters used in the Paris–Erdogan laws are as follows: C = 6 × 108, B = 3, referring to the crack growth model in [39]. The other parameters including unbalance eccentricity ε and orientation angle β and damping factor ζ are analyzed in detail in the following sections.For different rotary shafts, physical rotor parameters such as shaft diameter, length, and mass are determined by the design scheme or can be directly measured in real applications. Young’s modulusE and the parameters B and C used in the Paris–Erdogan laws are material related, which can be chosen according to the prior knowledge about the shaft material. Besides, without any prior knowledge, they can be measured by the experiment [20, 34, 46]. The other parameters which are difficult to measure directly, such as unbalance eccentricity, damping factor, and initial crack size, can be estimated through the off-line experiment or on-line measured data by various parameter identification methods, such as the least squares method and the maximum likelihood method [47–49]. Since the parameter measurement and identification are not the key point in this study, they are not introduced in detail and interested readers can refer to the literatures listed above.It is assumed the time interval isn = 103 rotation periods, during which the crack depth is constant. Floquet multipliers are computed with the time interval K = 50 and the order of Taylor expansion Nl = 4. Given crack depth a, the crack depth ratio r=a/D. It is assumed the initial crack depth r0=0.1 (r0=a0/D), and the critical depth for rapid crack growth rc=0.5 (rc=ac/D).


### 5.2. Parametric Sensitivity Analysis

Parametric sensitivity analysis is important to the system model analysis. It helps to know which parameter is more sensitive to the system output. In this study, four parameters are of interest, namely, rotation speedΩ, unbalance eccentricity ε and orientation angle β, and damping factor ζ, and the concerned system output is fatigue degradation life TD. According to the research in [19], the concerned parameter ranges and corresponding unit are listed in Table 1.Table 1
The concerned parameters and their global sensitivities.
Ω
ε
β
ζ
Range
[1000,10000]
[0.5,2]
[0,π]
[0.01,0.05]
Unit
rpm
×10−5 m
rad
—
Sensitivity
1.56
0.98
0.46
0.05Given a parameter spaceHn of n dimensions, the local sensitivity of the system to a parameter pi (i=1,⁡2,⁡…,⁡n) is defined as the partial derivative of the system function with respect to the parameter. According to [50], the local sensitivity of pi to TD about a nominal point p∗ in the parameter space is calculated as(23)SiLp∗=TDp1∗,…,⁡pi−1∗,⁡pi∗+Δpi∗,⁡pi+1∗,…,⁡pn∗−TDp∗Δpi∗×pi∗TDp∗,where Δpi∗ is the small increment in the parameter of interest.Since the global sensitivity analysis accounts for interactions between variables and does not depend on the choice of a nominal point, it is adopted in this study. The grid resolution of the parameter space is ten grid points in each parameter axis, and 104 points are used to calculate global sensitivities. Based on a finite summation technique [50], the global sensitivity of pi to TD is calculated as(24)SiG=1m∑k=1mSiLpk,where pk represents the kth point in Hn and m = 104 is the total number of points in Hn. The global sensitivities of the interested parameters are shown in Table 1. The most sensitive parameter is rotation speed, which means the rotation speed is more influential to the degradation life than other parameters. The second sensitive parameter is unbalance eccentricity, which is more sensitive than the corresponding unbalance orientation angle. The least sensitive parameter is the damping factor. Notably, small global sensitivity does not mean a parameter is not important because its local sensitivity may be large about some possible points in the parameter space.Since global sensitivity only presents the overall and general influence of a parameter to the system, we also need to know the specific and detailed influence of multiple parameters on the degradation behavior. The degradation behavior of the cracked rotor is studied by using degradation measureD, fatigue degradation life TD, and average degradation rate VM. Different influences of multiple parameters (rotation speed, unbalance eccentricity and orientation angle, and damping) are investigated in detail in Section 6.


## 5.1. Procedure and Parameters of Numerical Simulation

The degradation behavior of the cracked rotor described by the proposed coupling model is addressed by numerical simulation in this study, and the procedure is as follows:(1)
Given the initial crack depth and necessary parameters of the coupling model, the stiffness matrix of the cracked rotor is calculated using (4)–(8). Then, the Floquet multipliers under this parameter group are calculated according to Appendix B.(2)
Based on Floquet’s theory and the distribution of the Floquet multipliers on the complex plane, the vibration state of the cracked rotor is indentified. If the vibration state is stable or period-doubling bifurcation, go on to the next step. If the vibration state is the saddle-node bifurcation, it means degradation failure occurs and the computation stops. The crack depth at this time is the critical depth for unstable vibration (notated asau).(3)
The steady-state responses of the vibration differential equations are derived by the Runge–Kutta fourth order numerical integration method. Then, using (11)–(16), Ket of n rotation periods of time is calculated. The new crack depth after growth is derived by (19). If the vibration state is stable, (20) is used. If not, the rain-flow counting method is applied to get ΔKeMN.(4)
Judge whether the new crack depth reaches the given critical depth for rapid crack growth (notated asac). If so, it means degradation failure occurs and the computation stops. If not, go back to step 1, and continue the computation using the new crack depth until degradation failure occurs.In order to insure the veracity and reliability of the proposed model and simulation, the classical and validated models and theory are adopted, such as the harmonic model for the breathing effect, Paris’s law for the crack growth, Floquet’s theory for the stability analysis, and so on. Moreover, some intermediate results of the proposed simulation are compared to that of the published literatures, which is illustrated in Section6.Physical rotor parameters used for the simulations are as follows: shaft diameter,D = 0.045 m; length, L = 0.7 m; disk mass, m = 25 kg; and Young’s modulus, E = 2.1 × 105 MPa, referring the cracked rotor system in [16]. The empirical parameters used in the Paris–Erdogan laws are as follows: C = 6 × 108, B = 3, referring to the crack growth model in [39]. The other parameters including unbalance eccentricity ε and orientation angle β and damping factor ζ are analyzed in detail in the following sections.For different rotary shafts, physical rotor parameters such as shaft diameter, length, and mass are determined by the design scheme or can be directly measured in real applications. Young’s modulusE and the parameters B and C used in the Paris–Erdogan laws are material related, which can be chosen according to the prior knowledge about the shaft material. Besides, without any prior knowledge, they can be measured by the experiment [20, 34, 46]. The other parameters which are difficult to measure directly, such as unbalance eccentricity, damping factor, and initial crack size, can be estimated through the off-line experiment or on-line measured data by various parameter identification methods, such as the least squares method and the maximum likelihood method [47–49]. Since the parameter measurement and identification are not the key point in this study, they are not introduced in detail and interested readers can refer to the literatures listed above.It is assumed the time interval isn = 103 rotation periods, during which the crack depth is constant. Floquet multipliers are computed with the time interval K = 50 and the order of Taylor expansion Nl = 4. Given crack depth a, the crack depth ratio r=a/D. It is assumed the initial crack depth r0=0.1 (r0=a0/D), and the critical depth for rapid crack growth rc=0.5 (rc=ac/D).


## 5.2. Parametric Sensitivity Analysis

Parametric sensitivity analysis is important to the system model analysis. It helps to know which parameter is more sensitive to the system output. In this study, four parameters are of interest, namely, rotation speedΩ, unbalance eccentricity ε and orientation angle β, and damping factor ζ, and the concerned system output is fatigue degradation life TD. According to the research in [19], the concerned parameter ranges and corresponding unit are listed in Table 1.Table 1
The concerned parameters and their global sensitivities.
Ω
ε
β
ζ
Range
[1000,10000]
[0.5,2]
[0,π]
[0.01,0.05]
Unit
rpm
×10−5 m
rad
—
Sensitivity
1.56
0.98
0.46
0.05Given a parameter spaceHn of n dimensions, the local sensitivity of the system to a parameter pi (i=1,⁡2,⁡…,⁡n) is defined as the partial derivative of the system function with respect to the parameter. According to [50], the local sensitivity of pi to TD about a nominal point p∗ in the parameter space is calculated as(23)SiLp∗=TDp1∗,…,⁡pi−1∗,⁡pi∗+Δpi∗,⁡pi+1∗,…,⁡pn∗−TDp∗Δpi∗×pi∗TDp∗,where Δpi∗ is the small increment in the parameter of interest.Since the global sensitivity analysis accounts for interactions between variables and does not depend on the choice of a nominal point, it is adopted in this study. The grid resolution of the parameter space is ten grid points in each parameter axis, and 104 points are used to calculate global sensitivities. Based on a finite summation technique [50], the global sensitivity of pi to TD is calculated as(24)SiG=1m∑k=1mSiLpk,where pk represents the kth point in Hn and m = 104 is the total number of points in Hn. The global sensitivities of the interested parameters are shown in Table 1. The most sensitive parameter is rotation speed, which means the rotation speed is more influential to the degradation life than other parameters. The second sensitive parameter is unbalance eccentricity, which is more sensitive than the corresponding unbalance orientation angle. The least sensitive parameter is the damping factor. Notably, small global sensitivity does not mean a parameter is not important because its local sensitivity may be large about some possible points in the parameter space.Since global sensitivity only presents the overall and general influence of a parameter to the system, we also need to know the specific and detailed influence of multiple parameters on the degradation behavior. The degradation behavior of the cracked rotor is studied by using degradation measureD, fatigue degradation life TD, and average degradation rate VM. Different influences of multiple parameters (rotation speed, unbalance eccentricity and orientation angle, and damping) are investigated in detail in Section 6.


## 6. Results and Discussion

Given a specific rotation speed, the degradation evolution of the cracked rotor is computed, which is characterized by the degradation measureD. As comparison, the degradation evolution of the cracked rotor model omitted the coupling of the vibration response, and crack growth is also presented by the degradation measure D. Without accounting for the coupling effect, the vibration response of the initial state is used to calculate the crack growth. In particular, given the unbalance eccentricity ε = 1 × 10−5 m, unbalance orientation angle β=0, and damping factor ζ = 0.01 (ζ=c/2k0m), the degradation evolutions for rotation speed Ω = 1000 rpm, 3000 rpm, 5000 rpm, and 7000 rpm are shown in Figures 5(a)–5(d), respectively, where the degradation time is presented as number of stress cycles.Figure 5
Degradation evolution presented by different degradation measures for different rotation speeds: (a) 1000 rpm; (b) 3000 rpm; (c) 5000 rpm; (d) 7000 rpm.
(a)
(b)
(c)
(d)As shown in Figures5(a), 5(b), and 5(d), the degradation measure D with the coupling effect considered reaches 1 for less number of cycles than D without coupling, which means the degradation is slower if the coupling effect is omitted. This phenomenon is in agreement with the result of [21]. As to a linear system, magnitude of vibration response increases gradually as stiffness of the system decreases, which leads to the acceleration of degradation. However, as to the nonlinear cracked rotor system, magnitude of vibration response may decrease as stiffness of the system decreases for some specific parameter group. The nonlinear breathing effect of the crack may destroy or weaken the mechanism by which rotational energy is transformed into transverse vibrational energy, which slows down the degradation, as shown in Figure 5(c). This phenomenon of the decrease in vibration amplitude as crack propagation has been reported in [51]. In engineering applications, underestimating the degradation rate may cause risky accidents, while overestimating the degradation rate can lead to a waste of maintenance resources. Therefore, when performing the fatigue design or using physics-based method to prognostic or predict residual life, the coupling effect cannot be omitted.


### 6.1. Effect of Rotation Speed

For the uncracked rotor, the first order critical rotation speed can be calculated using the stiffness of the uncracked rotor, asΩ0=60k0/m/2π=4645⁡rpm. The simulations in this study run over the speed range 600 rpm–10000 rpm, which covers the double of the critical speed (2Ω0=9290⁡rpm). Given the parameters ε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the degradation evolutions for different rotation speeds (1000 rpm–10000 rpm with increment 500 rpm) are shown in Figure 6, where the degradation evolutions are characterized by the degradation measure D.Figure 6
Degradation evolutions for various rotation speeds.
(a)
(b)
(c)
(d)As Figure6 shows, the degradation rate increases as the crack growth and the degradation evolution can be approximately depicted by the exponential function. For some rotation speeds, the degradation is obviously faster than others, like the degradation evolution for 1500 rpm, 4500 rpm, and 9000 rpm. However, there is no intuitive monotonous regularity for degradation evolutions as rotation speed increases. In fact, since the degradation of the cracked rotor may be influenced by multiple factors, such as unbalance eccentricity ε and orientation angle β and damping factor ζ, a simple regularity controlled by a single factor is definitely not adequate to illustrate the degradation behavior.Givenε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the fatigue degradation lifes TD for different rotation speeds are shown in Figure 7. TD does not change monotonously with speed but has several valleys. Specifically, when the rotation speed Ω is around 4600 rpm, TD is of very few cycles, which means degradation failure occurs very fast. 4600 rpm is almost the critical speed noting that the critical speed of the cracked rotor is a little less than that of the uncracked rotor due to the decrease of stiffness. Therefore, the resonance causes the large amplitude of the vibration response and furthermore the fast degradation. When the rotation speed is around 1100 rpm, 1500 rpm, and 2300 rpm, which are the submultiples (1/2, 1/3, and 1/4) of the critical speed, TD also decreases. The reason is that subharmonic resonance occurs [16] and causes the acceleration of degradation. This phenomenon of fast degradation gradually vanishes as speed decreases.Figure 7
Fatigue degradation life versus rotation speed withε = 1 × 10−5 m, β=1/2π, and ζ = 0.01.However, 2900 rpm (almost 2/3 of the critical speed) is far away from the resonance and subresonance region, around whichTD also has a small valley. The reason is that unstable vibration occurs at this parameter group if damping is very small [52, 53]. And the bifurcation causes the fast degradation, which will be verified in the following part.When the rotation speed is around 9100 rpm (almost the double of the critical speed),TD decreases and the initial degradation measure D is not 0.2 (Figure 6(d)). In this situation, degradation failure occurs owing to unstable vibration (saddle-node bifurcation) before the crack grows to the critical depth for rapid crack growth (rc=0.5). As shown in Figure 8, when the rotation speed is 9000 rpm, once the crack depth ratio reaches 0.32, the crack grows rapidly and failure occurs due to unstable vibration. The crack depth ratio 0.32 responds to the critical crack depth for unstable vibration (notated as au). Similarly, when the rotation speed is 8500 rpm, failure occurs at r=0.43.Figure 8
Crack growth evolutions for various rotation speeds.Givenε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the 3D and 2D color graphs of the average degradation rate VM versus crack depth ratio and rotation speed are presented in Figure 9. It is observed that the rotation speeds around which VM is relatively large are coincided with the valleys of TD shown in Figure 7. These rotation speeds are the critical speed (4600 rpm), double (9100 rpm) and submultiples (1100 rpm, 1500 rpm, and 2300 rpm) of the critical speed, and 2/3 of the critical speed (2900 rpm), respectively. The regions of relatively large VM around the critical speed and double of it are well coincided with the unstable vibration regions of the cracked rotor with small damping [27, 37, 47, 48]. However, it is worth noting that, around the critical speed, VM increases gradually while around double of the critical speed, VM increases abruptly when the crack depth reaches some threshold. The reason is that, around the critical speed, the cracked rotor reaches the unstable vibration region through period-doubling bifurcation, in which case the response amplitude increases gradually as the crack grows. However, around double of the critical speed, the cracked rotor reaches the unstable vibration region through saddle-node bifurcation, in which case the divergence of vibration response occurs abruptly when the crack grows to some threshold depth. The Floquet multipliers for 4000 rpm–4600 rpm and 8000 rpm–9100 rpm with the crack depth rate changing from 0 to 0.5 are plotted in Figures 10(a) and 10(b), respectively, which illustrate the two different types of bifurcation. Besides, around 2900 rpm and r=0.45, VM is also relatively large, which is responsible for the valley of TD at 2900 rpm in Figure 7. The bifurcation occurs at this parameter group as shown in Figure 11. The response amplitude of the y coordinate direction presents a larger variation around rotation speed 2900 rpm (Figure 11(a)) and around crack depth rate 0.45 (Figure 11(b)).Figure 9
Average degrading rate versus crack depth ratio and rotation speed withε⁡=⁡1⁡×⁡10−5⁡m, β=1/2π, and ζ⁡=⁡0.01.
(a)
(b)Figure 10
Floquet multipliers withr changing from 0 to 0.5 for (a) 4000 rpm–4600 rpm (period doubling) and (b) 8000 rpm–9100 rpm (saddle-node).
(a)
(b)Figure 11
Bifurcation diagram for (a)⁡Ω⁡=⁡2500 to 3500⁡rpm with r=0.45, β=1/2π, ε⁡=⁡1⁡×⁡10−5⁡m, and ζ⁡=⁡0.01 and (b) ⁡r⁡=⁡0 to 0.5 with Ω⁡=⁡2900⁡rpm, β=1/2π, ε⁡=⁡1⁡×⁡10−5⁡m, and ζ⁡=⁡0.01.
(a)
(b)


### 6.2. Effect of Unbalance

In engineering practice, unbalance eccentricity is inevitable due to manufacturing and assembly errors. The unbalance eccentricity and orientation angle are two crucial parameters for the cracked rotor, which are analyzed in the following two sections, respectively.


#### 6.2.1. Effect of Unbalance Eccentricity

The unbalance eccentricity is usually very small, but the effect on vibration behavior cannot be ignored [16, 29]. Then, the vibration behavior affects the crack growth, so it is significant to study the effect of unbalance eccentricity on degradation behavior. Given β=1/2π and ζ = 0.01, TD for various rotation speeds and different unbalance eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 12. For different unbalance eccentricities, the regularity between TD and rotation speed is almost the same. When the rotation speed is very slow (Ω < 1800 rpm), TD for different unbalance eccentricities are nearly identical. However, as rotation speed increases, the larger the unbalance eccentricity is, the shorter the TD is, owing to the larger amplitude of the vibration responses. Besides, around the critical speed and submultiples of it, the effect of unbalance eccentricity is relatively slight. However, around the rotation speeds far away from the resonance and subresonance regions, unbalance eccentricity has considerable influence on TD.Figure 12
Fatigue degradation life versus rotation speed for various eccentricities withβ=1/2π and ζ⁡=⁡0.01.Givenβ=1/2π and ζ = 0.01, the average degradation rates of the cracked rotor for different eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 13. The regions of relatively fast average degradation rate are almost the same for different eccentricities. As the eccentricity increases, the fast degradation region around the critical speed becomes larger, while that around submultiples and double of the critical speed barely changes. It is illustrated that eccentricity has little effect on the subresonance motion compared to the resonance. Around double of the critical speed, degradation failure occurs abruptly due to unstable vibration when crack depth reaches the threshold depth. Moreover, the judgment of this unstable vibration region has nothing to do with eccentricity with only the homogeneous form of the vibration differential equations considered according to Floquet’s theory.Figure 13
Average degrading rate versus crack depth ratio and rotation speed withβ=1/2π and ζ = 0.01 for various eccentricities: (a) ⁡ε⁡=⁡0.5⁡×⁡10−5⁡m; (b) ⁡ε⁡=⁡1⁡×⁡10−5⁡m; (c) ⁡ε⁡=⁡1.5⁡×⁡10−5⁡m; (d) ⁡ε⁡=⁡2⁡×⁡10−5⁡m.
(a)
(b)
(c)
(d)


#### 6.2.2. Effect of Unbalance Orientation Angle

The unbalance orientation angle is a crucial parameter having inescapable effect on vibration behavior of the cracked rotor, which has been studied by some researchers [16, 54–57]. Since the vibration is coupled with the crack growth, unbalance orientation angle affects the degradation behavior consequently. Given ε = 1 × 10−5 m and ζ = 0.01, TD for various rotation speeds and different unbalance orientation angles (0-π, with increment 1/4π) are shown in Figure 14, noting that the case for β∈π,⁡2π is omitted on the account of geometric symmetry. As shown in Figure 14, the regularity between TD and rotation speed is almost the same for different β. Around the critical speed and submultiples of it, the effect of β is relatively slight. However, around the other rotation speeds, β has a considerable effect on TD. It is worth noting that the speed range affected by the unbalance orientation angle is consistent with that affected by unbalance eccentricity.Figure 14
Fatigue degradation life versus rotation speed for various unbalance orientation angles withε = 1 × 10−5 m and ζ = 0.01.WhenΩ > 2300 rpm (1/2 of the critical speed), TD is larger for β=3/4π and is smaller for β=1/4π. However, when Ω < 2300 rpm, TD is larger for β=π. Take 1800 rpm, 3500 rpm, and 7000 rpm as representative rotation speeds to investigate the different effects of β for different speed ranges. Given r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the average degradation rates for different β are shown in Figure 15 and the amplitude of vibration response are shown in Figure 16. When Ω = 1800 rpm, VM decreases as β increases. However when Ω = 3500 rpm and 7000 rpm, VM decreases firstly and then increases as β increases. As to the vibration amplitude shown in Figure 16, the amplitude decreases as β increases for all the three speeds, and the amplitude is the smallest at β=π, which is in agreement with the result of [56, 57]. It is worth noting that when Ω = 3500 rpm and 7000 rpm, the VM increases as the amplitude decreases for relatively large β, which is different than anticipated. The reason is that different β causes different variation of the vibration phase [58] and then leads to the offset or superposition of the tensile or compressive stresses of the directions ξ and η. Therefore, the range of the equivalent stress may be inconsistent with the magnitude of the vibration amplitude. Given Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the variation of σξ, ση, and σe during one rotation period for β=3/4π and β=π is shown in Figure 17. The different effects of β=3/4π and β=π on the offset or superposition of the stresses are not obviously presented in Figure 17 due to the very small variation of the stress during only one rotation period. However, for relatively large numbers of rotation periods, the effect becomes remarkable due to accumulation. For example, when Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the crack grows 4.1168 × 10−4 m after 1000 rotation periods for β=3/4π and 4.7291 × 10−4 m for β=π.Figure 15
Average degradation rate versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 16
Vibration amplitude versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 17
Stress variation during a rotation period withΩ = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01: (a) ⁡β=3/4π; (b) ⁡β=π.
(a)
(b)


### 6.3. Effect of Damping

Damping is unavoidable in mechanical systems in practice, which is also an important parameter of the cracked rotor system. Givenε = 1 × 10−5 m and β=1/2π, TD for different damping ratios (0.01–0.05 with increment 0.01) and various rotation speeds are shown in Figure 18. For different damping ratios, the regularity of TD versus Ω is almost the same. The larger the damping ratio is, the larger the TD is. The influence of ζ on TD is mainly reflected around the submultiples and double of the critical speed. At 2/3 of the critical speed (2900 rpm), larger damping (ζ≥0.02) prevents the unstable bifurcation motion and TD does not decrease for large damping, which is in agreement with the result reported by Gasch [52]. It is worth noting that the rotation speed ranges affected by damping is quite different from that by unbalance. Therefore, in order to increase TD of a rotor running at a specific speed in practice, the effectiveness of reducing unbalance or increasing damping should be carefully considered. Besides, the valley of TD is vanished around double of the critical speed when damping is large (ζ = 0.05). The reason is that degradation failure caused by unstable vibration response does not occur since the large damping prevents the bifurcation motion of the cracked rotor system.Figure 18
Fatigue degradation life versus rotation speed for various damping withε = 1 × 10−5 m and β=1/2π.Givenε = 1 × 10−5 m and β=1/2π, the average degradation rates for different damping ratios (ζ = 0.02, 0.03, 0.04, and 0.05) are shown in Figure 19. It is revealed that damping has obvious influence on the region of relatively fast degradation rate. As damping increases, each region of relatively fast degradation rate shrinks and the value of VM decreases. Especially for large damping (ζ = 0.05), the region of the fast degradation rate completely vanishes around double of the critical speed. The reason is that damping ratio ζ has great influence on the stability of the cracked rotor system. Specifically, the area of the unstable region decreases as damping ratio increases, which means damping is beneficial for the stability of the cracked rotor system [27, 37, 52, 53]. Given ε = 1 × 10−5 m and β=1/2π, for different damping ratios (ζ = 0.01 and 0.05), the 3D graphs of the maximum Floquet multipliers with the change of crack depth and rotation speed are shown in Figure 20. Moreover, the Floquet multipliers are plotted on the complex plane as also shown in Figure 20. For small damping ratio (ζ = 0.01), the maximum Floquet multipliers exceed 1 around the critical speed and double of it with a deep crack (r>0.2). The way in which the Floquet multipliers exceed the unit circle as shown in Figure 20(a) indicates the types of bifurcation for unstable vibration responses are period doubling and saddle-node. For large damping ratio (ζ = 0.05), the maximum Floquet multipliers are always smaller than 1 as shown in Figure 20(b), which means the vibration responses are stable for any crack depth and rotation speed labeled on the image.Figure 19
Average degrading rate versus crack depth ratio and rotation speed withε = 1 × 10−5 m and β=1/2π for various damping ratios: (a) ⁡ζ⁡=⁡0.02; (b) ⁡ζ⁡=⁡0.03; (c) ⁡ζ⁡=⁡0.04; (d) ⁡ζ⁡=⁡0.05.
(a)
(b)
(c)
(d)Figure 20
Maximum Floquet multipliers versus crack depth ratio and rotation speed and Floquet multipliers on the complex plane for (a)⁡ζ⁡=⁡0.01 and (b) ⁡ζ⁡=⁡0.05.
(a)
(b)


## 6.1. Effect of Rotation Speed

For the uncracked rotor, the first order critical rotation speed can be calculated using the stiffness of the uncracked rotor, asΩ0=60k0/m/2π=4645⁡rpm. The simulations in this study run over the speed range 600 rpm–10000 rpm, which covers the double of the critical speed (2Ω0=9290⁡rpm). Given the parameters ε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the degradation evolutions for different rotation speeds (1000 rpm–10000 rpm with increment 500 rpm) are shown in Figure 6, where the degradation evolutions are characterized by the degradation measure D.Figure 6
Degradation evolutions for various rotation speeds.
(a)
(b)
(c)
(d)As Figure6 shows, the degradation rate increases as the crack growth and the degradation evolution can be approximately depicted by the exponential function. For some rotation speeds, the degradation is obviously faster than others, like the degradation evolution for 1500 rpm, 4500 rpm, and 9000 rpm. However, there is no intuitive monotonous regularity for degradation evolutions as rotation speed increases. In fact, since the degradation of the cracked rotor may be influenced by multiple factors, such as unbalance eccentricity ε and orientation angle β and damping factor ζ, a simple regularity controlled by a single factor is definitely not adequate to illustrate the degradation behavior.Givenε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the fatigue degradation lifes TD for different rotation speeds are shown in Figure 7. TD does not change monotonously with speed but has several valleys. Specifically, when the rotation speed Ω is around 4600 rpm, TD is of very few cycles, which means degradation failure occurs very fast. 4600 rpm is almost the critical speed noting that the critical speed of the cracked rotor is a little less than that of the uncracked rotor due to the decrease of stiffness. Therefore, the resonance causes the large amplitude of the vibration response and furthermore the fast degradation. When the rotation speed is around 1100 rpm, 1500 rpm, and 2300 rpm, which are the submultiples (1/2, 1/3, and 1/4) of the critical speed, TD also decreases. The reason is that subharmonic resonance occurs [16] and causes the acceleration of degradation. This phenomenon of fast degradation gradually vanishes as speed decreases.Figure 7
Fatigue degradation life versus rotation speed withε = 1 × 10−5 m, β=1/2π, and ζ = 0.01.However, 2900 rpm (almost 2/3 of the critical speed) is far away from the resonance and subresonance region, around whichTD also has a small valley. The reason is that unstable vibration occurs at this parameter group if damping is very small [52, 53]. And the bifurcation causes the fast degradation, which will be verified in the following part.When the rotation speed is around 9100 rpm (almost the double of the critical speed),TD decreases and the initial degradation measure D is not 0.2 (Figure 6(d)). In this situation, degradation failure occurs owing to unstable vibration (saddle-node bifurcation) before the crack grows to the critical depth for rapid crack growth (rc=0.5). As shown in Figure 8, when the rotation speed is 9000 rpm, once the crack depth ratio reaches 0.32, the crack grows rapidly and failure occurs due to unstable vibration. The crack depth ratio 0.32 responds to the critical crack depth for unstable vibration (notated as au). Similarly, when the rotation speed is 8500 rpm, failure occurs at r=0.43.Figure 8
Crack growth evolutions for various rotation speeds.Givenε = 1 × 10−5 m, β=1/2π, and ζ = 0.01, the 3D and 2D color graphs of the average degradation rate VM versus crack depth ratio and rotation speed are presented in Figure 9. It is observed that the rotation speeds around which VM is relatively large are coincided with the valleys of TD shown in Figure 7. These rotation speeds are the critical speed (4600 rpm), double (9100 rpm) and submultiples (1100 rpm, 1500 rpm, and 2300 rpm) of the critical speed, and 2/3 of the critical speed (2900 rpm), respectively. The regions of relatively large VM around the critical speed and double of it are well coincided with the unstable vibration regions of the cracked rotor with small damping [27, 37, 47, 48]. However, it is worth noting that, around the critical speed, VM increases gradually while around double of the critical speed, VM increases abruptly when the crack depth reaches some threshold. The reason is that, around the critical speed, the cracked rotor reaches the unstable vibration region through period-doubling bifurcation, in which case the response amplitude increases gradually as the crack grows. However, around double of the critical speed, the cracked rotor reaches the unstable vibration region through saddle-node bifurcation, in which case the divergence of vibration response occurs abruptly when the crack grows to some threshold depth. The Floquet multipliers for 4000 rpm–4600 rpm and 8000 rpm–9100 rpm with the crack depth rate changing from 0 to 0.5 are plotted in Figures 10(a) and 10(b), respectively, which illustrate the two different types of bifurcation. Besides, around 2900 rpm and r=0.45, VM is also relatively large, which is responsible for the valley of TD at 2900 rpm in Figure 7. The bifurcation occurs at this parameter group as shown in Figure 11. The response amplitude of the y coordinate direction presents a larger variation around rotation speed 2900 rpm (Figure 11(a)) and around crack depth rate 0.45 (Figure 11(b)).Figure 9
Average degrading rate versus crack depth ratio and rotation speed withε⁡=⁡1⁡×⁡10−5⁡m, β=1/2π, and ζ⁡=⁡0.01.
(a)
(b)Figure 10
Floquet multipliers withr changing from 0 to 0.5 for (a) 4000 rpm–4600 rpm (period doubling) and (b) 8000 rpm–9100 rpm (saddle-node).
(a)
(b)Figure 11
Bifurcation diagram for (a)⁡Ω⁡=⁡2500 to 3500⁡rpm with r=0.45, β=1/2π, ε⁡=⁡1⁡×⁡10−5⁡m, and ζ⁡=⁡0.01 and (b) ⁡r⁡=⁡0 to 0.5 with Ω⁡=⁡2900⁡rpm, β=1/2π, ε⁡=⁡1⁡×⁡10−5⁡m, and ζ⁡=⁡0.01.
(a)
(b)


## 6.2. Effect of Unbalance

In engineering practice, unbalance eccentricity is inevitable due to manufacturing and assembly errors. The unbalance eccentricity and orientation angle are two crucial parameters for the cracked rotor, which are analyzed in the following two sections, respectively.


### 6.2.1. Effect of Unbalance Eccentricity

The unbalance eccentricity is usually very small, but the effect on vibration behavior cannot be ignored [16, 29]. Then, the vibration behavior affects the crack growth, so it is significant to study the effect of unbalance eccentricity on degradation behavior. Given β=1/2π and ζ = 0.01, TD for various rotation speeds and different unbalance eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 12. For different unbalance eccentricities, the regularity between TD and rotation speed is almost the same. When the rotation speed is very slow (Ω < 1800 rpm), TD for different unbalance eccentricities are nearly identical. However, as rotation speed increases, the larger the unbalance eccentricity is, the shorter the TD is, owing to the larger amplitude of the vibration responses. Besides, around the critical speed and submultiples of it, the effect of unbalance eccentricity is relatively slight. However, around the rotation speeds far away from the resonance and subresonance regions, unbalance eccentricity has considerable influence on TD.Figure 12
Fatigue degradation life versus rotation speed for various eccentricities withβ=1/2π and ζ⁡=⁡0.01.Givenβ=1/2π and ζ = 0.01, the average degradation rates of the cracked rotor for different eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 13. The regions of relatively fast average degradation rate are almost the same for different eccentricities. As the eccentricity increases, the fast degradation region around the critical speed becomes larger, while that around submultiples and double of the critical speed barely changes. It is illustrated that eccentricity has little effect on the subresonance motion compared to the resonance. Around double of the critical speed, degradation failure occurs abruptly due to unstable vibration when crack depth reaches the threshold depth. Moreover, the judgment of this unstable vibration region has nothing to do with eccentricity with only the homogeneous form of the vibration differential equations considered according to Floquet’s theory.Figure 13
Average degrading rate versus crack depth ratio and rotation speed withβ=1/2π and ζ = 0.01 for various eccentricities: (a) ⁡ε⁡=⁡0.5⁡×⁡10−5⁡m; (b) ⁡ε⁡=⁡1⁡×⁡10−5⁡m; (c) ⁡ε⁡=⁡1.5⁡×⁡10−5⁡m; (d) ⁡ε⁡=⁡2⁡×⁡10−5⁡m.
(a)
(b)
(c)
(d)


### 6.2.2. Effect of Unbalance Orientation Angle

The unbalance orientation angle is a crucial parameter having inescapable effect on vibration behavior of the cracked rotor, which has been studied by some researchers [16, 54–57]. Since the vibration is coupled with the crack growth, unbalance orientation angle affects the degradation behavior consequently. Given ε = 1 × 10−5 m and ζ = 0.01, TD for various rotation speeds and different unbalance orientation angles (0-π, with increment 1/4π) are shown in Figure 14, noting that the case for β∈π,⁡2π is omitted on the account of geometric symmetry. As shown in Figure 14, the regularity between TD and rotation speed is almost the same for different β. Around the critical speed and submultiples of it, the effect of β is relatively slight. However, around the other rotation speeds, β has a considerable effect on TD. It is worth noting that the speed range affected by the unbalance orientation angle is consistent with that affected by unbalance eccentricity.Figure 14
Fatigue degradation life versus rotation speed for various unbalance orientation angles withε = 1 × 10−5 m and ζ = 0.01.WhenΩ > 2300 rpm (1/2 of the critical speed), TD is larger for β=3/4π and is smaller for β=1/4π. However, when Ω < 2300 rpm, TD is larger for β=π. Take 1800 rpm, 3500 rpm, and 7000 rpm as representative rotation speeds to investigate the different effects of β for different speed ranges. Given r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the average degradation rates for different β are shown in Figure 15 and the amplitude of vibration response are shown in Figure 16. When Ω = 1800 rpm, VM decreases as β increases. However when Ω = 3500 rpm and 7000 rpm, VM decreases firstly and then increases as β increases. As to the vibration amplitude shown in Figure 16, the amplitude decreases as β increases for all the three speeds, and the amplitude is the smallest at β=π, which is in agreement with the result of [56, 57]. It is worth noting that when Ω = 3500 rpm and 7000 rpm, the VM increases as the amplitude decreases for relatively large β, which is different than anticipated. The reason is that different β causes different variation of the vibration phase [58] and then leads to the offset or superposition of the tensile or compressive stresses of the directions ξ and η. Therefore, the range of the equivalent stress may be inconsistent with the magnitude of the vibration amplitude. Given Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the variation of σξ, ση, and σe during one rotation period for β=3/4π and β=π is shown in Figure 17. The different effects of β=3/4π and β=π on the offset or superposition of the stresses are not obviously presented in Figure 17 due to the very small variation of the stress during only one rotation period. However, for relatively large numbers of rotation periods, the effect becomes remarkable due to accumulation. For example, when Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the crack grows 4.1168 × 10−4 m after 1000 rotation periods for β=3/4π and 4.7291 × 10−4 m for β=π.Figure 15
Average degradation rate versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 16
Vibration amplitude versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 17
Stress variation during a rotation period withΩ = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01: (a) ⁡β=3/4π; (b) ⁡β=π.
(a)
(b)


## 6.2.1. Effect of Unbalance Eccentricity

The unbalance eccentricity is usually very small, but the effect on vibration behavior cannot be ignored [16, 29]. Then, the vibration behavior affects the crack growth, so it is significant to study the effect of unbalance eccentricity on degradation behavior. Given β=1/2π and ζ = 0.01, TD for various rotation speeds and different unbalance eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 12. For different unbalance eccentricities, the regularity between TD and rotation speed is almost the same. When the rotation speed is very slow (Ω < 1800 rpm), TD for different unbalance eccentricities are nearly identical. However, as rotation speed increases, the larger the unbalance eccentricity is, the shorter the TD is, owing to the larger amplitude of the vibration responses. Besides, around the critical speed and submultiples of it, the effect of unbalance eccentricity is relatively slight. However, around the rotation speeds far away from the resonance and subresonance regions, unbalance eccentricity has considerable influence on TD.Figure 12
Fatigue degradation life versus rotation speed for various eccentricities withβ=1/2π and ζ⁡=⁡0.01.Givenβ=1/2π and ζ = 0.01, the average degradation rates of the cracked rotor for different eccentricities (ε = 0.5 × 10−5 m, 1 × 10−5 m, 1.5 × 10−5 m, and 2 × 10−5 m) are shown in Figure 13. The regions of relatively fast average degradation rate are almost the same for different eccentricities. As the eccentricity increases, the fast degradation region around the critical speed becomes larger, while that around submultiples and double of the critical speed barely changes. It is illustrated that eccentricity has little effect on the subresonance motion compared to the resonance. Around double of the critical speed, degradation failure occurs abruptly due to unstable vibration when crack depth reaches the threshold depth. Moreover, the judgment of this unstable vibration region has nothing to do with eccentricity with only the homogeneous form of the vibration differential equations considered according to Floquet’s theory.Figure 13
Average degrading rate versus crack depth ratio and rotation speed withβ=1/2π and ζ = 0.01 for various eccentricities: (a) ⁡ε⁡=⁡0.5⁡×⁡10−5⁡m; (b) ⁡ε⁡=⁡1⁡×⁡10−5⁡m; (c) ⁡ε⁡=⁡1.5⁡×⁡10−5⁡m; (d) ⁡ε⁡=⁡2⁡×⁡10−5⁡m.
(a)
(b)
(c)
(d)


## 6.2.2. Effect of Unbalance Orientation Angle

The unbalance orientation angle is a crucial parameter having inescapable effect on vibration behavior of the cracked rotor, which has been studied by some researchers [16, 54–57]. Since the vibration is coupled with the crack growth, unbalance orientation angle affects the degradation behavior consequently. Given ε = 1 × 10−5 m and ζ = 0.01, TD for various rotation speeds and different unbalance orientation angles (0-π, with increment 1/4π) are shown in Figure 14, noting that the case for β∈π,⁡2π is omitted on the account of geometric symmetry. As shown in Figure 14, the regularity between TD and rotation speed is almost the same for different β. Around the critical speed and submultiples of it, the effect of β is relatively slight. However, around the other rotation speeds, β has a considerable effect on TD. It is worth noting that the speed range affected by the unbalance orientation angle is consistent with that affected by unbalance eccentricity.Figure 14
Fatigue degradation life versus rotation speed for various unbalance orientation angles withε = 1 × 10−5 m and ζ = 0.01.WhenΩ > 2300 rpm (1/2 of the critical speed), TD is larger for β=3/4π and is smaller for β=1/4π. However, when Ω < 2300 rpm, TD is larger for β=π. Take 1800 rpm, 3500 rpm, and 7000 rpm as representative rotation speeds to investigate the different effects of β for different speed ranges. Given r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the average degradation rates for different β are shown in Figure 15 and the amplitude of vibration response are shown in Figure 16. When Ω = 1800 rpm, VM decreases as β increases. However when Ω = 3500 rpm and 7000 rpm, VM decreases firstly and then increases as β increases. As to the vibration amplitude shown in Figure 16, the amplitude decreases as β increases for all the three speeds, and the amplitude is the smallest at β=π, which is in agreement with the result of [56, 57]. It is worth noting that when Ω = 3500 rpm and 7000 rpm, the VM increases as the amplitude decreases for relatively large β, which is different than anticipated. The reason is that different β causes different variation of the vibration phase [58] and then leads to the offset or superposition of the tensile or compressive stresses of the directions ξ and η. Therefore, the range of the equivalent stress may be inconsistent with the magnitude of the vibration amplitude. Given Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the variation of σξ, ση, and σe during one rotation period for β=3/4π and β=π is shown in Figure 17. The different effects of β=3/4π and β=π on the offset or superposition of the stresses are not obviously presented in Figure 17 due to the very small variation of the stress during only one rotation period. However, for relatively large numbers of rotation periods, the effect becomes remarkable due to accumulation. For example, when Ω = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01, the crack grows 4.1168 × 10−4 m after 1000 rotation periods for β=3/4π and 4.7291 × 10−4 m for β=π.Figure 15
Average degradation rate versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 16
Vibration amplitude versus unbalance orientation angle for different rotation speeds withr=0.4, ε = 1 × 10−5 m, and ζ = 0.01.Figure 17
Stress variation during a rotation period withΩ = 7000 rpm, r=0.4, ε = 1 × 10−5 m, and ζ = 0.01: (a) ⁡β=3/4π; (b) ⁡β=π.
(a)
(b)


## 6.3. Effect of Damping

Damping is unavoidable in mechanical systems in practice, which is also an important parameter of the cracked rotor system. Givenε = 1 × 10−5 m and β=1/2π, TD for different damping ratios (0.01–0.05 with increment 0.01) and various rotation speeds are shown in Figure 18. For different damping ratios, the regularity of TD versus Ω is almost the same. The larger the damping ratio is, the larger the TD is. The influence of ζ on TD is mainly reflected around the submultiples and double of the critical speed. At 2/3 of the critical speed (2900 rpm), larger damping (ζ≥0.02) prevents the unstable bifurcation motion and TD does not decrease for large damping, which is in agreement with the result reported by Gasch [52]. It is worth noting that the rotation speed ranges affected by damping is quite different from that by unbalance. Therefore, in order to increase TD of a rotor running at a specific speed in practice, the effectiveness of reducing unbalance or increasing damping should be carefully considered. Besides, the valley of TD is vanished around double of the critical speed when damping is large (ζ = 0.05). The reason is that degradation failure caused by unstable vibration response does not occur since the large damping prevents the bifurcation motion of the cracked rotor system.Figure 18
Fatigue degradation life versus rotation speed for various damping withε = 1 × 10−5 m and β=1/2π.Givenε = 1 × 10−5 m and β=1/2π, the average degradation rates for different damping ratios (ζ = 0.02, 0.03, 0.04, and 0.05) are shown in Figure 19. It is revealed that damping has obvious influence on the region of relatively fast degradation rate. As damping increases, each region of relatively fast degradation rate shrinks and the value of VM decreases. Especially for large damping (ζ = 0.05), the region of the fast degradation rate completely vanishes around double of the critical speed. The reason is that damping ratio ζ has great influence on the stability of the cracked rotor system. Specifically, the area of the unstable region decreases as damping ratio increases, which means damping is beneficial for the stability of the cracked rotor system [27, 37, 52, 53]. Given ε = 1 × 10−5 m and β=1/2π, for different damping ratios (ζ = 0.01 and 0.05), the 3D graphs of the maximum Floquet multipliers with the change of crack depth and rotation speed are shown in Figure 20. Moreover, the Floquet multipliers are plotted on the complex plane as also shown in Figure 20. For small damping ratio (ζ = 0.01), the maximum Floquet multipliers exceed 1 around the critical speed and double of it with a deep crack (r>0.2). The way in which the Floquet multipliers exceed the unit circle as shown in Figure 20(a) indicates the types of bifurcation for unstable vibration responses are period doubling and saddle-node. For large damping ratio (ζ = 0.05), the maximum Floquet multipliers are always smaller than 1 as shown in Figure 20(b), which means the vibration responses are stable for any crack depth and rotation speed labeled on the image.Figure 19
Average degrading rate versus crack depth ratio and rotation speed withε = 1 × 10−5 m and β=1/2π for various damping ratios: (a) ⁡ζ⁡=⁡0.02; (b) ⁡ζ⁡=⁡0.03; (c) ⁡ζ⁡=⁡0.04; (d) ⁡ζ⁡=⁡0.05.
(a)
(b)
(c)
(d)Figure 20
Maximum Floquet multipliers versus crack depth ratio and rotation speed and Floquet multipliers on the complex plane for (a)⁡ζ⁡=⁡0.01 and (b) ⁡ζ⁡=⁡0.05.
(a)
(b)


## 7. Conclusions

The fast-slow coupling model of the rotor vibration and the crack growth is established, and an approximate numerical procedure based on sequential iteration is proposed to solve the coupling model. A competing degradation failure criterion is proposed through the analysis of the two degradation failure modes of the cracked rotor system, namely, crack rapid growth failure and unstable vibration failure. Numerical simulations show that if the coupling effect is omitted, one may obtain wrong degradation evolutions, which may cause risky accidents or waste of maintenance resources.The effects of multiple parameters on the degradation behavior, including rotation speed, unbalance eccentricity and orientation angle, and damping, are investigated in detail. Results show that there is no simple intuitive regularity of degradation evolution, since the degradation behavior is affected by multiple parameters of the coupling model. When the rotation speed is around the critical speed, the vibration becomes unstable through period-doubling bifurcation, as the crack propagates. Due to the increase of vibration amplitude, the cracked rotor degrades very fast. When the rotation speed is around submultiples of the critical speed, degradation is relatively fast due to subresonance. And this phenomenon of fast degradation gradually vanishes as speed decreases. When the rotation speed is far away from the resonance and subresonance region, the degradation may also be fast due to the unstable vibration responses under some parameter group, for example, 2/3 of the critical speed with relatively deep crack and small damping. When the rotation speed is around double of the critical speed, the vibration may become unstable abruptly due to saddle-node bifurcation. Once the crack grows to some threshold depth, degradation failure occurs abruptly.The degradation rate increases as unbalance eccentricity increases and as damping decreases. Unbalance orientation angle causes different variations of the vibration phase, which leads to the offset or superposition of the tensile or compressive stresses. Thus, when rotation speed is larger than 1/2 of the critical speed, degradation is slow forβ=3/4π; otherwise, degradation is slow for β=π. Moreover, the rotation speed ranges affected by damping are mainly the speed around the critical speed, submultiples, and double of it but that affected by unbalance eccentricity and orientation angle are just the opposite.In engineering practice, the speed ranges corresponding to the fast fatigue degradation should be avoided and reducing unbalance eccentricity, changing unbalance orientation, or increasing damping can effectively retard degradation at specific rotation speeds. The findings are of significance to guide the safety design of the rotor system for long time operation and help to the further research on prognostics and lifetime prediction.



---

*Source: 3239523-2018-08-29.xml*