The nonlinear vibrations of a vertical hard gyroscopic rotor with nonlinear characteristics

The paper considers an impact of viscous linear and cubic nonlinear damping of the elastic support on nonlinear vibrations of a vertical hard gyroscopic unbalanced rotor, taking into account nonlinear stiffness of the support material. Analyzing the research results shows that linear and cubic nonlinear damping can significantly suppress the resonance peak of the fundamental harmonic, eliminate the jumping phenomena of the nonlinear system. In non-resonance areas where the velocity is higher than the critical one, cubic nonlinear damping, unlike linear one, can slightly suppress amplitude of the rotor vibration. Therefore, in the high-velocity area, only nonlinear damping can maintain performance of a vibration isolator. In resonance area, an increase in linear or cubic nonlinear damping significantly suppresses the ability to absolute displacement. In non-resonance area, where the rotational velocity is lower than the critical one, they have almost no impact on ability to absolute displacement. In high velocity area, an increase in nonlinear damping may slightly increase the moment of force transmissibility, but linear damping has almost no impact on it. The obtained results can be successfully used to produce passive vibration isolators used for damping the vibrations of rotary machines, including gyroscopic ones.

1 Introduction

As you know, the rotary machines are widely used in many industries and have been studied for a long time. Despite this, there are many problems unsolved, in particular those associated with an action of mass imbalance on vibrations and stability, and subsequently with the stabilization of resonance vibrations of the rotary machines.

A simplified model with concentrated rotor system parameters is typically used to study the dynamics of a single rotor shaft on the bearing supports. It is very important to use the properties and characteristics of the support material for vibration attenuation and damping to stabilize a movement of the unbalanced rotor and the vibration systems. Supports are the mean for connecting devices between the rotor and the supporting structure, which have various shapes and designs depending on the specific assumptions. A convenient way is to introduce attenuation for supporting the bearings in the rotor system on viscoelastic flexible rubber supports (Zakaria et al., 2015). In parallel with the development of viscoelastic material simulation (Gil-Negrete et al., 2009; Richards and Singh, 1999), which helps to describe a complexity of material properties, the use of viscoelastic components in rotor dynamics and vibration systems (Ravindra and Mallik, 1994; Peng et al., 2012; Ho et al., 2012) as a whole also increased, including with nonlinear elastic characteristics and damping. For example, Ravindra and Mallick (1994) researched the parametric impact of different types of attenuation on the performance of nonlinear vibration isolators in case of harmonic excitation. The paper (Peng et al., 2012) considers an efficiency of passive vibration isolators with linear damping and cubic nonlinear damping in resonance and non-resonance vibration areas of the system with linear stiffness. It also provides an excellent overview of researches of the linear and nonlinear vibration isolation systems. In the paper (Ho et al., 2012), an impact of cubic nonlinear rigidity of the material on the performance of the isolator is additionally taken into account in the studies. The papers (Iskakov, 2015, 2017a, b) researched impacts of quadratic nonlinear damping on the resonance vibrations and stability of the vertical gyroscopic unbalanced rotor with quadratic and cubic nonlinear stiffness of an elastic support. Magnetorheological materials can be used as damping elements. The model is based on presentation of magnetorheological oil by bilinear material (Zapomel and Ferfecki, 2015). In the paper (Fujiwara et al., 2015), an experimental single-disc rotor system is prepared that is supported by ball bearings at both ends, and vibration is compared with a flexible support containing springs or rubber blankets and with a rigid support base, through simulation and experiment.

This paper is aimed at researching an impact of the material of the elastic bearing support with linear and cubic nonlinear damping on the main resonance curve and the moment of force transmissibility of the gyroscopic vertical hard unbalanced rotor, taking into account the nonlinear stiffness of the support.

2 Motion equations

The rotor flowchart shown in Fig. 1 is considered. Shaft with the length L is installed vertically, using a lower hinged and an upper elastic support at a distance l₀ from it. A disc having a mass m, a polar moment of inertia I_P and a transverse moment of inertia I_T, being the same for any direction, is fixed at the free end of the shaft. The shaft speed ω is so high that the rotor can be considered as a gyroscope, which fixed point is the lower shaft support. The geometric center of the rotor S is determined by the Cartesian coordinates x and y. At this time, the shaft deflections are determined by the angular coordinates θx and θy. The Cartesian and angular coordinates of the geometric center S are related by the relations x=Lθx and y=Lθy. The rotor movement towards the coordinate axis z is neglected. Next, we denote the coordinates of the center of masses m of the disc through x_m and y_m.

Figure 1 Rotor geometry.

Also assume that linear eccentricity e lies on the X-axis of the SXYZ coordinate system which rotates with the rotor. Limit its elves to small deviations of the rotor axis, so will take into account in the calculations only terms that are linear with respect to small quantities e, θ_x, θ_y, x_m, y_m.

Given the above, the projections of the angular velocity on the coordinate axes will be written as

$$\begin{array}{}\text{(1)}& {\mathit{\omega }}_x\approx -{\dot{\mathit{\theta }}}_y+\mathit{\omega }{\mathit{\theta }}_x,{\mathit{\omega }}_y\approx {\dot{\mathit{\theta }}}_x-\mathit{\omega }{\mathit{\theta }}_y,{\mathit{\omega }}_z\approx \mathit{\omega }+{\dot{\mathit{\theta }}}_x\cdot {\mathit{\theta }}_y,\end{array}$$

and kinetic energy of the system as

