Vibration and Instability of Rotating Composite Thin-Walled Shafts with Internal Damping

The dynamical analysis of a rotating thin-walled composite shaft with internal damping is carried out analytically.The equations of motion are derived using the thin-walled composite beam theory and the principle of virtual work. The internal damping of shafts is introduced by adopting the multiscale damping analysis method. Galerkin’s method is used to discretize and solve the governing equations. Numerical study shows the effect of design parameters on the natural frequencies, critical rotating speeds, and instability thresholds of shafts.


Introduction
The rotating composite material shafts are being used as structural elements in many application areas involving the rotating machinery systems.This is likely to contribute to the high strength to weight ratio, lower vibration level, and longer service life of composite materials.A significant weight saving can be achieved by the use of composite materials.Also by appropriate design of the composite layup configuration, orientation, and number of plies the improved performance of the shaft system can be obtained.Furthermore, the use of composite would permit the use of longer shafts in the supercritical range than what is possible with conventional metallic shafts.In the last few years, there exist numerous researches related to predicting critical speeds and natural frequencies of composite shaft.Zinberg and Symonds [1] investigated the critical speeds of rotating anisotropic cylindrical shafts based on an equivalent modulus beam theory (EMBT), and dos Reis et al. [2] evaluated the shaft of Zinberg and Symonds [1] by the finite element method.Kim and Bert [3] adopted the thin-and thick-shell theories of firstorder approximation to derive the motion equations of the rotating composite thin-walled shafts.They used this model to obtain a closed form solution for a simply supported drive shaft and to analyze the critical speeds of composite shafts.Singh and Gupta [4] developed two composite spinning shaft models employing EMBT and layerwise beam theory (LBT), respectively.It was shown that a discrepancy exists between the critical speeds obtained from both models for the unsymmetric laminated composite shaft.Chang et al. [5] presented a simple spinning composite shaft model based on a first-order shear deformable beam theory.The finite element method is used here to find the approximate solution of the system.The model was used to analyze the critical speeds, frequencies, mode shapes, and transient response of a particular composite shaft system.Gubran and Gupta [6] presented a modified EMBT model to account for the effects of a stacking sequence and different coupling mechanisms.Song et al. [7] used Rehfield's thin-walled beam theory [8] that presented a composite thin-walled shaft model.The effects of rotatory inertias, axial edge load, and boundary conditions on the natural frequencies and stability of the system were investigated.Ren et al. [9] proposed another composite thin-walled shaft model by means of the composite thin-walled beam theory, an asymptotically correct theory referred to as variational asymptotically method (VAM).The flexible composite shaft is assumed to support on bearings which are modeled as springs and dampers and containing of the rigid disks mounted on it.The natural frequencies and critical rotating speeds of the rotating composite shaft with the variation of the lamination angle, ratios of length over radius, ratios of radius over thickness, and shear deformation are then analyzed.
On the other hand, composite material shafts have higher internal damping than conventional metallic shafts.However, as it has been shown in previous research [10] that internal damping in rotating assemblies may lead to whirl instabilities in high speed rotors.Therefore, accurate prediction of effects of internal damping in composite material rotors is essential.So far, there have been less work related to the stability analysis, particularly effects of internal damping on the stability of a rotating composite shaft [4,[11][12][13][14].Singh and Gupta [4,11] introduced discrete viscous damping coefficients to account for effect of internal damping.In a similar approach, dynamic instability analysis of composite shaft has been performed by Montagnier and Hochard [13], Mazzei and Scott [12], Kim et al. [14].In above cases, internal damping terms were included simply in equations of motion of rotating shaft, no internal damping modeling has been described.Sino et al. [15] investigated the stability of an internally damped rotating composite shaft.Internal damping was introduced by the complex constitutive relation of a viscoelastic composite.The shaft was modeled by finite element method.However, only the mechanical coupling effects induced by symmetrical stacking are taken into account in their mode.
In the present work, an analytical model applicable to the dynamical analysis of rotating composite thin-walled shafts with internal damping is proposed.This model is based on the composite thin-walled beam theory, referred to as variational asymptotically method (VAM) by Berdichevsky et al. [16].The internal damping of shaft is introduced via the multiscale damping mechanics [17].The equations of motion of the composite shafts are derived by the principle of virtual work.Galerkin's method is used then to discretize and solve the governing equations.The natural frequencies, critical rotating speeds, and instability thresholds are obtained through numerical simulations.The effect of the ply angle and aspect ratio of cross-section are then assessed.The validity of the model is proved by comparing the results with those in literature and convergence examination.

Model of Composite Shaft
Consider the slender thin-walled composite shaft given in Figure 1.The length of the shaft is denoted by , its thickness is denoted by ℎ, the radius of curvature of the middle is denoted by , and the maximum cross-sectional dimension is denoted by .It is assumed that  ≪ , ℎ ≪ , and the coordinate  is measured along the normal to the middle surface within the limits −ℎ/2 ⩽  ⩽ ℎ/2.The shaft rotates about its longitudinal -axis at a constant rate Ω.

The Displacement and Strain.
The components of the displacement along the Cartesian coordinate (, , ) are expressed as follows [16]: where the primes in (1) denote differentiation with respect to .The tangential and normal displacements V 2 and V  can be expressed as follows: ( Based on the expressions shown in ( 1) and ( 2), the inplane strain components can be written in terms of the displacement variables as follows: (3)

Equations of Motion.
The equations of motion of the shaft free vibration can be described by a variational form where   ,   , and   are the variation of the strain energy, the kinetic energy, and the dissipated energy of the crosssection, respectively.The variation of the strain energy of the cross-section can be expressed as follows: where { , } T = { 11  22 2 12 },  is the cross-sectional area of the shaft, and [  ] is equivalent off-axis stiffness matrix.The variation of the dissipated energy of the cross-section can be expressed as follows [17]: where [  ] is equivalent off-axis damping matrix.

Shock and Vibration 3
Further, the variation of the kinetic energy of the crosssection with rotating motion is where [diag()] is a diagonal matrix with components equal to the mass density  of a ply.The position, velocity, and acceleration vectors for the deformed shaft are described as follows:

Equivalent Cross-Section Stiffness Matrix.
For the case of no internal pressure acting on the shaft, ( 5) can be simplified by using free hoop stress resultant ( 22 = 0) assumption as follows: where ∮(⋅) denotes the integral around the loop of the midline cross-section, and the reduced axial, coupling, and shear stiffness , , and  can be written as follows: In (9),   can also be expressed with respect to , V, , and  by combining (1)-( 3) and ( 9).Thus, one has where {Δ} is 4 × 1 column matrix of kinematic variables defined as {Δ} T = (      V  ) and [K] is 4 × 4 symmetric stiffness matrix.Its components   are given as follows: 2.4.Equivalent Cross-Section Damping Matrix.Similar to the derivation of the previous cross-section stiffness formulations, the variation of the dissipated energy of the crosssection in terms of the strains  11 and  12 can be modeled as follows: where The variation of the dissipated energy can be also expressed in terms of the kinematic variables as follows: where [C] is 4 × 4 symmetric damping matrix.The formulation of its components   is analogous to stiffness components   as shown in (12), but the terms , , and  in (10) should be replaced by the terms   ,   , and   , respectively.
2.5.Equivalent Cross-Section Mass.Substituting (1) into (8) and in view of ( 7), the variation of the kinetic energy of the cross-section can be obtained as follows: where 2.6.Approximate Solution Method.In order to find the approximate solution of the rotating composite shaft, the quantities (, ), V(, ), (, ), and (, ) are assumed in the form where   (),   (), and   ()are mode shape functions which fulfill all the boundary conditions of the composite shaft;  is complex eigenvalues of the system;   ,   ,   , and   are undetermined constants; and  = √ −1.Substituting ( 18) into the governing equations of motion equations ( 4)-( 6) and applying Galerkin's procedure, the following governing equations in matrix form can be found: where {} T = ( 1 ,  2 , . . .,   ,  1 ,  2 , . . .,   ,  1 ,  2 , . . .,   ,  1 ,  2 , . . .,   ) is a constant vector, [M] is the mass matrix; [G] is the gyroscopic matrix, [C] is the damping matrix, and [K] is the stiffness matrix which also includes contribution from the centrifugal forces.The detailed expressions of these matrices are as follows: where From ( 19), one can obtain the following complex eigenvalue problem: Complex eigenvalue  can be expressed in the form The damping natural frequency or whirl frequency of the system is the imaginary part , whereas its real part  gives the decay or growth of the amplitude of vibration.A negative value of  indicates a stable motion, whereas a positive value indicates an unstable motion, growing exponentially in time.

Numerical Results
The numerical calculations are performed by considering the shaft made of graphite-epoxy whose elastic characteristics are listed in Table 1.The shaft has rectangular cross-section of fixed geometrical characteristics width  = 0.32 m, length  = 4.5952 m, and wall thickness ℎ = 0.01016 m, whose layup is [] 16 with clamped-free boundary conditions.In order to examine the influence of the number of mode shape functions used in the solution of the equation on the accuracy of the results, the numerical results of natural frequency are shown in Tables 2 and 3 and modal damping in Tables 4 and 5 for an increasing number of mode shape functions, where , , and  are the length, width, and height, respectively.From these tables, it can be seen that to obtain accurate results of the first two natural frequencies and dampings, no more than five mode shape functions are required.This indicates clearly that the convergence of the present model is quite good.
A comparison of predictions using the present model with those obtained in [18] is also shown in Tables 2, 3, 4, and 5.A perfect agreement of numerical results with those in [18] can be seen.
Figure 2 shows the variation of the first two flexural natural frequencies versus rotating speed for various ply angles.As it can be seen, because of the nonsymmetry of the shaft crosssection (/ ̸ = 1), the standstill flexural frequencies about the two principal axes (flapping and sweeping denote bending about the and -axis, resp.) are unequal.The behaviors of the flapping and sweeping bending frequencies versus rotating speed are very different.In fact, due to the existence of the Coriolis effect the coupling between flapping and sweeping bending is induced, the first decreases until it becomes zero, while the second continues to increase.It is observed that instead of a rotating speed, there is a whole domain of rotating speed in which the first flapping frequency does not exist.In this domain, the flapping frequency becomes imaginary value, implying that the shaft becomes unstable.When the ply angle is decreased, in addition to shift of instability domain towards larger rotating speeds, it is also observed that the domain of instability is enlarged.
Figure 3 shows the variation of the first two flexural dampings versus rotating speed for various ply angles.It can be seen clearly that as the rotating speed is increased, the damping of flapping bending mode decreases and remains negative for all rotating speed, so the flapping bending mode is stable.From the results of Figure 3, it can be also observed that the dampings corresponding to sweeping bending mode are negative at low rotating speed and increase with increasing rotating speed, and at certain value of rotating speed the dampings vanish and then become positive.Transformation of damping from a negative to a positive value marks the onset of unstable motion.The rotating speed corresponding to zero damping is the threshold of instability of the shaft.The enclosed curves located nearby the threshold of instability represent that the real parts are conjugate.Figure 3  that the enclosed curves are shifted toward larger rotating speed but the extent of the enclosed curve is increased with the decrease of the ply angle.Figure 4 shows the variation of the first two extensiontwist natural frequencies versus rotating speed for various ply angles.As seen in Figure 4, due to the absence of the Coriolis effect, the first extension-twist natural frequency (the twist is dominant) decreases, while the second (the extension is dominant) remains constant at all rotating speeds.From Figure 4 it is seen that the effect of ply angle on the first extension-twist natural frequency is significant and is quite different from the case of flexural mode.The maximum critical speed can be reached when the ply angle  = 45 ∘ .Figure 5 shows the variation of the first two extensiontwist dampings versus rotating speed for various ply angles.It can be observed that the second extension-twist mode is stable whereas the first extension-twist mode is unstable as the rotating speed increases above certain value.The selfexcited range is easily identified from the figure by the sign of damping.It may also be noted that the effect of ply angle on the threshold of instability is similar to that previously described for the flexural mode.Figure 6 shows the variation of the first two flexural natural frequencies versus rotating speed for various aspect ratios.The results show that the critical rotating speed increases with the decrease of aspect ratios.
Figure 7 shows the variation of the first two flexural dampings versus rotating speed for various aspect ratios.From the results it can be seen that the threshold of instability increases as aspect ratio decreases.
Figures 8 and 9 present the effect of aspect ratio on the natural frequency-rotating speed curves and dampingrotating speed curves for the extension-twist mode, respectively.The results show that the decrease of aspect ratio yields

Figure 1 :
Figure 1: Geometry and coordinate systems of composite thinwalled shaft.

Figure 2 :Figure 3 :
Figure 2: The variation of natural frequencies with the rotating speed for various ply angles (/ = 3.6, first two flexural modes).

Figure 4 :Figure 5 :
Figure 4: The variation of natural frequencies with the rotating speed for various ply angles (/ = 3.6, first two extension-twist modes).

Figure 6 :
Figure 6: The variation of natural frequencies with the rotating speed for various aspect ratios ( = 45 ∘ first two flexural modes).

Figure 7 :
Figure 7: The variation of dampings with the rotating speed for various aspect ratios ( = 45 ∘ first two flexural modes).