# The global responses characteristics of a rotor/stator rubbing system with dry friction effects

Zhiyong Shang, Jun Jiang, Ling Hong

► To cite this version:

Zhiyong Shang, Jun Jiang, Ling Hong. The global responses characteristics of a rotor/stator rubbing system with dry friction effects. Journal of Sound and Vibration, 2011, 330 (10), pp.2150-2160. 10.1016/j.jsv.2010.06.004 . hal-01537669

HAL Id: hal-01537669

https://hal.science/hal-01537669

Submitted on 12 Jun 2017

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The global responses characteristics of a rotor/stator rubbing system with dry friction effects

Zhiyong Shang, Jun Jiang, Ling Hong

MOE Key Laboratory of Strength and Vibration, School of Aerospace Engineering, Xi'an Jiaotong University, 710049 Xi'an, China

Rotor/stator rubbing systems may undertake a number of quite different responses. Recent experiments on rotor/stator rubbing have revealed that two or three different responses may coexist. In this paper the global response characteristics of a general rotor/stator rubbing system, which takes into account the dominant factors in the process of rotor/stator rubbing, especially, the dry friction effect that is mostly neglected in the previous works and is the main factor for the self-excited dry friction backward whirl, are studied. The different solutions of the piecewise nonlinear system are derived and their stability are analyzed to get the existence boundaries of the different responses. An overall picture of the global response characteristics of this model is then obtained by drawing the existence boundaries in a same parameter space. The present results provide good understanding on the coexistence of different rubbing responses observed in tests. Moreover, deeper insight into the types of coexistence of different rubbing responses and their relationship with the system parameters is gained.

## 1. Introduction

With the need to enhance the efficiency of rotating machinery, the gap between the rotor and the stator of a rotating machine has been significantly reduced; meanwhile the rotor/stator rubbing is becoming not uncommon. Rotor/stator rubbing is a serious malfunction during the operation of a rotating machine, which can seriously degrade the machine performance and can even lead to disastrous consequences of the machine. There are a large amount of works on the rub-related phenomena in the rotor/stator systems in order to get deep insights into the dynamical behaviors and their relationship with system parameters. Some papers studied the possible types of rubbing responses induced through the rotor/stator rubbing and tried to derive the relationship between the rubbing responses and the system parameters, as summarized in [1]. Through the numerical and experimental studies, many rubbing behavior in rotor/stator contact systems have now become well known, for instance, rubbing can induce periodic synchronous full annular rubs [2-4], sub-and super-harmonic motions [5,6], the partial rubs in quasi-periodic whirl [7], the chaotic motion [8,9], as well as the destructive self-excited dry friction backward whirl [10-13]. The recent work aiming at the inclusion of more realistic effects in the rotor/stator rubbing systems was presented in [14] where the impact damping as well as the offset between the geometrical centers of the rotor and the stator are taken into account.

For the piecewise nonlinear systems of rotor/stator rubbing, the coexistence of multiple responses for a given set of parameters should be not seldom. A common nonlinear phenomenon — jump phenomenon — was observed
experimentally in a rotor/stator system around the mid-1950s last century by Yamamoto [15]. This phenomenon indicates the coexistence of the synchronous no-rub response and the synchronous full annular rub response. It was observed in tests [4] that when an additional external disturbance (e.g., hitting on the shaft) was applied on the rotor which is undertaking a normal operation at a very low rotating speed, dry friction backward whirl would be triggered. This indicates that the synchronous no-rub response coexists with the dry friction backward whirl. In Choi's experiments [18], it was observed that there are jump phenomena between two backward whirl frequencies during run-up and run-down of the rotor. This implies that two dry friction backward whirl responses with different whirl frequencies may coexist. In fact, there are no theoretical models in the literature that has successfully predicted this phenomenon. The results in this paper will well demonstrate the appearance of such phenomenon and also show the dependence of the phenomenon on the system parameters. Most recently, Ehehalt et al. have experimentally detected the coexistence of the no-rub motion, stable forward whirl and stable backward whirl within some speed ranges in their tests [18]. From above discussion it is found that different types of coexistence exist in the rotor/stator rubbing systems that may cause, under some disturbance, quite different response cascades with the change of the rotating speed. So it is of great necessity to reveal the global response characteristics of the rotor/stator rubbing systems in order to well understand the response behavior of the rotor/stator rubbing systems.

This paper aims to depict the global response characteristics of a general rotor/stator rubbing system, which takes into account both the dynamics of the rotor and the stator as well as the friction and the deformation on the contact surfaces. Most importantly, the dry friction effect, which is mostly neglected in the published works on the rotor/stator rubbing, is also considered in the present model in order to grasp the self-excited dry friction backward whirl. In obtaining different solutions of the piecewise nonlinear system analytically, different techniques are adopted in this paper. For instance, to get the synchronous full annular rub solutions analytically, the resultant contact force is expressed in a form which is different from the usual way to use direction of displacement to represent the direction of contact force. By deriving the solution for the self-excited dry friction backward whirl, the dry friction effects are properly dealt with based on the experimental and numerical observation that the amplitude for the backward whirl component is much larger than that of the force forward whirl component. To determine the boundaries between the quasi-periodic partial rub and the dry friction backward whirl, we start from the physical reason for the loss of stability of the partial rub, namely, the resonance at the nonlinear natural frequencies.

The paper is organized in the following ways. In Section 2, the model of the rotor/stator system is introduced. Then, the boundary of the synchronous no-rub motion and that of the synchronous full annular rub are analyzed in Section 3. The boundary between the partial rub and the dry friction backward whirl are predicted in Section 4. In Section 5, the existence boundaries of the backward whirl are determined analytically by using a proper criterion for the start of the solution. Section 6 illustrates the global response characteristics of rotor/stator system in the parameter space through examples. Finally, Section 7 is the conclusions of this work.

## 2. Mathematical model

The rotor/stator model studied in this paper is shown in Fig. 1. A weightless shaft is supported in ideal bearings having the effective transverse stiffness $k_r$ and rotates at an angular speed $\omega$. A rigid disk of mass $m_r$ is mounted at the midpoint of the shaft and possesses a mass eccentricity of $e$. Concentric with the disk is an annular stator of mass $m_s$. The stator is elastically supported by a symmetrical set of springs with isotropic radial stiffness $k_s$. There exists a clearance $\delta$ between the rotor and the stator. During rotor/stator rubbing, a friction force with friction coefficient $\mu$ applies tangentially at the contact point. The pointing direction of the friction force (forward or backward) depends on the direction of the relative

<figure><img src="image_2.png" /><figcaption>Fig. 1. The model of rotor/stator rubbing system.</figcaption></figure>
velocity between the circumferential speed and the whirl speed of the rotor with respect to the stator at the contact point. To consider the deformation at the contact surfaces, a symmetrical set of fictive springs with isotropic radial stiffness $k_c$ is assumed being laid in the inner ring of the stator to model the contact stiffness.

The equations that govern the motion of the rotor/stator system can be written in the complex form as