$$T={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}m\left({\dot{x}}_m^{\mathrm{2}}+{\dot{y}}_m^{\mathrm{2}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\left(I_x{\mathit{\omega }}_x^{\mathrm{2}}+I_y{\mathit{\omega }}_y^{\mathrm{2}}+I_z{\mathit{\omega }}_z^{\mathrm{2}}\right),$$

given that $I_x=I_y=I_{\mathrm{T}}$, I_z=I_P and Eq. (1), obtain in a form of

$$\begin{array}{}\text{(2)}& \begin{array}{rl}T& ={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}m\left({\dot{x}}_m^{\mathrm{2}}+{\dot{y}}_m^{\mathrm{2}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\left\{I_{\mathrm{T}}[\left({\dot{\mathit{\theta }}}_x^{\mathrm{2}}+{\dot{\mathit{\theta }}}_y^{\mathrm{2}}\right)\right.\\ & \left.-\mathrm{2}\mathit{\omega }\left({\mathit{\theta }}_x{\dot{\mathit{\theta }}}_y+{\mathit{\theta }}_y{\dot{\mathit{\theta }}}_x\right)]+I_{\mathrm{P}}\left({\mathit{\omega }}^{\mathrm{2}}+\mathrm{2}\mathit{\omega }{\dot{\mathit{\theta }}}_x{\mathit{\theta }}_y\right)\right\},\end{array}\end{array}$$

where

$$\begin{array}{}\text{(3)}& \begin{array}{rl}& x_m=x+e\cos \mathit{\phi }=L{\mathit{\theta }}_x+e\cos \mathit{\omega }t,\\ & y_m=y+e\sin \mathit{\phi }=L{\mathit{\theta }}_y+e\sin \mathit{\omega }t.\end{array}\end{array}$$

Moments of external forces are as follows

$$\begin{array}{}\text{(4)}& M_x=\left(L{\mathit{\theta }}_x+e\cos \mathit{\omega }t\right)G,M_y=\left(L{\mathit{\theta }}_y+e\sin \mathit{\omega }t\right)G,\end{array}$$

where L is the shaft length, G=mg is the disc weight.

Mallick et al. (1999) experimentally confirmed that the restoring and damping forces in elastomeric isolators should be described, using a nonlinear model. Richards and Singh (1999) found that rubber dampers have both nonlinear damping and nonlinear stiffness. To achieve higher performance, any nonlinearities in the design should be taken into account. Consequently, the elastic support of the upper bearing of the gyroscopic rotor can be made of non-linear materials, such as rubber, resin and other polymers widely used as the damper of the arising vibrations. Given all this, we set the dissipative forces in the elastic support as

$$\begin{array}{}\text{(5)}& \mathit{\varphi }={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\mathit{\mu }}_{d\mathrm{1}}\left({\dot{\mathit{\theta }}}_x^{\mathrm{2}}+{\dot{\mathit{\theta }}}_y^{\mathrm{2}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{\mathit{\mu }}_{d\mathrm{3}}\left({\dot{\mathit{\theta }}}_x^{\mathrm{4}}+{\dot{\mathit{\theta }}}_y^{\mathrm{4}}\right),\end{array}$$

where μ_d1 is a coefficient of linear viscous damping of the support material, μ_d3 is a coefficient of nonlinear cubic viscous damping of the support material. If consider that the rotor shaft is hard and only its upper support has elasticity, the potential energy of the system can be represented as

$$\begin{array}{}\text{(6)}& V={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{k}_{\mathrm{1}}{l}_{\mathrm{0}}^{\mathrm{2}}\left({\mathit{\theta }}_x^{\mathrm{2}}+{\mathit{\theta }}_y^{\mathrm{2}}\right)+{\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{k}_{\mathrm{3}}{l}_{\mathrm{0}}^{\mathrm{4}}\left({\mathit{\theta }}_x^{\mathrm{4}}+{\mathit{\theta }}_y^{\mathrm{4}}\right),\end{array}$$

where k₁ is a coefficient of linear stiffness of the elastic support, k₃ is a coefficient of nonlinear stiffness of the elastic support.

The Lagrange's Equations of Second Kind (Yablonsky, 2007) for the rotor system are written as

$$\begin{array}{}\text{(7)}& {\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}\left({\displaystyle \frac{\partial T}{\partial {\dot{q}}_i}}\right)-{\displaystyle \frac{\partial T}{\partial q_i}}+{\displaystyle \frac{\partial V}{\partial q_i}}=-{\displaystyle \frac{\partial \mathit{\varphi }}{\partial {\dot{q}}_i}}+Q_i.\end{array}$$

Here, q_i:θ_x, θ_y is the joint coordinates; Q_i:M_x, M_y are the generalized forces.

Substituting the expressions (2)–(6) in (7), obtain the equations of motion of the rotor

$$\begin{array}{}\text{(8)}& \begin{array}{rl}& \left(I_{\mathrm{T}}+m{L}^{\mathrm{2}}\right){\ddot{\mathit{\theta }}}_x+I_{\mathrm{P}}\mathit{\omega }{\dot{\mathit{\theta }}}_y+{\mathit{\mu }}_{d\mathrm{1}}{\dot{\mathit{\theta }}}_x+{\mathit{\mu }}_{d\mathrm{3}}{\dot{\mathit{\theta }}}_x^{\mathrm{3}}\\ & +\left(k_{\mathrm{1}}{l}_{\mathrm{0}}^{\mathrm{2}}-GL\right){\mathit{\theta }}_x+k_{\mathrm{3}}{l}_{\mathrm{0}}^{\mathrm{4}}{\mathit{\theta }}_x^{\mathrm{3}}\\ & =(me{\mathit{\omega }}^{\mathrm{2}}L+mge)\cos \mathit{\omega }t,\\ & \left(I_{\mathrm{T}}+m{L}^{\mathrm{2}}\right){\ddot{\mathit{\theta }}}_y-I_{\mathrm{P}}\mathit{\omega }{\dot{\mathit{\theta }}}_x+{\mathit{\mu }}_{d\mathrm{1}}{\dot{\mathit{\theta }}}_y+{\mathit{\mu }}_{d\mathrm{3}}{\dot{\mathit{\theta }}}_y^{\mathrm{3}}\\ & +\left(k_{\mathrm{1}}{l}_{\mathrm{0}}^{\mathrm{2}}-GL\right){\mathit{\theta }}_y+k_{\mathrm{3}}{l}_{\mathrm{0}}^{\mathrm{4}}{\mathit{\theta }}_y^{\mathrm{3}}\\ & =(me{\mathit{\omega }}^{\mathrm{2}}L+mge)\sin \mathit{\omega }t,\end{array}\end{array}$$

where

$$\begin{array}{}\text{(9)}& me{\mathit{\omega }}^{\mathrm{2}}L+Ge=M_f\end{array}$$

is an amplitude of the resultant moment of the external forces.

By entering the following dimensionless parameters

$$\begin{array}{}\text{(10)}& \begin{array}{rl}\mathit{\epsilon }& =e/L;l=l_{\mathrm{0}}/L;\stackrel{\mathrm{‾}}{t}=t{\mathit{\omega }}_{\mathrm{0}};\mathrm{\Omega }=\mathit{\omega }/{\mathit{\omega }}_{\mathrm{0}};\\ {\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}& =I_{\mathrm{P}}/\left(m{L}^{\mathrm{2}}\right);{\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}=I_{\mathrm{T}}/\left(m{L}^{\mathrm{2}}\right);{\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}=k_{\mathrm{1}}/\left(m{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\right);\\ {\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}& =k_{\mathrm{3}}{L}^{\mathrm{2}}/\left(m{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\right);P=G/\left(mL{\mathit{\omega }}_{\mathrm{0}}^{\mathrm{2}}\right);\\ {\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}& ={\mathit{\mu }}_{d\mathrm{1}}/\left(m{L}^{\mathrm{2}}{\mathit{\omega }}_{\mathrm{0}}\right);{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}={\mathit{\mu }}_{d\mathrm{3}}{\mathit{\omega }}_{\mathrm{0}}/\left(m{L}^{\mathrm{2}}\right),\end{array}\end{array}$$

where

$$\begin{array}{}\text{(11)}& {\mathit{\omega }}_{\mathrm{0}}=\sqrt{{\displaystyle \frac{k_{\mathrm{1}}{l}_{\mathrm{0}}^{\mathrm{2}}-GL}{m{L}^{\mathrm{2}}-\left(I_{\mathrm{P}}-I_{\mathrm{T}}\right)}}}\end{array}$$

is a critical velocity of a damping-free linear system, using a designation of the expression of the dimensionless amplitude of the resultant moment of external forces

$$\begin{array}{}\text{(12)}& \stackrel{\mathrm{‾}}{M}=\mathit{\epsilon }{\mathrm{\Omega }}^{\mathrm{2}}+\mathit{\epsilon }P\end{array}$$

can provide the equations of motion (8) with a compact dimensionless view

$$\begin{array}{}\text{(13)}& {\displaystyle }\begin{array}{rl}\left(\mathrm{1}+{\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}\right){\mathit{\theta }}_x^{\prime \prime }& +{\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}\mathrm{\Omega }{\mathit{\theta }}_y^{\prime }+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}{\mathit{\theta }}_x^{\prime }+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{{\mathit{\theta }}^{\prime }}_x^{\mathrm{3}}\\ & +\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right){\mathit{\theta }}_x+{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{\mathit{\theta }}_x^{\mathrm{3}}=\stackrel{\mathrm{‾}}{M}\cos \mathrm{\Omega }\stackrel{\mathrm{‾}}{t},\end{array}\text{(14)}& {\displaystyle }\begin{array}{rl}\left(\mathrm{1}+{\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}\right){\mathit{\theta }}_y^{\prime \prime }& -{\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}\mathrm{\Omega }{\mathit{\theta }}_x^{\prime }+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}{\mathit{\theta }}_y^{\prime }+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{{\mathit{\theta }}^{\prime }}_y^{\mathrm{3}}\\ & +\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right){\mathit{\theta }}_y+{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{\mathit{\theta }}_y^{\mathrm{3}}=\stackrel{\mathrm{‾}}{M}\sin \mathrm{\Omega }\stackrel{\mathrm{‾}}{t}.\end{array}\end{array}$$

Here the strokes denote derivatives in dimensionless time $\stackrel{\mathrm{‾}}{t}$.

Thus, it turns out that the steady-state motion of the considered rotor is described by the Duffing type differential equation system (13) and (14). To determine a periodic solution with a period equal to the period of external action, a method of decomposition of solutions (13) and (14) into Fourier series with uncertain coefficients is usually used. The coefficients can be found by the harmonic balance method (Hayashi, 1964; Szemplinska-Stupnicka, 1968; Kydyrbekuly, 2006), taking into account the finite and usually small number of terms.

3 Amplitude- and phase-frequency characteristic and moment of force transmissibility of the system

Expanding the solutions of the equations of motion (13) and (14) into the Fourier series with indefinite coefficients, can verify that the approximation of the solutions with a simple harmonic and oscillation frequency equal to the disturbing moment frequency is quite satisfactory in the case of the main resonance. Taking into account the following notation of the vibration parameters A=A₁ and α=α₁, the solutions of the equations of motio (13) and (14) at a first approximation can be written as

$$\begin{array}{}\text{(15)}& {\displaystyle }{\mathit{\theta }}_x=A\cos (\mathrm{\Omega }\stackrel{\mathrm{‾}}{\mathrm{t}}-\mathit{\alpha }),\text{(16)}& {\displaystyle }{\mathit{\theta }}_y=A\sin (\mathrm{\Omega }\stackrel{\mathrm{‾}}{\mathrm{t}}-\mathit{\alpha }).\end{array}$$

After substitution the expressions (15)–(16) in (13)–(14) and applying the harmonic balance method (Hayashi, 1964; Szemplinska-Stupnicka, 1968; Kydyrbekuly, 2006), obtain the amplitude-frequency and the phase- frequency dependences of the fundamental harmonic

$$\begin{array}{}\text{(17)}& {\displaystyle }\begin{array}{rl}& \left\{{\left[\left(\mathrm{1}-H\right){\mathrm{\Omega }}^{\mathrm{2}}-\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right)-\mathrm{0.75}{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{A}^{\mathrm{2}}\right]}^{\mathrm{2}}\right.\\ & \left.+{\left({\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}+\mathrm{0.75}{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{A}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{2}}\right)}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{2}}\right\}{A}^{\mathrm{2}}={\stackrel{\mathrm{‾}}{M}}^{\mathrm{2}},\end{array}\text{(18)}& {\displaystyle }\tan \mathit{\alpha }={\displaystyle \frac{\left({\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}+\mathrm{0.75}{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{A}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{2}}\right)\mathrm{\Omega }}{-(\mathrm{1}-H){\mathrm{\Omega }}^{\mathrm{2}}+({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P)+\mathrm{0.75}{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{A}^{\mathrm{2}}}},\end{array}$$

where $H={\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}-{\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}$ is the nominal thickness of the disc.

The formula for the amplitude-frequency characteristic (17) taking into account the expressions of dimensionless parameters (10), the critical velocity formula (11), and the formula for the resulting moment amplitude (12), may be needed when experimentally determining the values of the linear damping coefficients μ₁ and nonlinear cubic damping μ₃ of the support under the values of the amplitude A and rotation speeds Ω of the experimental resonance curve A=A(Ω).

The results of numerical solutions of the equations of motion (13) and (14) of a rotor system with nonlinear rigidity and the dependency graphs (15) and (16) taking into account expressions (17) and (18) are presented in Figures 2 and 3 for the following parameters: l=0.88, P=0.012, ${\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}=\mathrm{0.19},{\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}=\mathrm{0.1}$, H=0.09, ε=0.91, K₁=1.19, Ω=1.5.

Figure 2 Oscillogram θ_x=θ_x(t) and θ_y=θ_y(t) when ${\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}=\mathrm{0.075},{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}=\mathrm{0.09},{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}=\mathrm{0}$. Solid lines are the results of a numerical solution, points are the results of an analytical solution.

Figure 3 Oscillogram θ_x=θ_x(t) and θ_y=θ_y(t) when ${\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}=\mathrm{0.075},{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}=\mathrm{0.09},{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}=\mathrm{0.18}$. Solid lines are the results of a numerical solution, points are the results of an analytical solution.

The graphs in Figs. 2 and 3 show a sufficient approximation of the results of numerical modeling and analytically obtained results.

In the absence of nonlinear terms in Eqs. (13) and (14), the results for the linear rotor model are obtained from expressions (17) and (18) (Iskakov and Kalybaeva, 2010).

With the introduction of additional value notations

$$\begin{array}{}\text{(19)}& \begin{array}{rl}\mathrm{1}& =\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right)/\left(\mathrm{1}-H\right);K_{\mathrm{3}}={\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}/\left(\mathrm{1}-H\right);\\ {\mathit{\mu }}_{\mathrm{1}}& ={\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}/\left(\mathrm{1}-H\right),{\mathit{\mu }}_{\mathrm{3}}={\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}/\left(\mathrm{1}-H\right);\\ M& =\stackrel{\mathrm{‾}}{M}/\left(\mathrm{1}-H\right)\end{array}\end{array}$$

obtain an expression for the amplitude-frequency and the phase-frequency characteristics

$$\begin{array}{}\text{(20)}& {\displaystyle }{\left({\mathit{\omega }}_{\ast }^{\mathrm{2}}-{\mathrm{\Omega }}^{\mathrm{2}}\right)}^{\mathrm{2}}+{\left({\mathit{\mu }}_{\mathrm{1}}\mathrm{\Omega }+\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}\mathrm{\Omega }{}^{\mathrm{3}}{A}^{\mathrm{2}}\right)}^{\mathrm{2}}=M^{\mathrm{2}}/A^{\mathrm{2}},\text{(21)}& {\displaystyle }\tan \mathit{\alpha }={\displaystyle \frac{\left({\mathit{\mu }}_{\mathrm{1}}+\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}{A}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{2}}\right)\mathrm{\Omega }}{\mathrm{1}-{\mathrm{\Omega }}^{\mathrm{2}}+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{2}}}},\end{array}$$

where the notations are introduced

$$\begin{array}{}\text{(22)}& {\mathit{\omega }}_{\ast }^{\mathrm{2}}=\mathrm{1}+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{2}}.\end{array}$$

Amplitude-frequency supporting curves typically describe a relationship between amplitude and frequency of free vibrations of the system without damping (Panovko, 1971). Assuming the expression M to be equal to zero due to the external moment, and the damping coefficients μ₁ and μ₃ in Eq. (20), obtain the equation of the supporting curve for vibrations at the fundamental resonance frequency

$$\begin{array}{}\text{(23)}& A={\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{0.75}{K}_{\mathrm{3}}}}}\sqrt{{\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}}.\end{array}$$

Here

$$\mathrm{\Omega }\ge \mathrm{1}.$$

It can be seen from formula (23) that the supporting curve is a parabola symmetric about the axis Ω, the larger the valueK₃ the greater the slope to the right, of the supporting curve Formula (23) and the reference curve constructed from it are important for the experimental determination of the nonlinear stiffness coefficient of the support material K₃ from the intersection point of the reference curve with the resonance curve A=A(Ω).

In case when H=0.99 and moment P of gravitational force is neglected the amplitude of the centrifugal couple is M=Ω². Solving Eq. (20) with respect to A depending on Ω and μ₁ at a constant value of μ₃= 0 and depending on Ω and μ₃ at a constant value of μ₁= 0 with K₃=0, get

$$\begin{array}{}\text{(24)}& {\displaystyle }A={\displaystyle \frac{\mathrm{\Omega }}{\sqrt{{\left(\mathrm{1}-{\mathrm{\Omega }}^{\mathrm{2}}\right)}^{\mathrm{2}}+{\left({\mathit{\mu }}_{\mathrm{1}}\mathrm{\Omega }\right)}^{\mathrm{2}}}}},\text{(25)}& {\displaystyle }A=\sqrt{\sqrt[\mathrm{3}]{-q/\mathrm{2}+\sqrt{Q}}+\sqrt[\mathrm{3}]{-q/\mathrm{2}-\sqrt{Q}}},\end{array}$$

where

$$Q={\left(p/\mathrm{3}\right)}^{\mathrm{3}}+{\left(q/\mathrm{2}\right)}^{\mathrm{2}},p={\left({\displaystyle \frac{\mathrm{1}-{\mathrm{\Omega }}^{\mathrm{2}}}{\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}}}\right)}^{\mathrm{2}},q={\displaystyle \frac{\mathrm{1}}{{\left(\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}\mathrm{\Omega }\right)}^{\mathrm{2}}}}.$$

If shaft speeds are close to zero, the action of the gravitation moment is more noticeable. Then, if Ω=0, the Eq. (20) takes the following form

$$\begin{array}{}\text{(26)}& (\mathrm{1}+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{2}})A=P.\end{array}$$

Solutions of Eq. (26) if K₃=0:

$$\begin{array}{}\text{(27)}& A=P\end{array}$$

and if K₃≠0

$$\begin{array}{}\text{(28)}& \begin{array}{rl}A& =\sqrt[\mathrm{3}]{{\displaystyle \frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}}+\sqrt{{\left({\displaystyle \frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}}\right)}^{\mathrm{3}}+{\left({\displaystyle \frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}}\right)}^{\mathrm{2}}}}\\ & +\sqrt[\mathrm{3}]{{\displaystyle \frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}}-\sqrt{{\left({\displaystyle \frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}}\right)}^{\mathrm{3}}+{\left({\displaystyle \frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}}\right)}^{\mathrm{2}}}}.\end{array}\end{array}$$

Using dimensionless parameters in accordance with formulas (10), it is possible to give the following form to expressions of the projections of moment of transmitted force (transmissive force)

$$\begin{array}{}\text{(29)}& \begin{array}{rl}& {\stackrel{\mathrm{‾}}{M}}_{\mathit{\tau }x}={\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}{{\mathit{\theta }}^{\prime }}_x+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{{\mathit{\theta }}^{\prime }}_x^{\mathrm{3}}+\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right){\mathit{\theta }}_x+{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{\mathit{\theta }}_x^{\mathrm{3}},\\ & {\stackrel{\mathrm{‾}}{M}}_{\mathit{\tau }y}={\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}{{\mathit{\theta }}^{\prime }}_y+{\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}{{\mathit{\theta }}^{\prime }}_y^{\mathrm{3}}+\left({\stackrel{\mathrm{‾}}{K}}_{\mathrm{1}}{l}^{\mathrm{2}}-P\right){\mathit{\theta }}_y+{\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}{l}^{\mathrm{4}}{\mathit{\theta }}_y^{\mathrm{3}}.\end{array}\end{array}$$

On the other hand, take the projection of the moment of transmitted force (29) in a form of the following harmonics

$$\begin{array}{}\text{(30)}& \begin{array}{rl}& {\stackrel{\mathrm{‾}}{M}}_{\mathit{\tau }x}=A_{\mathit{\tau }}\cos (\mathrm{\Omega }\stackrel{\mathrm{‾}}{t}-\mathit{\alpha }),\\ & {\stackrel{\mathrm{‾}}{M}}_{\mathit{\tau }y}=A_{\mathit{\tau }}\sin (\mathrm{\Omega }\stackrel{\mathrm{‾}}{t}-\mathit{\alpha })\end{array}\end{array}$$

Substituting solutions (15), (16) and (30) into expressions (29) and taking into account additional notations of values (19), obtain an expression for determining amplitude of the moment of transmitted force A_τ depending on the rotation velocity Ω, the linear damping coefficient μ_1, and the nonlinear damping coefficient μ₃:

$$\begin{array}{}\text{(31)}& A_{\mathit{\tau }}^{\mathrm{2}}={\left({\mathit{\mu }}_{\mathrm{1}}A\mathrm{\Omega }+\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}{A}^{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}\right)}^{\mathrm{2}}+{\left(A+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{3}}\right)}^{\mathrm{2}}.\end{array}$$