$$ 
\begin{aligned}
& m_r \ddot{r}_r + c_r \dot{r}_r + k_r r_r + \boldsymbol{\Theta}(1 + j\mu \operatorname{sign}(v_{\mathrm{rel}})) k_c \quad r_r - r_s - \delta \frac{r_r - r_s}{|r_r - r_s|} = m_r e^{j\omega^2 t} \\
& m_s \ddot{r}_s + c_s \dot{r}_s + k_s r_s - \boldsymbol{\Theta}(1 + j\mu \operatorname{sign}(v_{\mathrm{rel}})) k_c \quad r_r - r_s - \delta \frac{r_r - r_s}{|r_r - r_s|} = 0 \\
& v_{\mathrm{rel}} = r_{\mathrm{disk}} \omega + |r_r - r_s| \omega_w
\end{aligned}
\quad (1) $$

where $r_r=y_r+jz_r$ and $r_s=y_s+jz_s$ are, respectively, the complex deflections of the rotor and the stator, and $c_r$ and $c_s$ the damping constant. $v_{\mathrm{rel}}$ is the relative velocity at the contact point with $\omega_w$ being the whirl frequency of and $r_{\mathrm{disk}}$ the radius of the rotor. Additionally, we define $\Theta=1$, if $|r_r-r_s| \ge \delta$ and $\Theta=0$, if $|r_r-r_s| < \delta$.

Eq. (1) can be formulated into the non-dimensional form as