The moment of force transmissibility is defined as a ratio of amplitude of the moment of transmitted force to the amplitude of the disturbing moment (Ho et al., 2012):

$$\begin{array}{}\text{(32)}& \begin{array}{rl}\mathit{\gamma }& =A_{\mathit{\tau }}/M=A_{\mathit{\tau }}/{\mathrm{\Omega }}^{\mathrm{2}}\\ & =\sqrt{\left[{\left({\mathit{\mu }}_{\mathrm{1}}A\mathrm{\Omega }+\mathrm{0.75}{\mathit{\mu }}_{\mathrm{3}}{A}^{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}\right)}^{\mathrm{2}}+{\left(A+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{3}}\right)}^{\mathrm{2}}\right]/{\mathrm{\Omega }}^{\mathrm{4}}}.\end{array}\end{array}$$

Introduce a function

$$\begin{array}{}\text{(33)}& \begin{array}{rl}f\left(K_{\mathrm{3}},{\mathit{\mu }}_{\mathrm{1}},{\mathit{\mu }}_{\mathrm{3}},\mathrm{\Omega }\right)& =\left[{\left({\mathit{\mu }}_{\mathrm{1}}A\mathrm{\Omega }+\mathrm{0},\mathrm{75}{\mathit{\mu }}_{\mathrm{3}}{A}^{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}\right)}^{\mathrm{2}}\right.\\ & \left.+{\left(A+\mathrm{0},\mathrm{75}{K}_{\mathrm{3}}{A}^{\mathrm{3}}\right)}^{\mathrm{2}}\right]/{\mathrm{\Omega }}^{\mathrm{4}}.\end{array}\end{array}$$

Calculate a derivative of f at Ω and equate it to zero:

$$\begin{array}{}\text{(34)}& {\displaystyle \frac{\mathrm{d}f}{\mathrm{d}\mathrm{\Omega }}}=\mathrm{0}.\end{array}$$

The solution of Eq. (34) is a certain value Ω₁. If K₃=0, this value is equal Ω₁≈1. It can be shown that the second derivative of f at Ω at this point is negative, i.e. the function f(Ω) really has the maximum value, while

$$\begin{array}{}\text{(35)}& {\mathit{\gamma }}_{max}=\sqrt{f\left({\mathrm{\Omega }}_{\mathrm{1}}\right)}.\end{array}$$

If the shaft speed is Ω=0 and K₃=0 the moment of force transmissibility of the system is

$$\begin{array}{}\text{(36)}& \mathit{\gamma }={\displaystyle \frac{A}{P}}=\mathrm{1}\end{array}$$

and if Ω=0 and K₃≠0

$$\begin{array}{}\text{(37)}& \begin{array}{rl}& \mathit{\gamma }={\displaystyle \frac{\mathrm{1}}{P}}\left(A+\mathrm{0.75}{K}_{\mathrm{3}}{A}^{\mathrm{3}}\right)\\ & ={\displaystyle \frac{\mathrm{1}}{P}}\left\{\begin{array}{l}\sqrt[\mathrm{3}]{\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}+\sqrt{{\left(\frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}\right)}^{\mathrm{3}}+{\left(\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}\right)}^{\mathrm{2}}}}\\ +\sqrt[\mathrm{3}]{\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}-\sqrt{{\left(\frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}\right)}^{\mathrm{3}}+{\left(\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}\right)}^{\mathrm{2}}}}+\\ +\mathrm{0.75}{K}_{\mathrm{3}}\left[\sqrt[\mathrm{3}]{\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}+\sqrt{{\left(\frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}\right)}^{\mathrm{3}}+{\left(\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}\right)}^{\mathrm{2}}}}\right.\\ {\left.+\sqrt[\mathrm{3}]{\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}-\sqrt{{\left(\frac{\mathrm{1}}{\mathrm{2.25}{K}_{\mathrm{3}}}\right)}^{\mathrm{3}}+{\left(\frac{P}{\mathrm{1.5}{K}_{\mathrm{3}}}\right)}^{\mathrm{2}}}}\right]}^{\mathrm{3}}\end{array}\right\}.\end{array}\end{array}$$

4 Impact of damping on rotor dynamics

The numerical solution of the Eq. (20) and the calculation from the expression (32) were performed using MathLab to demonstrate an impact of linear and nonlinear viscous damping on vibration isolation of the Duffing type rotor system. The results are shown in Figs. 4–13. The rotary system with linear stiffness of the elastic support as described in researches (Peng et al., 2012; Iskakov, 2018a, b) is presented in Figs. 4–7. In Figures 4 and 5, the elastic support of the rotor has linear stiffness and linear damping with μ₁=0.1, 0.2, 0.4. When the parameter of linear damping μ₁ increases, the amplitude of vibrations A and the moment of force transmissibility of the system γ in the resonance area near Ω≈1 decrease, and in non-resonance areas, where Ω≻1 or Ω≺1, the impacts of μ₁ on A and γ are very weak and negligible. In Figs. 6 and 7, the μ₁ value is kept constant at the level of 0.1, while the cubic damping parameter μ₃ takes 0, 0.2, 0.4. In this case, the resonance peak of the fundamental harmonic at Ω≈1 is significantly suppressed with increasing μ_3, and in the non-resonance area, where Ω≺1 μ₃ has virtually no impact on A and γ. However, in the non-resonance area, where Ω≻1, it can be seen that an increase in μ₃ may slightly suppress the A vibration amplitude and slightly increase γ.

Figure 4 Dependence of the vibration amplitude A on the rotor velocity Ω with linear stiffness of the support and at different values of the coefficient of linear viscous damping μ₁.

Figure 5 Dependence of the moment of force transmissibility of the system γ on the rotor velocity Ω with linear stiffness of the support and at different values of the coefficient of linear viscous damping μ₁.

Figures 8–9 shows the effects of linear damping in case of nonlinear stiffness of the elastic support. A contrast between the effects of linear and nonlinear viscous damping can be observed by comparing Figs. 8–9 with Figs. 10–13, where the cubic viscous damping coefficient μ₃ is 0, 0.2, 0.4 for K₃=0.05 and K₃=0.1, respectively. Similar to the case of pure linear stiffness, the linear viscous damping μ₁ again modifies the resonance area without causing noticeable effects in non-resonance areas.

To ensure the stable movement, avoiding jumping should be an important feature of the vibration isolator. Researches of Ravindra and Mallick (1994) have shown that the jumping phenomena can be eliminated by linear damping. With pure linear viscous damping, the calculation results in Figs. 8 and 9 show jumps occurring at Ω=1.54 and Ω=1.52 when μ₁=0.1 and μ₁=0.2 for K₃=0.1 and μ₃=0, but when the level of linear damping increases to μ₁=0.4, jumps no longer exist. Nonlinear viscous damping can also take into account the known jumping phenomenon in the presence of nonlinear stiffness of the elastic support of the rotor. In Figs. 10–13, jumps disappear when μ₃ increases from 0 to 0.2. A compromise between the use of linear damping to eliminate the jumping phenomena is in the overwhelming impact of nonlinear damping on the vibration amplitude in high velocity areas.

Figure 6 Dependence of the vibration amplitude A on the rotor velocity Ω with linear stiffness of the support and at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Figure 7 Dependence of the moment of force transmissibility of the system γ on the rotor velocity Ω with linear stiffness of the support and at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Figure 8 Dependence of the vibration amplitude A on the rotor velocity Ω with nonlinear stiffness coefficient of the support K₃=0.1 and at different values of the coefficient of linear viscous damping μ₁.

Figure 9 Dependence of the moment of force transmissibility of the system γ on the rotor velocity Ω with nonlinear stiffness of the support K₃=0.1 and at different values of the coefficient of linear viscous damping μ₁.

In non-resonance area, where Ω≺1, linear damping μ₁ and cubic nonlinear damping μ₃ have virtually no impact on the ability to the isolator absolute displacement. Derivatives $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{1}}}\approx \mathrm{0}$ and $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{3}}}\approx \mathrm{0},$ moment of force transferability of the system γ weakly depends on μ₁ and μ₃. In the resonance area, where Ω≈1, an increase in either linear damping μ₁ or cubic nonlinear damping μ₃ can significantly suppress ability to absolute displacement. Derivatives $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{1}}}\prec \mathrm{0}$ and $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{3}}}\prec \mathrm{0}$, with an increase in μ₁ and μ₃, the moment of force transmissibility of the system γ decreases. In high velocity area of, where Ω≻1, an increase in the linear damping μ₁ does not affect the ability to the isolator absolute displacement. Derivative $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{1}}}\approx \mathrm{0},$ the moment of force transmissibility of the system γ remains almost unchanged, an increase in the cubic nonlinear damping μ₃ may slightly reduce the ability to relative displacement and an increase in ability to absolute displacement and the moment of force transmissibility of the system γ, i.e $\frac{\mathrm{d}\mathit{\gamma }}{\mathrm{d}{\mathit{\mu }}_{\mathrm{3}}}\succ \mathrm{0}$.

Figure 10 Dependence of the vibration amplitude A on the rotor velocity Ω with nonlinear stiffness coefficient of the support K₃=0.05 and at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Figure 11 Dependence of the moment of force transmissibility of the system γ on the rotor velocity Ω with nonlinear stiffness coefficient of the support K₃=0.05 at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Figure 12 Dependence of the vibration amplitude A on the rotor velocity Ω with nonlinear stiffness coefficient of the support K₃=0.1 at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Figure 13 Dependence of the moment of force transmissibility of the system γ on the rotor velocity Ω with nonlinear stiffness coefficient of the support K₃=0.1 at value of the coefficient of linear viscous damping μ₁=0.1 and various values of nonlinear viscous damping coefficient μ₃.