$$ 
\begin{aligned}
& \hat{r}'' + 2\zeta_r \hat{r}' + \hat{r}_r + \boldsymbol{\Theta} \beta_{cr} (1 + j\mu \operatorname{sign}(V_{\mathrm{rel}})) \left( 1 - \frac{\Delta}{|\hat{r}_r - \hat{r}_s|} \right) (\hat{r}_r - \hat{r}_s) = \boldsymbol{\Omega}^2 e^{j\Omega\tau} \\
& M_{sr} \hat{r}'' + 2\zeta_s \sqrt{M_{sr} \beta_{sr} \hat{r}' + \beta_{sr} \hat{r}_s - \boldsymbol{\Theta} \beta_{cr} (1 + j\mu \operatorname{sign}(V_{\mathrm{rel}}))} \left( 1 - \frac{\Delta}{|\hat{r}_r - \hat{r}_s|} \right) (\hat{r}_r - \hat{r}_s) = 0 \\
& V_{\mathrm{rel}} = R_{\mathrm{disk}} \boldsymbol{\Omega} + |\hat{r}_r - \hat{r}_s| \omega_b
\end{aligned}
\quad (2) $$

where $'$ represents the differentiation with respect to the non-dimensional time $\tau = \omega_0 t$ with $\omega_0 = \sqrt{k_r/m_r}$ being the natural frequency of the rotor. The other non-dimensional variables are defined as follows:

$$ 
\begin{aligned}
\hat{r}_r &= \frac{r_r}{e}, \quad \hat{r}_s = \frac{r_s}{e}, \quad \hat{F} = \frac{F}{ek_r}, \quad M_{sr} = \frac{m_s}{m_r}, \quad \beta_{sr} = \frac{k_s}{k_r}, \quad \beta_{cr} = \frac{k_c}{k_r}, \quad A = \frac{\delta}{e}, \\
\boldsymbol{\Omega} &= \frac{\omega}{\omega_0}, \quad \zeta_r = \frac{c_r}{2\sqrt{k_r m_r}}, \quad \zeta_s = \frac{c_s}{2\sqrt{k_s m_s}}, \quad R_{\mathrm{disk}} = \frac{r_{\mathrm{disk}}}{e}, \quad \omega_b = \frac{\omega_w}{\omega_0}, \quad \bar{\zeta}_s = \zeta_s \sqrt{M_{sr} \beta_{sr}}
\end{aligned}
$$

The dynamics of rotor/stator system can be divided into different phases of motion which are governed by different equations depending upon the motions of the rotor and stator deflection. For instance, when the relative deflection between the rotor and the stator is less than the clearance, the rotor and the stator are not in contact. In this case, $\Theta$ in Eq. (1) or (2) will be equal to 0. When the rotor and the stator are in full annular contact, namely, $\Theta=1$, there are still two cases that should be distinguished. Depending upon whether the friction changes direction or not, the governing equation with $\Theta=1$ may have two different forms. As will be shown in this paper, the most often observed responses in rotor/stator rubbing systems, such as the synchronous no-rub motion, the synchronous full annular rub motion and the dry friction backward whirl, can be determined from different forms of the governing equations. Although the quasi-periodic partial rub motion alternates between the phases without rubbing and with rubbing, its existence boundaries can be determined indirectly by analyzing the equation with rubbing and the finding in [13] that the loss of the stability of partial rub is due to the resonance at the natural frequencies of the nonlinear normal modes. With all the existence boundaries of the above mentioned responses being determined, the global response characteristics of the rotor/stator system in parameter space can be depicted.

### 3. Boundaries for the synchronous responses

For the rotor/stator system, there are two types of synchronous responses: one is the normal operation condition by which the rotor whirls at the same frequency as the frequency it spins and there is no contact between the rotor and stator; the other is the mild rubbing condition by which the rotor touches the stator all the time and whirls synchronously. Below the boundaries of these two responses will be determined.

#### 3.1. The boundary of no-rub synchronous motion

When the relative deflection between the rotor and the stator is less than the clearance, there is no rubbing in the system, in this case $\Theta$ in Eq. (2) will be equal to 0:

$$ 
\begin{aligned}
& \hat{r}'' + 2\zeta_r \hat{r}' + \hat{r}_r = \boldsymbol{\Omega}^2 e^{j\Omega\tau} \\
& M_{sr} \hat{r}'' + 2\zeta_s \sqrt{M_{sr} \beta_{sr} \hat{r}' + \beta_{sr} \hat{r}_s} = 0
\end{aligned}
$$
The synchronous no-rub motion is the only response and the periodic solution can be written as

$$ \dot{r}_r = A e^{j(\Omega\tau + \phi)} \quad (3) $$

$$ \text{where } A = \frac{\Omega^2}{\sqrt{(1-\Omega^2)^2 + (2\zeta_r\Omega)^2}}, \quad \text{tg}\phi = -\frac{2\zeta_r\Omega}{1-\Omega^2} $$

This period solution is unconditionally stable. The solution is physically valid only when the amplitude of the rotor is smaller than the clearance, namely $A < \Delta$. In order to find the parameter regions where the synchronous no-rub motion exists, let $A=\Delta$ in (3). The parameter equation defining the boundaries of no-rub motion is

$$ (\Delta^2 - 1)\Omega^4 + \Delta^2(4\zeta_r^2 - 2)\Omega^2 + \Delta^2 = 0 \quad (4) $$

This equation has two pairs of real roots, each with an equal magnitude but opposite sign. Since the rotating speed is supposed to be positive and the synchronous no-rub response always undertakes forward whirl, the two positive roots, denoted as $\Omega_L$ and $\Omega_U$, define the boundaries of the synchronous no-rub motion. The synchronous no-rub motion exists when the rotating speed is less than $\Omega_L$ or greater than $\Omega_U$ (see Figs.2 and 3).

## 3.2. Boundaries of synchronous full annular rub

To solve the explicit form of the synchronous full annular rub solution from the nonlinear equation of (2) with $\Theta=1$, an equivalent form of the resultant contact force should be adopted as done in [21], which is in the form

$$ \hat{F} = \beta_{cr}(1+j\mu) \quad \hat{r}_r - \hat{r}_s - \Delta \frac{1-j\mu}{\sqrt{1+\mu^2}} \frac{\hat{F}}{|\hat{F}|} \quad (5) $$

Suppose the synchronous full annular rub solution has the following form

$$ \hat{r}_r = \rho_r e^{j\Omega t}, \quad \hat{r}_s = \rho_s e^{j\Omega t}, \quad \hat{F} = \Phi e^{j\Omega t} \quad (6) $$

where $\rho_r$, $\rho_s$ and $\Phi$ are the complex constants because this solution, by its definition, has a constant amplitude of deflection and a whirl speed equal to $\Omega$. After substituting (6) into Eq. (2) with the terms representing contact force being replaced by (5), the equations governing the complex amplitudes of the rotor, the stator and the contact force are given as

$$ C_r \rho_r = \Omega^2 - \Phi $$

$$ C_s \rho_s = \Phi $$

$$ (1+j\mu)\beta_{cr} \quad \rho_r - \rho_s - \Delta \frac{1-j\mu}{\sqrt{1+\mu^2}} \frac{\Phi}{|\Phi|} = \Phi \quad (7) $$

where $C_r$ and $C_s$ denote the dynamical stiffness of the rotor and the stator

$$ C_r = 1 - \Omega^2 + \Omega_{\gamma r}^2 + j(2\Omega_{\gamma r}^2 - \kappa_r), \quad C_s = \beta_{sr} - \Omega^2 M_{sr} + \Omega_{\gamma s}^2 + j(2\Omega_{\gamma s}^2 - \kappa_s) $$

The synchronous full annular rub solutions are in the form

$$ \hat{r}_{r0} = \rho_{r0} e^{j\Omega\tau} = B_r e^{j(\Omega\tau + \alpha)} = y_{r0}(\tau) + jz_{r0}(\tau) $$

$$ \hat{r}_{s0} = \rho_{s0} e^{j\Omega\tau} = B_s e^{j(\Omega\tau + \beta)} = y_{s0}(\tau) + jz_{s0}(\tau) $$

<figure><img src="image_3.png" /><figcaption>Fig. 2. Global response characteristics on the parameter plane $\mu-\Omega$ when $\beta_{cr}=200.0$.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>Fig. 3. Global response characteristics on parameter plane $\mu-\Omega$ when $\beta_{cr}=2.0$.</figcaption></figure>

Defining the state vector in the inertial coordinate system

$$ \mathbf{x}(\tau) = \left\{ y_r \quad z_r \quad y_s \quad z_s \quad y'_r \quad z'_r \quad y'_s \quad z'_s \right\}^T $$

Eq. (2) can be written into the first-order differential equations in the state space

$$ \dot{\mathbf{x}} = A\mathbf{x} + g(\mathbf{x}, \tau) = G(\mathbf{x}, \tau) \quad (8) $$

where $A$ and $g(\mathbf{x}, \tau)$ are the coefficient matrix and the nonlinear vector (see [21]).

After linearizing Eq. (8) about the synchronous full annular rub solution, which is denoted by

$$ \mathbf{x}^0 = \left\{ B_r \cos \varphi_\alpha \quad B_r \sin \varphi_\alpha \quad B_s \cos \varphi_\beta \quad B_s \sin \varphi_\beta \quad -\Omega B_r \sin \varphi_\alpha \quad \Omega B_r \cos \varphi_\alpha \quad -\Omega B_s \sin \varphi_\beta \quad \Omega B_s \cos \varphi_\beta \right\}^T $$

with $\varphi_\alpha=\Omega\tau+\alpha$, $\varphi_\beta=\Omega\tau+\beta$, a set of time-variant linear differential equations is obtained:

$$ \delta \dot{\mathbf{x}} = [J(\mathbf{x}^0, \tau)] \delta \mathbf{x} $$

where $\delta \mathbf{x} = \mathbf{x} - \mathbf{x}^0$ is the perturbation to the synchronous full annular rub solution. The so-called Jacobian matrix is given as

$$ [J(\mathbf{x}^0, \tau)]_{8 \times 8} = \left. \frac{\partial \mathbf{G}}{\partial \mathbf{x}} \right|_{\mathbf{x} = \mathbf{x}^0} = \mathbf{A} + \left. \frac{\partial \mathbf{g}}{\partial \mathbf{x}} \right|_{\mathbf{x} = \mathbf{x}^0} $$

Using the Floquet theory, the stability of the synchronous full annular rub solution can be determined through the examination of the sign of the real parts of the Floquet exponents, which is governed by the following eigenvalue problem with $[J_n]$ is given in [21]:

$$ [J_n] - \lambda I_{8 \times 8} = 0 $$

The solution is stable if all the real parts of its characteristic roots are less than zero. Otherwise it is unstable. Since the stability boundaries are of interest, bifurcation theory is used as is done in [3]. The saddle node bifurcation, denoted by SN in Fig. 2, is one of the existence boundaries of the synchronous full annular rub. The Hopf bifurcation condition, denoted by HP in Fig. 2, gives the boundary between the synchronous full annular rub and the quasi-periodic partial rub. It is found that there is no synchronous full annular rub but only the partial rub when the rotor is in contact with the stator in case of Fig. 3 where the contact stiffness is relative small.

## 4. The boundary of partial rubs

Quasi-period partial rubbing is a response of the rotor/stator system by which the rotor rubs the stator intermittently. So Eq. (2) with piecewise smooth equations of $\Theta=0$ and 1 should be analyzed together. This is a difficult task to get the solution of the response. However, the existence boundaries of partial rub can be determined based on the study of other responses. For example, it is easy to determine the boundary between no-rub and partial rub. Since the synchronous full annular rub bifurcates into partial rub through Hopf bifurcation, the Hopf bifurcation condition of synchronous full annular rub serves as the boundary between synchronous full annular rub and partial rub. The partial rubs cease to exist after dry friction backward whirl is triggered by the imbalance. As known from [13], the rotor in resonance at the natural frequency of the nonlinear normal modes of the coupled rotor/stator system is the reason for onset of dry friction backward whirl. From this understanding, a method proposed in this paper to determine the approximate onset boundary.
## 4.1. Natural frequencies of nonlinear normal modes

To derive the formula for calculation of the natural frequencies of the nonlinear normal mode, we set right-hand side of Eq. (2) to zero and omit the dry friction effect to get the so-called unforced coupled nonlinear rotor/stator system:

$$ 
\begin{aligned}
& \dot{r}_r'' + 2\zeta_r \dot{r}_r' + \dot{r}_r + \beta_{cr}(1+j\mu) \left( 1 - \frac{\Delta}{|\dot{r}_r - \dot{r}_s|} \right) (\dot{r}_r - \dot{r}_s) = 0 \\
& M_{sr} \dot{r}_s'' + 2\zeta_s \sqrt{M_{sr} \beta_{sr} \dot{r}_s'} + \beta_{sr} \dot{r}_s - \beta_{cr}(1+j\mu) \left( 1 - \frac{\Delta}{|\dot{r}_r - \dot{r}_s|} \right) (\dot{r}_r - \dot{r}_s) = 0
\end{aligned}
\quad (9) $$

Supposed the solution has the following form:

$$ \dot{r}_r = H_r e^{j\omega_w \tau}, \quad \dot{r}_s = H_s e^{j(\omega_w \tau + \varphi)} $$

By substituting above solution into Eq. (9), it yields

$$ 
\begin{aligned}
& -H_r \omega_w^2 + 2\zeta_r H_r \omega_w j + H_r + \beta_{cr}(1+j\mu) \left( 1 - \frac{\Delta}{|H_r - H_s|} \right) (H_r - H_s e^{j\varphi}) = 0 \\
& -M_{sr} H_s \omega_w^2 e^{j\varphi} + 2\zeta_s \sqrt{M_{sr} \beta_{sr} H_s \omega_w} e^{j\varphi} j + \beta_{sr} H_s e^{j\varphi} - \beta_{cr}(1+j\mu) \left( 1 - \frac{\Delta}{|H_r - H_s|} \right) (H_r - H_s e^{j\varphi_s}) = 0
\end{aligned}
$$

Adding the two equations to solve $H_r$ as the function of $H_s$ and $\varphi$, we get

$$ H_r = \frac{\beta_{sr} - M_{sr} \omega_w^2 + 2j \zeta_s \sqrt{M_{sr} \beta_{sr} \omega_w}}{1 - \omega_w^2 + 2j \zeta_r \omega_w} H_s e^{j\varphi} \quad (10) $$

Through some manipulation to delete $H_r$, $H_s$ and $\varphi$, one obtains a real polynomial equations about the natural frequency $\omega_w$,

$$ A_6 \omega_w^6 + A_5 \omega_w^5 + A_4 \omega_w^4 + A_3 \omega_w^3 + A_2 \omega_w^2 + A_1 \omega_w + A_0 = 0 \quad (11) $$

where

$$ 
\begin{aligned}
A_0 &= -\mu\beta_{sr}(1 + \beta_{sr}), \quad A_1 = 2\zeta_r\beta_{sr}^2 + 2\zeta_s\sqrt{M_{sr}\beta_{sr}} \\
A_2 &= 2M_{sr}\mu\beta_{sr} + M_{sr}\mu + 2\mu\beta_{sr} + \mu\beta_{sr}^2 - 4\mu\zeta_r^2\beta_{sr} - 4\mu\zeta_s^2 M_{sr} \\
A_3 &= 8\zeta_r^2\zeta_s\sqrt{M_{sr}\beta_{sr}} - 4\zeta_r\beta_{sr}M_{sr} + 8\zeta_r\zeta_s^2 M_{sr}\beta_{sr} - 4\zeta_s\sqrt{M_{sr}\beta_{sr}} \\
A_4 &= -2\mu\beta_{sr}M_{sr} - \mu\beta_{sr} + 4\mu\zeta_r^2 M_{sr} + 4\mu\zeta_s^2 M_{sr}\beta_{sr} - M_{sr}^2\mu - 2M_{sr}\mu \\
A_5 &= 2\zeta_r M_{sr}^2 + 2\zeta_s\sqrt{M_{sr}\beta_{sr}}, \quad A_6 = M_{sr}^2\mu + M_{sr}\mu
\end{aligned}
$$

Below the mode shapes are explored, which can be written as

$$ \frac{H_s e^{j\varphi}}{H_r} = Y_s + Z_s j $$

From Eq. (10) we know

$$ Y_s = \frac{(1 - \omega_n^2)(\beta_{sr} - \omega_n^2 M_{sr}) + 4\zeta_r \zeta_s \omega_n^2}{(\beta_{sr} - \omega_n^2 M_{sr})^2 + 4\zeta_s^2 \omega_n^2}, \quad Z_s = \frac{2\zeta_r \omega_n (\beta_{sr} - \omega_n^2 M_{sr}) - 2\zeta_s \omega_n (1 - \omega_n^2)}{(\beta_{sr} - \omega_n^2 M_{sr})^2 + 4\zeta_s^2 \omega_n^2} $$

The modulus of $Y_s + jZ_s$ is the amplitude ratio of the stator to the rotor during backward whirl motion, and its argument gives the phase difference between the rotor and the stator.

It should be noted that not all roots of Eq. (11) are the natural frequencies of the nonlinear normal modes of the coupled nonlinear rotor/stator system except for those which are real and make $|H_{rs}| = |H_r - H_s| > \Delta$ so as to meet the condition of coupling between the rotor and the stator. There are generally three negative roots meeting the condition, that is, the natural frequencies of the nonlinear normal modes are negative, implying that the nonlinear normal modes correspond to the motions of backward whirl.

## 4.2. Critical rotating speed for the onset of dry friction backward whirl

To use the finding that the loss of stability of partial rub is due to the resonance of system at the nonlinear natural frequencies, Eq. (2) is transformed into a form that describes the dynamics of the rotor/stator system near the nonlinear
nature frequencies:

$$
\begin{equation}
\begin{aligned}
& \hat{r}'' + 2\zeta_r \hat{r}' + \hat{r} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) (\hat{r}_r - \hat{r}_s) = \Omega^2 e^{j\Omega\tau} + \Delta\beta_{cr}(1+j\mu) \left(\frac{1}{|\hat{r}_r - \hat{r}_s|} - \frac{1}{|H_{rs}|}\right) (\hat{r}_r - \hat{r}_s) \\
& M_{sr}\hat{r}''_s + 2\zeta_s \sqrt{M_{sr}\beta_{sr}\hat{r}'_s} + \beta_{sr}\hat{r}_s - \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) (\hat{r}_r - \hat{r}_s) = -\Delta\beta_{cr}(1+j\mu) \left(\frac{1}{|\hat{r}_r - \hat{r}_s|} - \frac{1}{|H_{rs}|}\right) (\hat{r}_r - \hat{r}_s)
\end{aligned}
\tag{12}
\end{equation}
$$

where the left-hand sides have the natural frequencies equal to $\omega_w$. An approximated linear equation to (12) can be set up by defining

$$
\epsilon = \frac{1}{|\hat{r}_r - \hat{r}_s|} - \frac{1}{|H_{rs}|} = o(\epsilon) \le 0 \quad (13)
$$

where $\epsilon$ is a small parameter. The reason for doing this is that $|H_{rs}|$ is the amplitude of the rotor/stator system in backward
whirl without imbalance excitation and $|\hat{r}_r - \hat{r}_s|$ is the amplitude of the system in backward whirl with imbalance
excitation. So $|\hat{r}_r - \hat{r}_s|$ should be only a small modification on $|H_{rs}|$ with $\Omega^2$ is not large. Thus the difference between the
corresponding inverse relations should be also small. By denoting $W_r = \hat{r}_r$, $W_s = \hat{r}_s$, Eq. (12) can be rewritten
approximately as

$$
\begin{align}
& W_r'' + 2\zeta_r W_r' + W_r + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|W_r - W_s|}\right) (W_r - W_s) = \Omega^2 e^{j\Omega\tau} + \Delta\beta_{cr}(1+j\mu)\varepsilon(W_r - W_s) \nonumber \\
& M_{sr}W_s'' + 2\zeta_s \sqrt{M_{sr}\beta_{sr}W_s'} + \beta_{sr}W_s - \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|W_r - W_s|}\right) (W_r - W_s) = -\Delta\beta_{cr}(1+j\mu)\varepsilon(W_r - W_s) \tag{14}
\end{align}
$$

The approximated equation can be solved by the multiple scale method [20], at the rotating speed near the magnitude
of the nonlinear natural frequency, $\Omega = |\omega_w| + \varepsilon\sigma$. $\sigma$ is a small detuning. Assume that the problem depends on many time
scales:

$$
T_0 = \tau, \quad T_1 = \varepsilon\tau, \quad T_2 = \varepsilon^2\tau, \dots
$$

Suppose

$$
W_r(\tau, \epsilon) = W_{r0} + \epsilon W_{r1} + \dots
$$

$$
W_s(\tau, \epsilon) = W_{s0} + \epsilon W_{s1} + \dots
$$

Substitute them into (14) to get the equations of order $\epsilon^0$ and $\epsilon^1$ are as follows:

Order $\epsilon^0$:

$$
\begin{align}
& D_0^2 W_{r0} + 2\zeta_r D_0 W_{r0} + W_{r0} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{r0} \nonumber \\
& \qquad = \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{s0} + \Omega^2 e^{j\Omega\tau} \tag{15}
\end{align}
$$

$$
\begin{align}
& M_{sr} D_0^2 W_{s0} + 2 \zeta_s D_0 W_{s0} + \beta_{sr} W_{s0} + \beta_{cr} (1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{s0} \nonumber \\
& \qquad = \beta_{cr} (1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{r0} \tag{16}
\end{align}
$$

Order $\epsilon^1$:

$$
\begin{align}
& D_0^2 W_{r1} + 2\zeta_r D_0 W_{r1} + W_{r1} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{r1} \nonumber \\
&= -2D_0 D_1 W_{r0} - 2\zeta_r D_1 W_{r0} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{s1} \nonumber \\
&\quad + \beta_{cr}\Delta(1+j\mu)(W_{r0} - W_{s0}) \tag{17}
\end{align}
$$

$$
\begin{align}
& M_{sr} D_0^2 W_{s1} + 2\zeta_s D_0 W_{s1} + \beta_{sr} W_{s1} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{s1} \nonumber \\
&= -2M_{sr} D_0 D_1 W_{s0} - 2\zeta_s D_1 W_{s0} + \beta_{cr}(1+j\mu) \left(1 - \frac{\Delta}{|H_{rs}|}\right) W_{r1} \nonumber \\
&\quad - \beta_{cr}\Delta(1+j\mu)(W_{r0} - W_{s0}) \tag{18}
\end{align}
$$
Assume $W_{r0}^* = A_r e^{j\Omega\tau}$ and $W_{s0}^* = A_s e^{j\Omega\tau}$, then substitute them into Eqs. (15) and (16) to get

$$(1+k-\Omega^2+2\zeta_r\Omega j)A_r = kA_s+\Omega^2$$

$$(\beta_{sr}+k-M_{sr}\Omega^2+2\zeta_s\Omega j)A_s = kA_r$$

By solving above equations, we obtain

$$A_s = \frac{k\Omega^2}{(1+k-\Omega^2+2\zeta_r\Omega j)(\beta_{sr}+k-M_{sr}\Omega^2+2\zeta_s\sqrt{M_{sr}\beta_{sr}\Omega j})-k^2}$$

$$A_r = \frac{\Omega^2(\beta_{sr}+k-M_{sr}\Omega^2+2\zeta_s\sqrt{M_{sr}\beta_{sr}\Omega j})}{(1+k-\Omega^2+2\zeta_r\Omega j)(\beta_{sr}+k-M_{sr}\Omega^2+2\zeta_s\sqrt{M_{sr}\beta_{sr}\Omega j})-k^2}$$

$$k = \beta_{cr}(1+j\mu) \left( 1 - \frac{\Delta}{|H_{rs}|} \right)$$

The solution of Eqs. (15) and (16) are in the form

$$\begin{aligned} W_{r0} &= a_r(T_1)e^{j(\omega_w T_0 + \theta_r(T_1))} + A_r e^{j\Omega T_0} \\ W_{s0} &= a_s(T_1)e^{j(\omega_w T_0 + \theta_s(T_1))} + A_s e^{j\Omega T_0} \end{aligned} \quad (19)$$

Substituting (19) into (17) and (18) leads to

$$\begin{aligned} D_0^2 W_{r1} + 2\zeta_r D_0 W_{r1} + (1+k)W_{r1} - kW_{s1} &= -2\omega_w [a'_r(T_1)j - a_r(T_1)\theta'_r(T_1)]e^{j(\omega_w T_0 + \theta_r(T_1))} \\ &\quad -2\zeta_r [a'_r(T_1) + a_r(T_1)\theta'_r(T_1)]e^{j(\omega_w T_0 + \theta_r(T_1))} + \beta_{cr}\Delta(1+j\mu)a_r(T_1)e^{j(\omega_w T_0 + \theta_r(T_1))} \\ &\quad -\beta_{cr}\Delta(1+j\mu)a_s(T_1)e^{j(\omega_w T_0 + \theta_s(T_1))} + \beta_{cr}\Delta(1+j\mu)(A_r-A_s)e^{j\Omega T_0} \end{aligned} \quad (20)$$

$$\begin{aligned} M_{sr} D_0^2 W_{s1} + 2\zeta_s D_0 W_{s1} + (\beta_{sr}+k)W_{s1} - kW_{r1} &= -2M_{sr}\omega_w [a'_s(T_1)j - a_s(T_1)\theta'_s(T_1)]e^{j(\omega_w T_0 + \theta_s(T_1))} \\ &\quad -2\zeta_s \omega_w [a'_s(T_1) + a_s(T_1)\theta'_s(T_1)]e^{j(\omega_w T_0 + \theta_s(T_1))} - \beta_{cr}\Delta(1+j\mu)a_r(T_1)e^{j(\omega_w T_0 + \theta_r(T_1))} \\ &\quad + \beta_{cr}\Delta(1+j\mu)a_s(T_1)e^{j(\omega_w T_0 + \theta_s(T_1))} - \beta_{cr}\Delta(1+j\mu)(A_r-A_s)e^{j\Omega T_0} \end{aligned} \quad (21)$$

Now substituting the above expression of $\Omega$ into (20) and (21) and setting the coefficients of terms with $\exp(j\omega_0 T_0)$ on the right-hand side to zero in order to eliminate the secular terms. By separating the real and the imaginary parts of the obtained equations, we get

$$a'_r = c_1(\zeta_r + \mu\omega_w)a_r + c_1 c_2 a_s \sin(\theta_r - \theta_s - \alpha_1 - \alpha_4) + c_1 c_2 c_3 \sin(\sigma T_1 + \alpha_1 + \alpha_2 - \theta_r)$$

$$a_r \theta'_r = c_1(\mu\zeta_r - \omega_w)a_r + c_1 c_2 a_s \cos(\theta_r - \theta_s - \alpha_1 - \alpha_4) - c_1 c_2 c_3 \cos(\sigma T_1 + \alpha_1 + \alpha_2 - \theta_r)$$

$$a'_s = c_4(\zeta_s + \mu M_{sr} \omega_w)a_s - c_4 c_5 a_r \sin(\theta_r - \theta_s + \alpha_3 + \alpha_4) - c_4 c_5 c_3 \sin(\sigma T_1 + \alpha_2 + \alpha_3 - \theta_s)$$

$$a_s \theta'_s = c_4(\mu\zeta_s - M_{sr}\omega_w)a_s + c_4 c_5 a_r \cos(\theta_r - \theta_s + \alpha_3 + \alpha_4) + c_4 c_5 c_3 \cos(\sigma T_1 + \alpha_2 + \alpha_3 - \theta_s) \quad (22)$$

where

$$A_r = A_{rr} + jA_{ri}, \quad A_s = A_{sr} + jA_{si}, \quad c_1 = \frac{\beta_{cr}\Delta}{2(\zeta_r^2 + \omega_w^2)}, \quad c_2 = \sqrt{\zeta_r^2 + \omega_w^2}, \quad c_3 = \sqrt{d_1^2 + d_2^2}, \quad c_4 = \frac{\beta_{cr}\Delta}{2(\zeta_s^2 + M_{sr}^2 \omega_w^2)}$$

$$c_5 = \sqrt{\zeta_s^2 + M_{sr}^2 \omega_w^2}, \quad \alpha_1 = \arctan \frac{\zeta_r}{\omega_w}, \quad \alpha_3 = \arctan \frac{\zeta_s}{M_{sr}\omega_w}$$

$$d_1 = A_{rr} - A_{sr} - \mu(A_{ri} - A_{si}), \quad d_2 = \mu(A_{rr} - A_{sr}) + A_{ri} - A_{si}$$

$$\alpha_2 = \arctan \frac{d_1}{d_2} = \arctan \frac{A_{rr} - A_{sr} - \mu(A_{ri} - A_{si})}{\mu(A_{rr} - A_{sr}) + A_{ri} - A_{si}}, \quad \alpha_4 = \arctan \mu$$

Next (22) is transformed into an autonomous system with $T_1$ not explicitly appearing. Let

$$\gamma_r = \sigma T_1 + \alpha_1 + \alpha_2 - \theta_r \quad \text{and} \quad \gamma_s = \sigma T_1 + \alpha_2 + \alpha_3 - \theta_s$$

and substitute them into (22) to yield

$$a'_r = c_1(\zeta_r + \mu\omega_w)a_r + c_1 c_2 a_s \sin(\gamma_s - \gamma_r - \alpha_3 - \alpha_4) + c_1 c_2 c_3 \sin(\gamma_r)$$

$$a_r \gamma'_r = a_r \sigma - c_1(\mu\zeta_r - \omega_w)a_r - c_1 c_2 a_s \cos(\gamma_s - \gamma_r - \alpha_3 - \alpha_4) + c_1 c_2 c_3 \cos(\gamma_r)$$

$$a'_s = c_4(\zeta_s \sqrt{M_{sr}\beta_{sr}} + \mu M_{sr}\omega_w)a_s + c_4 c_5 a_r \sin(\gamma_s - \gamma_r + \alpha_1 + \alpha_4) - c_4 c_5 c_3 \sin(\gamma_s)$$

$$a_s \gamma'_s = c_4(\mu\zeta_s - M_{sr}\omega_w)a_s + c_4 c_5 a_r \cos(\gamma_s - \gamma_r + \alpha_1 + \alpha_4) + c_4 c_5 c_3 \cos(\gamma_s) \quad (23)$$

Let right-hand side equal zero, $a_r$ and $a_s$ could be solved from Eq. (23) by using the phase difference between the rotor and the stator from the normal mode analysis in Section 4.1. Then, $|\dot{r}_r - \dot{r}_s| \approx |W_{r0} - W_{s0}|$ can be obtained. The critical rotating speed, $\Omega = |\omega_w| + \varepsilon\sigma$, for the onset of the dry friction backward whirl or the loss of stability of the partial rub can be
solved after the deviation $\epsilon\sigma$ to $|\omega_w|$ can be determined. To do this, we can first set $\sigma$ to a reasonable constant number, i.e., $\sigma=1$ as used in this paper, then a trial and error process one is carried out in order to find the most proper that make it a largest negative value, i.e., mostly satisfying (13). In this way we determine the boundary by which dry whip is triggered on from partial rub.

In Figs. 2 and 3, curve DW represents the boundary between the partial rub and the dry friction backward whirl. It is found that the critical rotating speeds are near the magnitude of the first nonlinear natural frequency in both cases. The points marked by triangular stand for the boundary points where the dry friction backward whirl is triggered on from partial rub by the imbalance. It shows good agreement between the predicted results and the simulated ones.

## 5. The existence boundaries of dry friction backward whirl

As observed in tests [4] and in simulations [12], the response of the dry friction backward whirl of the rotor/stator system is composed of a backward whirl motion with a negative super-synchronous whirl frequency and a forced forward whirl motion with a synchronous whirl frequency whose amplitude is usually negligible. The solutions of dry friction backward whirl can, thus, be assumed in the form

$$ \hat{r}_i = W_{bi} + W_{fi} \quad \text{with} \quad W_{bi} = H_i e^{(\alpha + j\omega_b)\tau + j\phi_i} \quad \text{and} \quad W_{fi} = A_i e^{i\Omega\tau}, \quad i=r,s \quad (24) $$

where $W_{bi}$ represents the solution part of the backward whirl motion with $H_i e^{\alpha\tau}$ the time-varying amplitude and $\omega_b$ the negative whirl frequency to be determined, and $W_{fi}$ the solution part of the forced synchronous forward whirl motion with $A_i$ to be determined. Without loss of generality, we take $\phi_r=0$, $\phi_s=\phi$:

$$ |W_r - W_s| = |(W_{br} + W_{fr}) - (W_{bs} + W_{fs})| = -\frac{R_{disk}\Omega}{\omega_b} \quad (25) $$

Substituting (24) and (25) into the third equation of Eq. (2), and forced forward whirl motion with a synchronous whirl frequency whose amplitude is usually negligible, by setting $\Theta=1$ and the right-hand term to zero:

$$ \begin{aligned} & [(\alpha^2 - \omega_b^2 + 2j\alpha\omega_b) + 2\bar{\zeta}_r(\alpha + j\omega_b) + 1]H_r + (1 + j\mu)\beta_{cr}\left(1 + \frac{\omega_b A}{R_{disk}\Omega}\right)(H_r - H_s e^{j\phi}) = 0 \\ & \left[M_{sr}(\alpha^2 - \omega_b^2 + 2j\alpha\omega_b) + 2\bar{\zeta}_s(\alpha + j\omega_b) + \beta_{sr}\right]H_s e^{j\phi} - (1 + j\mu)\beta_{cr}\left(1 + \frac{\omega_b A}{R_{disk}\Omega}\right)(H_r - H_s e^{j\phi}) = 0 \end{aligned} \quad (26) $$

With some manipulation as done in [22], a complex equation of $\alpha$ and $\omega_b$ as function of system parameters is obtained, which yields two real polynomial equations about the two unknown variables $\alpha$ and $\omega_b$ after separating the real and the imaginary parts:

$$ \begin{aligned} & [M_{sr}(\alpha^2 - \omega_b^2) + 2\bar{\zeta}_s\alpha + \beta_{sr}][(\alpha^2 - \omega_b^2) + 2\bar{\zeta}_r\alpha + 1] - 4\omega_b^2(\alpha M_{sr} + \bar{\zeta}_s)(\alpha + \bar{\zeta}_r) \\ & + \beta_{cr}(1 + C\omega_b)[((1 + M_{sr})(\alpha^2 - \omega_b^2) + 2\alpha(\bar{\zeta}_r + \bar{\zeta}_s) + (1 + \beta_{sr})) - 2\omega_b\mu[\alpha(1 + M_{sr}) + (\bar{\zeta}_r + \bar{\zeta}_s)]] = 0 \end{aligned} \quad (27) $$

$$ \begin{aligned} & 2\omega_b[(\alpha + \bar{\zeta}_r)[M_{sr}(\alpha^2 - \omega_b^2) + 2\bar{\zeta}_s\alpha + \beta_{sr}] + (\alpha M_{sr} + \bar{\zeta}_s)[(\alpha^2 - \omega_b^2) + 2\bar{\zeta}_r\alpha + 1]] \\ & + \beta_{cr}(1 + C\omega_b)[2\omega_b[\alpha(1 + M_{sr}) + (\bar{\zeta}_r + \bar{\zeta}_s)] + \mu[(1 + M_{sr})(\alpha^2 - \omega_b^2) + 2\alpha(\bar{\zeta}_r + \bar{\zeta}_s) + (1 + \beta_{sr})]] = 0 \end{aligned} \quad (28) $$

where $C=A/R_{disk}\Omega$ is the sole parameter containing the rotating speed.

As discussed above, the solution of dry friction backward whirl should have an exponent with a positive real part, namely $\alpha > 0$. The existence condition of dry friction backward whirl is the critical case that dry friction backward whirl starts to exist, that is, the real part of the exponent in the solution should be zero, namely, $\alpha=0$. Applying this condition to Eqs. (27) and (28) yields two equations of $\omega_b$

$$ \begin{aligned} & \beta_{cr}\omega_b[-(1 + M_{sr})\omega_b^2 - 2\omega_b\mu(\bar{\zeta}_r + \bar{\zeta}_s) + (\beta_{sr} + 1)]C \\ & = 4\omega_b^2\bar{\zeta}_r\bar{\zeta}_s - \beta_{cr}[-(1 + M_{sr})\omega_b^2 - 2\omega_b\mu(\bar{\zeta}_r + \bar{\zeta}_s) + (\beta_{sr} + 1)] - (M_{sr}\omega_b^2 - \beta_{sr})(\omega_b^2 - 1) \end{aligned} \quad (29) $$

$$ \begin{aligned} & \{\mu\beta_{cr}\omega_b[-(1 + M_{sr})\omega_b^2 + (\beta_{sr} + 1)] + 2\omega_b^2(\bar{\zeta}_r + \bar{\zeta}_s)\}C \\ & = 2\omega_b[\bar{\zeta}_r(M_{sr}\omega_b^2 - \beta_{sr}) + \bar{\zeta}_s(\omega_b^2 - 1)] - \mu\beta_{cr}[-(1 + M_{sr})\omega_b^2 + (\beta_{sr} + 1)] - 2\omega_b\beta_{cr}(\bar{\zeta}_r + \bar{\zeta}_s) \end{aligned} \quad (30) $$

To solve the critical rotating speed for the start of the existence of dry friction backward whirl, we first delete $C$ from Eqs. (29) and (30) to get a polynomial equation containing only $\omega_b$ as the unknown, then solve the equation for $\omega_b$, after that substitute them into either of Eq. (29) or (30) to yield a equation containing $C$ as the unknown, finally solve the equation for $C$ and thus the critical rotating speeds for the existence of dry friction backward whirl. By changing the parameters, we can determine the boundary for dry friction backward whirl in the parameter space.

As pointed out in [22], there are two existence regions of the dry friction backward whirl with different backward whirl frequencies. Fig. 2 shows that the two existence regions are isolated with their boundaries marked by DF1 and DF2 for the case of relatively rigid contact surface. For the case with soft contact surface as shown in Fig. 3, the second existence region
is completely contained in the first one. So there is coexistence of two types of dry friction backward whirls in the parameter region of the second existence region.

## 6. Global response characteristics

From the above analysis the boundaries of the synchronous no-rub motion, the synchronous full annular rub, the partial rub and the dry friction backward whirl are determined. With these results the characteristics of the global dynamics of the rotor/stator system can be depicted. In Figs. 2 and 3, the characteristics of the global dynamics of the rotor/stator system on the parameter plane $\mu-\Omega$ with the system parameters $M_{sr}=0.2$, $\zeta_r=0.05$, $\zeta_s=0.05$, $\beta_{sr}=2.0$, $\Delta=2.0$ and $R_{disk}=20$ are demonstrated as examples.

For assistance of the reading, boundaries of different responses in Figs. 2 and 3 are summarized here: curves $\Omega_L$ and $\Omega_U$ are boundaries for the synchronous no-rub motion. Curve SN and HP represent, respectively, the boundaries of saddle-node bifurcation and Hopf bifurcation of the synchronous full annular rub solutions. Curve DW is the boundary by which dry friction backward whirl is triggered by imbalance and the partial rub ceases to exist. Curve DFi stands for the boundary of the existence region of the ith dry friction backward whirl.

The parameter regions in $\mu-\Omega$ plane of Figs. 2 and 3 where only one kind of response exists are denoted with the corresponding names. Especially, the synchronous full annular rub is abbreviated as 'sfar' and the dry friction backward whirl is denoted shortly as dry whip. The meanings of abbreviations in Figs. 2 and 3 are summarized in Table 1. For the parameter regions in $\mu-\Omega$ plane where coexistence of different responses occurs are marked with numbers and the types of coexistence are listed in Table 2.

As can be seen from Table 2, there are at least ten types of coexistence of different responses. The experimental observed coexistence phenomena in rotor/stator rubbing tests reported in the literature are all included in the results. For instance, the coexistence types of 1-3 were reported in the experimental works [4,16,17]. The coexistence type 10 should correspond to what was observed experimentally in [18]. Coexistence type 7 should be the one shown in the tests of Ref. [19].

There are significant differences in the global response characteristics of the rotor/stator rubbing system with rigid or soft contact surfaces. In the case of rigid contact surfaces as shown in Fig. 2, the two existence regions of dry friction backward whirl are isolated. The first existence region bounded by DF1 is very small and generally hard to be detected in the real operation condition. In contrast, the second one with higher backward whirl frequency and larger response amplitude occupies very large region and coexists with almost all other types of responses. In this case, the synchronous full annular rub as mild rubbing is observable. In the case of soft contact surfaces as shown in Fig. 3, the existence region of first dry friction backward whirl bounded by DF1 deforms from the U-shape and overlaps with the second existence region edged by DF2. In this case, the partial rub rather than the synchronous full annular rub appears when the rotor initiates contact with the stator. The synchronous full annular rub occurs in the right-hand side of Hopf bifurcation boundary.

<table><caption>**Table 1**
Meaning of some abbreviations in Figs. 2 and 3.</caption><thead><tr><td>Abbreviation</td><td>Meaning</td></tr></thead><tbody><tr><td>sfar</td><td>Synchronous full annular rub</td></tr><tr><td>SN</td><td>Saddle-node bifurcation boundary of synchronous full annular rub solution</td></tr><tr><td>HP</td><td>Hopf bifurcation boundary of synchronous full annular rub solution</td></tr><tr><td>DF</td><td>Existence boundary of dry friction backward whirls</td></tr><tr><td>DW</td><td>Onset boundary of dry friction backward whirl</td></tr></tbody></table>

<table><caption>**Table 2**
List of different types of coexistences.</caption><thead><tr><td>Mark</td><td>Different responses in coexistence</td></tr></thead><tbody><tr><td>"1"</td><td>No-rub motion+1st dry friction backward whirl</td></tr><tr><td>"2"</td><td>No-rub motion+2nd dry friction backward whirl</td></tr><tr><td>"3"</td><td>No-rub motion+synchronous full annular rub</td></tr><tr><td>"4"</td><td>Synchronous full annular rub+2nd dry friction backward whirl</td></tr><tr><td>"5"</td><td>Partial rub+1st dry friction backward whirl</td></tr><tr><td>"6"</td><td>Partial rub+2nd dry friction backward whirl</td></tr><tr><td>"7"</td><td>No-rub motion+synchronous full annular rub+2nd dry friction backward whirl</td></tr><tr><td>"8"</td><td>No-rub motion+1st dry friction backward whirl+2nd dry friction backward whirl</td></tr><tr><td>"9"</td><td>Partial rub+1st dry friction backward whirl+2nd dry friction backward whirl</td></tr><tr><td>"10"</td><td>1st dry friction backward whirl+2nd dry friction backward whirl</td></tr></tbody></table>
## 7. Conclusions

In this paper, the global response characteristics of a general rotor/stator system with dry friction effects are investigated. The rotor/stator rubbing model takes into account both the dynamics of the rotor and the stator as well as the friction and the deformation on the contact surfaces. Most importantly, the dry friction effects, which are mostly neglected in the published works on the rotor/stator rubbing, are also considered in the present model in order to study the self-excited dry friction backward whirl.

In order to get a deep insight into the global response characteristics of this model has demonstrated of this system, the exact and approximated analytical existence boundaries of the typical responses, such as the synchronous no-rub motion, the synchronous full annular rub, the partial rub and the dry friction backward whirl, are derived by analyzing different pieces of the piecewise smooth nonlinear governing equations one by one.

Through our analysis, two regions of synchronous no-rub response with two different ranges of whirl frequencies, one region of the synchronous full annular rub and one region of the quasi-periodic partial rub as well as two regions of dry friction backward whirl with two different ranges of backward whirl frequencies are determined. It is found that the contact stiffness plays a dominant role in deciding the characteristics of the dry friction backward whirl. In the cases of a relatively large, there are two separate existence regions with the first region being U-shaped. In the case of relatively small contact stiffness, the first existence region may deform away from U-shape and overlap with the second existence region so that the coexistence of two dry friction backward whirl responses with two different backward whirl frequencies occurs. For some friction coefficients the jump phenomenon in the whirl frequencies during the run-up or run-down of the rotor may be detected as observed in [18]. There are totally ten types of coexistence of different responses being found in this paper that cover all the types of coexistence observed in the tests of rotor/stator rubbing until now. Moreover, the present results provide an overall picture on the global response characteristics of the rotor/stator rubbing system that shows how the coexistences in dependence on the system parameters.

## Acknowledgement

This work was supported by the National Natural Science Foundation of China (NSFC) under the Grant nos. 10872155, 10472086 and 10772140.

## References

[1] A. Muszynska, Rotor-to-stationary element rub-related vibration phenomena in rotating machinery—literature survey, *Sound and Vibration Digest* 21 (1989) 3–11.

[2] H.F. Black, Interaction of a whirling rotor with a vibrating stator across a clearance annulus, *International Journal of Mechanical Engineering Science* 10 (1968) 1–12.

[3] J. Jiang, H. Ulbrich, Stability analysis of sliding whirl in a nonlinear jeffcott rotor with cross-coupling stiffness coefficients, *Nonlinear Dynamics* 24 (2001) 269–283.

[4] J.J. Yu, P. Goldman, D.E. Bently, A. Muszynska, Rotor/seal experimental and analytical study on full annular rub, *Journal of Engineering for Gas Turbines and Power—Transactions of the ASME* 124 (2002) 340–350.

[5] D.W. Childs, Rub-induced parametric excitation in rotors, *Journal of Mechanical Design* 101 (1979) 640–644.

[6] F.F. Ehrich, High order sub-harmonic response of high speed rotors in bearing clearance, *Journal of Vibration Acoustics, Stress and Reliability in Design—Transactions of the ASME* 110 (1988) 9–16.

[7] Y.B. Kim, S.T. Noah, Quasi-periodic response and stability analysis for a nonlinear jeffcott rotor, *Journal of Sound and Vibration* 190 (1996) 239–253.

[8] P. Goldman, A. Muszynska, Chaotic behavior of rotor/stator system with rubs, *Journal of Engineering for Gas Turbines and Power* 116 (1994) 693–701.

[9] E.V. Karpenko, M. Wiercigroch, M.P. Cartmell, Regular and chaotic dynamics of a discontinuously nonlinear rotor system, *Chaos, Solitons and Fractals* 13 (2002) 1231–1242.

[10] W. Zhang, Dynamic instability of multi-degree-of-freedom flexible rotor systems due to full annular rub, *ImechE C252/88*, 1988, pp. 305–308.

[11] S. Crandall, From whirl to whip in rotor dynamics, *Proceedings of the IFToMM Third International Conference on Rotor Dynamics*, Lyon, France, 1990, pp. 19–26.

[12] J. Jiang, H. Ulbrich, The physical reason and the analytical condition for the onset of dry whip in rotor-to-stator contact systems, *Journal of Acoustics and Vibration—Transactions of the ASME* 127 (2005) 594–603.

[13] J. Jiang, The analytical solution and existence condition of dry friction backward whirl in rotor-to-stator contact systems, *Journal of Acoustics and Vibrations—Transactions of the ASME* 129 (2007) 260–264.

[14] S. Popprath, H. Ecker, Nonlinear dynamics of a rotor contacting an elastically suspended stator, *Journal of Sound and Vibration* 308 (2007) 767–784.

[15] T.T. Yamamoto, On critical speeds of a shaft, *Memoris of Faculty of Engineering, Nagoya University* 6 (2) (1954) 1–20.

[16] A. Muszynska, Synchronous self-excited rotor vibration caused by a full annular rub, *Machinery Dynamics Eighth Seminar*, Halifax, Nova Scotia, Canada, 1984.

[17] A. Lingener, Experimental investigation of reverse whirl of a flexible rotor, *Transactions of IFToMM Third International Conference on Rotordynamics*, Lyon, France, 1990, pp. 13–18.

[18] Y.S. Choi, Investigation on the whirl motion of full annular rotor rub, *Journal of Sound and Vibration* 258 (1) (2002) 191–198.

[19] U. Ehehalt, E. Hahn, R. Markert, Experimental validation of various motion patterns at rotor stator contact, *Proceedings of the ISROMAC-11*, Hawaii, USA, 2006.

[20] A.H. Nayfeh, D.T. Mook, in: *Nonlinear Oscillations*, John Wiley & Sons, 1979.

[21] Z.Y. Shang, J. Jiang, L. Hong, The influence of the cross-coupling effects on the dynamics of rotor/stator rubbing, *Proceedings of the Second Conference on Dynamics, Vibration and Control*, Chendu, China, 2009.

[22] J. Jiang, Z. Shang, L. Hong, Characteristics of dry friction backward whirl—a self-excited oscillation in rotor-to-stator contact system, *Science in China Series E: Technological Sciences* 53 (3) (2010) 674–683.