Diagrams for study of the phase-frequency characteristic of the rotor are shown in Figs. 14–17. When plotting the diagrams, the properties of the arctangent function were used. These diagrams show that in case of linear stiffness of the support, a rotation of the vibration phase occurs when passing through the critical velocity, i.e. when Ω=1 (Figs. 14–17), and in case of a nonlinear component of stiffness of the support, a rotation of the vibration phase is observed in a phenomenon of the amplitude jump, in Figs. 16–17 this corresponds to the shaft velocity Ω=1.6. The joint impact of the linear and nonlinear components of the support damping essentially eliminates the jump phenomenon, which indicates the rotation of the vibration phase when the shaft velocity is close to the critical one, i.e. if Ω≈1 (Fig. 17), as in the case of only linear stiffness of the support (Figs. 14–15).

Figure 14 Dependence of the vibration phase α on the rotor speed Ω with linear stiffness of the support and at different values of the coefficient of linear viscous damping μ₁.

Figure 15 Dependence of the vibration phase α on the rotor speed Ω with linear stiffness of the support and at value of the coefficient of linear viscous damping μ₁=0.1 and different values of the nonlinear viscous damping coefficient μ₃.

Figure 16 Dependence of the vibration phase α on the rotor speed Ω with coefficient of non-linear stiffness of the support K₃=0.1 and at different values of the coefficient of linear viscous damping μ₁.

Figure 17 Dependence of the vibration phase α on the rotor speed Ω with coefficient of non-linear stiffness of the support K₃=0.1 and at value of the coefficient of linear viscous damping μ₁=0.1 and different values of the nonlinear viscous damping coefficient μ₃.

The results of numerical solutions of the equations of motio (13) and (14) of a rotor system with nonlinear rigidity with the output of the A=A(Ω) dependence for the parameters: l=0.88, P=0.012, ${\stackrel{\mathrm{‾}}{I}}_p=\mathrm{0.19}$, ${\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}=\mathrm{0.1}$, H=0.09, ε=0.91, K₁=1.19 are shown in Figs. 18 and 19. A comparison of the graphs in Figs. 18 and 19 shows that nonlinear damping of the elastic support significantly suppresses the main resonant amplitude at Ω≈1.57 and completely eliminates the height of the vibration amplitude in the region of the rotation speed Ω∼1.178.

Figure 18 Amplitude response of oscillations at ${\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}=\mathrm{0.075}$, ${\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}=\mathrm{0.09}$, ${\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}=\mathrm{0}$.

Figure 19 Amplitude response of oscillations at ${\stackrel{\mathrm{‾}}{K}}_{\mathrm{3}}=\mathrm{0.075}$, ${\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{1}}=\mathrm{0.09}$, ${\stackrel{\mathrm{‾}}{\mathit{\mu }}}_{\mathrm{3}}=\mathrm{0},\mathrm{18}$.

5 Motion stability

Consider motion stability using an approximate theory (Iskakov, 2018b; Iskakov, 2019; Van Dooren, 1971). The geometric location of the points at which the amplitude curves for the oscillations of the principal resonance have vertical tangents is determined by the equation

$$\begin{array}{}\text{(38)}& {\displaystyle \frac{\mathrm{d}\mathrm{\Omega }}{\mathrm{d}A}}=\mathrm{0}.\end{array}$$

In accordance with Eq. (20), equality takes place

$$\begin{array}{}\text{(39)}& \begin{array}{rl}{f}_{\mathrm{1}}(\mathrm{\Omega },a)& =\left[{\left(\mathrm{1}+\mathrm{0},\mathrm{75}{K}_{\mathrm{3}}{A}^{\mathrm{2}}-{\mathrm{\Omega }}^{\mathrm{2}}\right)}^{\mathrm{2}}\right.\\ & +\left.{\left({\mathit{\mu }}_{\mathrm{1}}\mathrm{\Omega }+\mathrm{0},\mathrm{75}{\mathit{\mu }}_{\mathrm{3}}\mathrm{\Omega }{}^{\mathrm{3}}{A}^{\mathrm{2}}\right)}^{\mathrm{2}}\right]A^{\mathrm{2}}-{\mathrm{\Omega }}^{\mathrm{2}}=\mathrm{0}.\end{array}\end{array}$$

Differentiated the last Eq. (39) in frequency, obtaining

$$\begin{array}{}\text{(40)}& {\displaystyle \frac{\partial f_{\mathrm{1}}}{\partial \mathrm{\Omega }}}+{\displaystyle \frac{\partial f_{\mathrm{1}}}{\partial A}}{\displaystyle \frac{\partial A}{\partial \mathrm{\Omega }}}=\mathrm{0}.\end{array}$$

It follows that

$$\begin{array}{}\text{(41)}& {\displaystyle \frac{\partial \mathrm{\Omega }}{\partial A}}=-{\displaystyle \frac{\partial f_{\mathrm{1}}/\partial A}{\partial f_{\mathrm{1}}/\partial \mathrm{\Omega }}}.\end{array}$$

Condition (36) is approximately satisfied subject to the equality

$$\begin{array}{}\text{(42)}& {\displaystyle \frac{\partial f_{\mathrm{1}}}{\partial A}}=\mathrm{0}.\end{array}$$

Equality (42) with allowance for (39) leads to the equation

$$\begin{array}{}\text{(43)}& \begin{array}{rl}& \left(\mathrm{1}+{\displaystyle \frac{\mathrm{3}}{\mathrm{4}}}{K}_{\mathrm{3}}{A}^{\mathrm{2}}-{\mathrm{\Omega }}^{\mathrm{2}}\right)\left(\mathrm{1}+{\displaystyle \frac{\mathrm{9}}{\mathrm{4}}}{K}_{\mathrm{3}}{A}^{\mathrm{2}}-{\mathrm{\Omega }}^{\mathrm{2}}\right)\\ & +\left({\mathit{\mu }}_{\mathrm{1}}\mathrm{\Omega }+{\displaystyle \frac{\mathrm{3}}{\mathrm{4}}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}{A}^{\mathrm{2}}\right)\left({\mathit{\mu }}_{\mathrm{1}}\mathrm{\Omega }+{\displaystyle \frac{\mathrm{9}}{\mathrm{4}}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{3}}{A}^{\mathrm{2}}\right)=\mathrm{0}.\end{array}\end{array}$$

Equation (43) describes the boundary curves of the stability region for oscillations at the principal resonance frequency. Therefore, the geometrical location of the points at which the amplitude curves have vertical tangents determine the boundaries of the stability region (Iskakov, 2018b, 2019; Van Dooren, 1971).

The solutions of Eq. (43) have the form

$$\begin{array}{}\text{(44)}& A_{\mathrm{1},\mathrm{2}}={\left\{\begin{array}{l}\frac{\mathrm{8}}{\mathrm{9}}\frac{K_{\mathrm{3}}\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)-{\mathit{\mu }}_{\mathrm{1}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{4}}}{K_{\mathrm{3}}^{\mathrm{2}}+{\mathit{\mu }}_{\mathrm{3}}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{6}}}\\ \mp {\left\{\begin{array}{l}{\left[\frac{\mathrm{8}}{\mathrm{9}}\frac{K_{\mathrm{3}}\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)-{\mathit{\mu }}_{\mathrm{1}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{4}}}{K_{\mathrm{3}}^{\mathrm{2}}+{\mathit{\mu }}_{\mathrm{3}}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{6}}}\right]}^{\mathrm{2}}\\ -\frac{\mathrm{16}}{\mathrm{27}}\frac{{\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)}^{\mathrm{2}}+{\left(\mathit{\mu }{\mathrm{\Omega }}_{\mathrm{1}}\right)}^{\mathrm{2}}}{K_{\mathrm{3}}^{\mathrm{2}}+{\mathit{\mu }}_{\mathrm{3}}^{\mathrm{2}}{\mathrm{\Omega }}^{\mathrm{6}}}\end{array}\right\}}^{\frac{\mathrm{1}}{\mathrm{2}}}\end{array}\right\}}^{\frac{\mathrm{1}}{\mathrm{2}}}.\end{array}$$

Equation (43) have real roots (44) under the condition

$$\begin{array}{rl}& {\left\{{\displaystyle \frac{\mathrm{8}}{\mathrm{9}}}\left[K_{\mathrm{3}}\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)-{\mathit{\mu }}_{\mathrm{1}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{4}}\right]\right\}}^{\mathrm{2}}\\ & \succ {\displaystyle \frac{\mathrm{16}}{\mathrm{27}}}\left[{\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)}^{\mathrm{2}}+{\left(\mathit{\mu }{\mathrm{\Omega }}_{\mathrm{1}}\right)}^{\mathrm{2}}\right],\end{array}$$

where

$$K_{\mathrm{3}}\left({\mathrm{\Omega }}^{\mathrm{2}}-\mathrm{1}\right)\succ {\mathit{\mu }}_{\mathrm{1}}{\mathit{\mu }}_{\mathrm{3}}{\mathrm{\Omega }}^{\mathrm{4}}.$$

The results of calculations by the formula (44) for a rotor system with the support of cubic non-linear stiffness under the influence of linear and non-linear damping of the support are shown in Figs. 20 to 21. Figure 20 were given the boundaries of the instability region of an oscillatory system with non-linear stiffness K₃=0.2 depending on the rotation velocity Ω for different values of the linear damping coefficient μ₁=0.1; 0.2; 0.4. A slight narrowing of the height of the instability region is observed as the linear damping coefficient increases from 0.1 to 0.4.

Figure 21 shows the boundaries of the instability region of an oscillatory system with non-linear stiffness K₃=0.2 depending on the rotation velocity Ω for different values of the non-linear damping coefficient μ₃=0; 0.01; 0.02. Here it was take that μ₁=0.1. This shows that with an increase the coefficient of nonlinear damping μ₃ from 0 to 0.02, not only the height of the instability region decreases, but also its width.

Figure 20 The boundaries of the instability region for a system with a coefficient of non-linear stiffness of the support K₃=0.2 and different values of linear viscous damping coefficient μ₁.

Figure 21 The boundaries of the instability region for a system with a coefficient of non-linear stiffness of the support K₃=0.2 with a linear viscous damping coefficient μ₁=0.1 and different values of non-linear viscous damping coefficient μ₃.

6 The dynamics of the experimental rotor

Experimental studies of an effect of the support material on dynamics of the centrifugal gyroscopic rotor were carried out on a test unit which is generally shown in Fig. 22, and a flow chart with all measuring instruments is shown in Fig. 23.

Figure 22 Test Unit Overview: 1 – centrifuge with rotor, 2 – F 5034 frequency meter, 3 – H 043.1 electromagnetic oscillograph, 4 – Topaz-3 DC amplifier, 5 – amplifier UT-4 AC amplifier, 6 – Agate power supply.

Figure 23 Experimental installation flowchart.

The unit consists of an electrically powered rotor, a channel for measuring and recording the precession motion performance, a channel for measuring and recording the rotor shaft speed, a channel for measuring and recording current in driven electric motor.

The rotor to be studied is made in form of a cylindrical cup made of duralumin D 16 T with a top transparent PMMA cover. Distance between supports l₀=0.33 m, distance from the bottom support to the cup center L=0.52 m, eccentricity $e=\mathrm{0.2}\times {\mathrm{10}}^{-\mathrm{3}}$ m, cup mass m=3 kg, cup weight G=29.4 N, the rubber support rigidity coefficient $k_{\mathrm{1}}=\mathrm{2}\times {\mathrm{10}}^{\mathrm{4}}$ N m⁻¹. The rotor is mounted on the upper cone part of the shaft inside the centrifuge (1).

The unit uses an electromagnetic speed sensor, voltage pulses from which are fed to an input of frequency meter F 5034 (2) and/or to one of the channels of the electromagnetic oscillograph H 043.1 (3).

Strain-gauge and variable inductance transducers of displacement were used to record precession motion performance. Transducers are coupled by a half-bridge circuit and are connected to an input of Topaz-3 (4) amplifier. Offset voltage of the bridge between a movable rotor and variable inductance transducer of displacement is fed to the input of UT-4 (5) amplifier and from its output to H 043.1 (3) oscilloscope. Supply voltage of Topaz-3 and UT-4 is produced by Agate (6) power supply unit.

To measure current intensity in circuit of the driving motor, compensation method was used since the capacity variations to be determined are much less than the primary engine capacity.

Dynamic performance of the rotor motion and parameters of interaction with the driving motor were read according to the following procedure.

According to the results of measuring the amplitude of signals received in H 043.1 (3) oscilloscope, from the strain resistor and inductive sensors through Topaz-3 (4) and UT-4 amplifiers (5), respectively, and according to the readings of F 5034 (2) frequency meter, the amplitude-frequency performance of oscillations of the experimental rotor was read.

Sensitivity and stability of amplifier channels was monitored and precession performance was recorded before and after each experiment by displacing the rotor shaft from an equilibrium position towards the sensor using a pusher with a distance meter. Value of voltage applied to the motor from LATR-1M laboratory autotransformer was recorded at all shaft velocities.

Along with determination of the rotor oscillation amplitude and the shaft velocity, the current in the circuit of the motor was recorded by H 043.1 oscilloscope.

The above operations were repeated at different damping of upper elastic support and values of radial clearance of bearing. Changing in damping was achieved by installing the PTFE support material in one case, and the soft rubber support material in another one. Changing in radial clearance was carried out by changing the bearings. Dynamic values of the rotor motion were compared for two limiting options of these parameters – a fluoroplastic support and bearings with 0.1 mm clearance in one case and a soft rubber elastic support and bearings with 0.01 mm clearance in another one.

As can be seen from the above, the measured values were amplitude and frequency of the rotor precession, shaft speed, voltage and current in the driving motor circuit.

Figure 24 Amplitude-frequency characteristic of the centrifugal rotor in different versions of the support material.

Amplitude-frequency characteristics of the rotor measured for cases of fluoroplastic support, bearings with a clearance of 0.1 mm and a rubber support, bearings with a clearance of 0.01 mm are shown in Fig. 24. Similar curves were previously obtained respectively in the papers by Kelzon and Pryadilov (1965), Tuleshov (1987). Comparison of the amplitude-frequency characteristic curves shows that the replacement of the fluoroplastic material of the support with rubber damps a resonance curve of the main harmonic in the area of the critical velocity ω≈68 s⁻¹ and in the high-velocity area up to the rotation velocity ω=284 s⁻¹, completely suppresses a significant elevation of the vibration amplitude in the area of the rotation velocity ω∼600 s⁻¹. It describes a previously undetected peak of much smaller amplitude (relative to the above one) in the area of the rotation velocity of ω∼300 s⁻¹.

Zones of higher vibration amplitudes are not related to a coincidence of the shaft velocity with one of the rotor eigenfrequencies, but are related to the ripple frequency of power in the electric motor because of the limited energy source power. The presence of modes of coupled vibrations of power and precession amplitude (Felix et al., 2015; Tuleshov, 1987), at frequencies that are multiples of the frequency of the alternating component of the actuating moment (628 s⁻¹) can be, along with the nonlinearity of the system itself, also due to the presence of subharmonics in the time law of the supply voltage change. The results of the experimental works qualitatively correspond to the results of theoretical studies of the impact of nonlinear damping on the amplitude-frequency characteristic of the rotor, carried out in the absence of nonlinear stiffness of the elastic support.

Figure 25 presents a dependence of the angular shaft velocity on the tension applied during acceleration and braking of the engine. The obtained dependence indicates the quadratic character of an increase in losses in the engine with an increase in its power, which is consistent with (Vishnevsky, 1955). Significant deviation from this law when accelerating the rotor in the range of angular velocities 550 s${}^{-\mathrm{1}}\prec \mathit{\omega }\prec \mathrm{700}$ s⁻¹ and when braking in the range 400 s${}^{-\mathrm{1}}\prec \mathit{\omega }\prec \mathrm{600}$ s⁻¹ indicates the development of a process of consuming additional power from the engine (Kelzon and Pryadilov, 1965). This explains the significant elevation of the vibration amplitude in the range of angular velocities 550 s${}^{-\mathrm{1}}\prec \mathit{\omega }\prec \mathrm{700}$ s⁻¹ (Fig. 24).

Swinging of the centrifugal rotor by torsional vibrations occurs due to the cyclicity of the instantaneous power in the engine when powered from the AC network, i.e., there is a phenomenon of capture on a ripple frequency of the engine power. When replacing the support material with rubber, the movement of an area of increased amplitudes to the area of the shaft velocity value ω∼300 s⁻¹ is associated with the development of parametric vibrations in which the most noticeable response is at a frequency twice as low as a frequency of the exciting force (Tondl, 1971).

Figure 25 Dependence of the shaft velocity on the tension applied to the engine in case of the fluoroplastic material of the rotor support.

Thus, selecting the values of the stiffness coefficient k, the linear damping coefficient μ_d1 and the coefficient of nonlinear cubic damping μ_d3 of the support, it is possible to create such an elastic support that would have simultaneously a restoring and damping impact on the resonance amplitude and enable the rotor to pass painlessly through the critical velocity. The restoring and damping properties of the support are the principle of operation of the centrifuge elastic support, which operates based on a gyroscopic rotor, for which the Kazakhstan patent of invention was obtained (Iskakov and Kunelbayev, 2018).

The invention is illustrated by drawings in Fig. 26. In Fig. 26, the driving electric motor (4) with the power take-off shaft (5) is strengthened on the platform (3) by means of the cylindrical ring (1) and the brackets (2). The electric motor shaft is connected to the power take-off shaft by a transient coupling (6). The driving electric motor with a power take-off shaft is attached to the frame by a lower hinge (7) and an upper spring (8) and rubber (9) elastic supports. A cylindrical vessel (10) is mounted on the upper conical part of the shaft (5) and is held on the shaft (5) by a lock nut (11). The cylindrical vessel has an internal lining (12) and is closed on top with a transparent cover (13). To measure characteristics of the precessional movement, the displacement sensors (14) are installed on the bracket (2), connected to the control device by wires (15), (16) and (17). To measure the shaft velocity (5), a magnetoelectric velocity sensor (18) is used, connected by wires (19) and (20) to the recording device in the control apparatus. The rotor velocity in the operating modes is set by the power source voltage regulator.

Figure 26 Gyroscopic rotor based centrifuge.

The electric motor is powered by wires (21) and (22). The shaft precession amplitude for the operating velocity and at full fill with liquid is regulated by tightening the spring of the upper elastic support (8) with the help of a nut (23), depending on the production process and product quality requirements. In this case, a rubber nozzle of the lower elastic support (9) is selected in such a way that the value of the vibration isolation coefficient (moment of force transmissibility of the system) of a special resin or rubber material is acceptable for safe passage through the critical velocity in case of determining the operating velocity beyond the critical velocity.

7 Conclusions

The following conclusions can be drawn based on the analysis and discussion of the research results:

The impact of an elastic bearing material with linear viscous damping and with cubic non-linear viscous damping on non-linear vibrations of a vertical hard gyroscopic rotor was studied, taking into account the non-linear stiffness of the elastic support.

The equations of motion of the rotor were composed based on the Lagrange's Equations of Second Kind, which were reduced to a dimensionless form.

The equations of motion of the rotor were solved by the harmonic balance method, and the expressions of amplitude-frequency and phase-frequency characteristics of the fundamental harmonic and the expression for determining the moment of force transmissibility of the system were obtained.

Researches of the amplitude-frequency characteristic depending on the coefficient of linear viscous damping and nonlinear viscous cubic damping of the elastic support showed that both linear and nonlinear cubic damping significantly suppresses the resonance amplitude of the fundamental harmonic, eliminate jumping effect of the nonlinear system.

In non-resonance areas, an impact of linear damping on the vibration amplitude is very weak and negligible, and nonlinear cubic damping in an area where the velocity is many times greater than its critical value can slightly suppress the amplitude of the rotor vibrations and, therefore, only nonlinear damping can maintain performance of the vibration isolator throughout the entire range of the rotor velocity.

The research results were used when preparing an application for a patent for invention of a gyroscopic rotor based centrifuge to develop a description, abstract and claims. The invented centrifuge can be used in pharmaceutical and food industries to effectively intensify the mixing of suspensions of medicinal herbs, dairy and cultured milk products, other liquid food products throughout the volume of a container due to controlled precession.

Appendix A:  

Nomenclature
A= vibration amplitude, rad
Aτ= amplitude of the moment of transmitted force, dimensionless
e= linear eccentricity, m
G= disc weight, N
H= nominal thickness of the disc, dimensionless
IP= polar moment of inertia, kg m2
${\stackrel{\mathrm{‾}}{I}}_{\mathrm{P}}=$ polar moment of inertia, dimensionless
IT= transverse moment of inertia, kg m2
${\stackrel{\mathrm{‾}}{I}}_{\mathrm{T}}=$ transverse moment of inertia, dimensionless
k1= coefficient of linear stiffness, N m−1
k3= coefficient of nonlinear stiffness, N m−3
K3= coefficient of nonlinear stiffness, dimensionless
l0= distance between supports, m
L= shaft length, m
m= disc mass, kg
M= amplitude of the resultant moment, dimensionless
α= vibration phase, ∘
γ = moment of force transmissibility, dimensionless
${\mathit{\theta }}_x,{\mathit{\theta }}_y=$ angles, rad
μd1= coefficient of linear viscous damping, N ms rad−1
μ1= coefficient of linear viscous damping, dimensionless
μd3= coefficient of nonlinear cubic viscous damping, N ms3 rad−3
μ3= coefficient of nonlinear cubic viscous damping, dimensionless
ω= shaft speed, rad s−1
Ω= shaft speed, dimensionless
ω0= critical speed, rad s−1