Dynamic Modeling and Stability Analysis for a Spur Gear System Considering Gear Backlash and Bearing Clearance

: In practice, gear backlash and bearing clearance usually exist together in a gear system. T hey may induce complicated dynamic responses and degrade transmission performance. Up to now, although each of them has been researched, little attention has been paid to the coupling dynamic characteristics of gear backlash and bearing clearance. In a limited number of relevant studies, since the linear collision models they adopted are difficult to realistically depict actual collision behaviors caused by bearing clearance, these studies cannot accurately reveal the coupling dynamic characteristics of gear backlash and bearing clearance. Furthermore, system stability of the gear sys-tem considering gear backlash and bearing clearance has not been thoroughly investigated. In view of this, this paper contributes to the research on dynamic modeling and stability analysis for the spur gear system considering gear backlash and bearing clearance. A nonlinear collision model with time-varying contact stiffness/damping is suggested for describing the bearing collision behaviors. Based on the geometrical relationship of dynamic center distance, dynamic working pressure angle, and dynamic backlash, the coupling motion model of gear backlash and bearing clearance is developed. On this basis, the dynamic model of the spur gear system considering gear backlash and bearing clearance is established and verified by numerical simulations, virtual prototyping simulations and experiments. Afterwards, to thoroughly explore the complicated dynamic characteristics of the gear system considering gear backlash and bearing clearance, several important parameters, i.e., rotational speed, gear backlash and bearing clearance, are chosen as bifurcation parameters to study their influences on system stability via bifurcation diagrams, time-domain waveforms, FFT spectra, Poincaré maps, and phase diagrams. Various complicated nonlinear behaviors, such as hopping, multiple periodic motion, quasi-periodic motion, and chaotic motion, are revealed. This study can provide useful reference for the multi-clearance coupling research of the gear system in complicated working environments.


Introduction
Gear systems are widely used in many industrial applications, such as wind turbines, ships, robots, and so on.However, various internal and external excitations may induce complicated dynamic responses and degrade transmission performance.To effectively deal with these problems, it is necessary to establish a relatively accurate dynamic model to gain an understanding of the important dynamic characteristics of gear systems.
In past decades, various dynamic models of gear systems have been put forward for different purposes.Several factors, such as time-varying meshing stiffness [1][2][3], transmission error [4][5][6], and tooth worn [7,8] have been taken into account in these dynamic models.Accordingly, the influence of these factors on nonlinear dynamic characteristics of gear systems have been analyzed based on these dynamic models.Kank Recently, the clearance problems of gear systems have attracted extensive attention.To implement good lubrication and prevent gear teeth from sticking, a certain amount of backlash needs to be reserved in the gear pair.This gear backlash is referred to as initial backlash or static backlash.Inevitably, it may induce the collision and vibration of the gear system.Kahraman and Singh [9] firstly established the dynamic equations considering gear backlash and used piecewise functions to describe the contact and separation states of gear teeth.Kahraman [10] then developed an experimental setup with gear backlash and found that gear backlash can cause the separation and impact of the gear pair in the case of light loads.Based on Kahraman's work, Sun [11] used the harmonic balance method to calculate the dynamic response of a planetary gear system considering gear backlash and analyzed the influence of meshing stiffness on nonlinear dynamic characteristics of the gear system.Zhao [12] developed a pure torsional model of a planetary gear system with gear backlash and analyzed the influence of internal and external excitations on transmission performance.Yang [13] and Li [14], respectively, analyzed the influence of gear backlash on transmission error and system stability for different types of planetary gear systems.Wang [15] established a gear meshing model using the finite element method and analyzed the influence of gear backlash on dynamic transmission error.Qi [16] suggested a gear backlash model based on fractal theory and studied the influence of gear backlash on the dynamic performance of the gear system.Nevertheless, most of them only pay attention to static backlash.In fact, gear geometric eccentricity or shaft elastic deformation can result in time-varying gear backlash and time-varying pressure angle.Chen [17] first studied dynamic gear backlash and dynamic working pressure angle in 2011.Based on Chen's work, Yi [18] suggested a nonlinear dynamic model of the spur gear system considering dynamic gear backlash and dynamic working pressure angle.Liu [19] analyzed the influence of gear center distance variation on gear meshing states and studied the influence of speed and torque on the dynamic characteristics of the gear system.
In addition to gear backlash, bearing clearance can also seriously degrade the transmission performance of gear systems.Kahraman [20] first introduced bearing clearance into the dynamic model of the gear system and analyzed the influence of bearing stiffness on dynamic characteristics of the gear system.Huang [21] established a three-degree-offreedom dynamic model with bearing clearance and studied the influence of bearing clearance on the system stability of the gear system.Zhang [8] established a six-degreeof-freedom bending-torsional gear-rotor-bearing model with bearing clearance and analyzed the influence of different bearing clearance on the radial vibration of the gear system.Suzuki [22] developed a planetary gear test bench and measured the load distribution of the planetary gear system in the case of different bearing clearance.Gu [7] developed an experimental setup of the gear system and studied the influence of different bearing clearance on the transmission performance of the gear system by changing gears and bearings.However, these studies used the linear collision model to depict the bearing collision behaviors caused by bearing clearance.Although this widely accepted linear model is convenient for calculation, it is difficult to accurately describe actual collision process characterized by strong nonlinearity.
In practice, gear backlash and bearing clearance usually exist together in a gear system.Since bearing clearance may induce the radial vibration of gears, gear center distance and pressure angle are time-varying, thereby causing the change of gear backlash.This dynamic gear backlash, in turn, further influence the radial vibration of gears by altering the meshing force.Therefore, there exists the complicated coupling of gear backlash and bearing clearance.Up to now, although gear backlash and bearing clearance have been studied separately, few studies have been conducted on their coupling dynamic characteristics.Gu [7] established a six-degree-of-freedom gear pair model considering dynamic backlash and analyzed the influence of gear backlash and bearing clearance on the steady responses of the gear system.Zhang [8] established a gear-rotor-bearing model considering multiple clearances and analyzed the influence of different clearance on vibration amplitude.However, these studies were limited in number and still adopted the linear collision model.Therefore, their dynamic models cannot reveal the coupling effect of dynamic working pressure angle and nonlinear bearing collision.Furthermore, the system stability problem of the gear system considering gear backlash and bearing clearance has not been thoroughly investigated.
In view of this, this paper contributes to the research on dynamic modeling and stability analysis for the spur gear system considering gear backlash and bearing clearance.A nonlinear collision model with time-varying contact stiffness/damping is suggested for describing the bearing collision behaviors.Based on the geometrical relationship of dynamic center distance, dynamic working pressure angle, and dynamic backlash, the coupling motion model of gear backlash and bearing clearance is developed.On this basis, the dynamic model of the spur gear system considering gear backlash and bearing clearance is established and verified by numerical simulations, virtual prototyping simulations, and experiments.Afterwards, to thoroughly explore the complicated dynamic characteristics of the gear system considering gear backlash and bearing clearance, several important parameters, i.e., rotational speed, gear backlash, and bearing clearance, are chosen as bifurcation parameters to study their influences on system stability.Various complicated nonlinear behaviors, such as hopping, multiple periodic motion, quasi-periodic motion, and chaotic motion are revealed.This study can provide useful reference for the further multi-clearance coupling research of the gear system in complicated working environments.

Meshing Model of a Spur Gear Pair with Gear Backlash and Bearing Clearance
A lumped parameter model for a pair of parallel shaft involute spur gears with gear backlash and bearing clearance is established, as shown in Figure 1.The pinion and the gear are supported by two rolling bearings.The blue dotted circles indicate their positions without radial bearing clearance, and the black circles indicate their positions with radial bearing clearance.The meshing between the pinion and the gear is characterized by a nonlinear spring and a viscous damper.Dynamic analysis is conducted in the plane of the gear pair.Here, m1, m2, I1, I2, rb1 and rb2 denote the mass, the moment of inertia, and the base circle radius of the pinion and the gear, respectively; 1 T and 2 T denote the input and the output torques, respectively; r1 f and r2 f denote the bear- ing stiffness, the bearing damping and the bearing clearance of the pinion and the gear, respectively; ( ) t K t , t D and t f denote the meshing stiffness, the meshing damping and the gear backlash, respectively; r1 F and r2 F denote the radial collision forces of the pin- ion and the gear, respectively; t1 F and t2 F denote the dynamic meshing forces of the pinion and the gear, respectively; e(t) denotes the transmission error of the gear system.Figure 1.System description of an involute spur gear pair with gear backlash and bearing clearance.

Nonlinear Bearing Collision Model Considering Bearing Clearance
Due to the existence of bearing clearance, various complicated dynamic phenomena inevitably arise.Up to now, most relevant studies simplified the bearing as two mutually perpendicular linear springs with viscous damping.
Although this widely accepted model is convenient for calculation, it is difficult to accurately depict actual situations.Firstly, the bearing collision force should not be simplified as a linear function of the embedding amount.In fact, it is a nonlinear function of the embedding amount.Secondly, since the embedding amount is decomposed in two mutually perpendicular directions, the bearing collision force in one direction relies only on the embedding amount in this direction.Therefore, the coupling effects between two mutually perpendicular directions are ignored.Thirdly, the bearing clearance area is viewed as a rectangle.This assumption is unreasonable and not applicable to actual situations.Finally, the bearing stiffness and the bearing damping are viewed as constant.However, this assumption is only applicable to large clearance and small load.Strictly speaking, both of these are influenced by the embedding amount, and thus are time-varying.
In recent years, the nonlinear contact theory has been used to describe the collision process of the mechanical system with joint clearance.However, it is rarely used in dynamic modelling of the gear system with gear backlash and bearing clearance.In view of the above difficulties that the linear collision model faces, a new bearing collision model is suggested based on the nonlinear contact theory.As shown in Figure 2, bearing collision is modeled as contact between a shaft and a sleeve.Due to bearing clearance, there exist two states between the shaft and the sleeve, i.e., contact and separation.To accurately describe these states, the embedding amount function δ is modeled as , the shaft and the sleeve are in contact state, and the collision force is modeled as where r K and r D are the contact stiffness and the contact damping, respectively; n is the coefficient determined by material and is usually set to 1.5 for the metallic material; v is the relative velocity between the shaft and the sleeve, and is calculated by From Equation ( 2), it can be seen that the suggested collision model is nonlinear.Since the collision force is not a linear function of the embedding amount, the collision force in the X direction no longer relies only on the embedding amount in the X direction.The same is true for the collision force in the Y direction.Therefore, the coupling effects that have been ignored by the linear collision model can be taken into account in this nonlinear collision model.Besides, the assumption of a rectangular bearing clearance area is discarded in this nonlinear collision model.
In the suggested nonlinear collision model, the contact stiffness r K and the contact damping r D are two important parameters used to depict bearing collision behaviors.
The former indicates the relationship between the collision force and the embedding amount.The latter indicates energy dissipation during contact process.Based on the Hertz contact theory and the collision experiments, they can be obtained by where L is the axial length of the bearing; 1 υ and 2 υ are the Poisson's ratio of the shaft and the sleeve, respectively; 1 E and 2 E are the elastic modulus of the shaft and the sleeve, respectively. is the initial collision velocity.As shown in Equation ( 4), the contact stiffness r K is not a constant value.It relies on the structural dimension, material, bearing clearance, and embedding amount.Even though the structural dimension, material, and bearing clearance remain unchanged, the contact stiffness still varies with the time-varying embedding amount.Therefore, the assumption of constant contact stiffness in the linear collision model is not applicable to actual situations.
Similarly, the contact damping r D is not a constant value either, as shown in Equa- tion (5).It depends on the recovery coefficient, embedding amount, initial collision velocity, and contact stiffness.Obviously, it varies with the time-varying embedding amount and the time-varying contact stiffness during the contact process.Therefore, the assumption of constant contact damping in the linear collision model is not applicable to actual situations.
In conclusion, compared to the linear collision model, the suggested nonlinear collision model can better describe actual bearing collision behaviors.

Coupling Motion Model of Gear Backlash and Bearing Clearance
In practice, gear backlash and bearing clearance exist together in a gear system.Although each of them has been researched, their coupling dynamic characteristics are seldom investigated.In view of this, the coupling motion model of gear backlash and bearing clearance is developed in this section.
As shown in Figure 3, O1 and O2 are the centers of the pinion and the gear, respectively; DH is the common tangent of the base circles of the pinion and the gear, and is referred to as the line of action (LOA); b1 r and b2 r are the radius of the base circle of the pinion and the gear, respectively; Since bearing clearance may induce the radial vibration of gears, both the actual center distance ' A and the pressure angle '  α are time-varying, and can be obtained by ( ) ( ) where A0 is the initial center distance; xi and yi (i = 1, 2) are the radial vibration displacement in the X and Y directions of the pinion and the gear, respectively.In the ideal meshing process, the point B of the pinion should coincide with the point C of the gear.However, due to the existence of both the initial backlash and dynamic center distance, the dynamic backlash 2bt between the point B and the point C becomes where ' t is the tooth pitch on the pitch circle of the pinion, ( ) ( ) α α ; t is the tooth pitch on the standard pitch circle of the gear; where 1 s and 2 s are the tooth thickness on the standard pitch circle of the pinion and the gear, respectively; 1 r and 2 r are the radius of the standard pitch circle of the pinion and the gear, respectively.Based on Equations ( 8)-( 11), the dynamic backlash t 2b is obtained by It can be seen that the dynamic backlash is a nonlinear function of the actual center distance and the actual pressure angle.Since radial bearing vibration may change the center distance and the pressure angle, the actual gear backlash varies with the time-varying center distance and the time-varying pressure angle.

Meshing Force Model Considering Dynamic Backlash
The dynamic backlash induced by bearing clearance may lead to a complicated change in the meshing force.As a result, gear meshing states are divided into three types, i.e., positive meshing, separation, and tooth back meshing.As shown in Figure 4, these meshing states can be obtained from the dynamic transmission error (DTE) along the line of action DH.
where b1 1 b2 2 -r r θ θ is angle margin; Based on the dynamic backlash and the dynamic transmission error, the embedding amount is obtained by According to the embedding amount, the corresponding meshing state can be determined.When t g is greater than or equal to t b , the gear pair is in positive meshing state.When t g is less than t b and greater than t b − , the gear pair is in separation state and the meshing force is zero.When t g is less than or equal to t b − , the gear pair is in tooth back meshing state.
Based on the dynamic backlash and the meshing state, the time-varying meshing force t F is expressed by where t K is the meshing stiffness and can be obtained via the finite element analysis or physical experiments, t D is the meshing damping and can be expressed by where ξ is the meshing damping ratio of the gear and is usually between 0.03-0.17.

Derivation of Equations of Motion
Based on the above meshing model of a spur gear pair with gear backlash and bearing clearance, the differential equations of motion for the gear system are derived using Lagrange's equation.They have six generalized coordinates, i.e., [ ] , , , , , q x y x y θ θ = (17) The equations of motion of the gear system are obtained as follows Compared with previous models, this dynamic model adopts the nonlinear collision model, and thus can better reveal nonlinear dynamic characteristics of the gear system.

Comparison and Validation of Dynamic Models
The previous dynamic model is established based on the linear collision model and the constant contact stiffness/damping, while the new dynamic model is established based on the nonlinear collision model and the time-varying contact stiffness/damping.In this section, a pair of parallel shaft involute spur gears are used as an example to compare and verify these two models.The bearing clearance is 20 μm and the initial gear backlash is 50 μm.Other parameters are listed in Table 1.

Model Comparison
Although both the previous and new dynamic models consider gear backlash and bearing clearance, the latter can depict the dynamic behaviors of the gear system more realistically than the former does, because it takes into account the coupling effect of dynamic working pressure angle and nonlinear bearing collision.
For the convenience of understanding, three models are compared in this section with the help of numerical simulations based on the Runge-Kutta method.They are the dynamic nonlinear model (i.e., the new dynamic model characterized by the dynamic working pressure angle and the nonlinear bearing collision model), the dynamic linear model (i.e., the previous dynamic model characterized by the dynamic working pressure angle and the linear bearing collision model), and the static nonlinear model, characterized by the static working pressure angle and nonlinear bearing collision model.
As shown in Figure 5, the area enclosed by the blue line represents the bearing clearance area, the centre of which indicates the center of the sleeve; the red point represents the trajectory of the center of the driven gear.In order to depict the trajectory of the center more clearly, it is enlarged in the Figure .When the center of the driven gear falls into the bearing clearance area, the shaft and the sleeve are in separation state.Otherwise, the shaft and the sleeve are in contact state.
In the dynamic nonlinear model, since the working pressure angle is dynamic, the direction of the LOA is also time-varying, thereby exerting the excitations that do not coincide with the direction of the initial LOA to the driven gear.Under the action of these excitations with time-varying direction, the bearing collision may occur.Since the collision force is a nonlinear function of the embedding amount in the nonlinear bearing collision model, the resultant collision response may further alter the direction of the meshing force.Therefore, there exists the complicated coupling of the dynamic working pressure angle and the nonlinear bearing collision.In this case, the motion of the center of the driven gear includes not only the motion along the LOA, but also the circular motion perpendicular to the initial LOA, as shown in Figure 5a.This phenomenon exhibits the coupling effect of dynamic working pressure angle and nonlinear bearing collision on the dynamic behaviors of the gear system.
In the static nonlinear model, since the working pressure angle is constant, the direction of the LOA is almost constant.Therefore, the direction of the meshing force is almost constant too.Although the collision force is a nonlinear function of the embedding amount, the resultant collision response cannot significantly alter the direction of the meshing force due to constant pressure angle.As a result, the trajectory of the center of the driven gear is almost a straight line, as shown in Figure 5b.
In the dynamic linear model, since the working pressure angle is dynamic, the direction of the LOA is also time-varying, thereby exerting the excitations that do not coincide with the direction of the initial LOA to the driven gear.However, since the bearing colli-sion in the X and Y directions are viewed as independent of each other in the linear bearing collision model, the coupling of the dynamic working pressure angle and the bearing collision is weakened due to the linear collision assumption.Therefore, the trajectory of the center of the driven gear is almost a straight line, as shown in Figure 5c.
In conclusion, the dynamic working pressure angle exerts the time-varying excitations that do not coincide with the direction of the initial LOA to the driven gear.The nonlinear bearing collision further amplifies this effect, thereby generating a circular motion perpendicular to the initial LOA.This phenomenon can be found in the following virtual prototyping simulations and physical experiments.Therefore, the new dynamic model can depict the dynamic behaviors of the gear system more realistically than the previous dynamic model does.At first, the vibration responses of these models are analyzed.Since the vibration behaviors in the X and Y directions are similar and the latter is larger than the former, the vibration response in the Y direction is used as an example to compare them.As shown in Figure 6a, the vibration amplitude of the dynamic nonlinear model is larger than those of other two models, due to the coupling effect of dynamic working pressure angle and nonlinear bearing collision.In addition, the linear bearing collision model separately calculates the bearing clearances in the X and Y directions.Therefore, the vibration displacement of the static nonlinear model is larger than those of other two models.In the frequency-domain, as shown in Figure 6b, the vibration response of both the static nonlinear model and the dynamic linear model mainly converges at the meshing frequency e f ( e f = z1 r f = 133 Hz).However, in addition to the meshing frequency, the vibration re- sponse of the dynamic nonlinear model is distributed at a circular vibration frequency a f , i.e., 345 Hz.This circular vibration response results from the coupling effect of dynamic working pressure angle and nonlinear bearing collision on the dynamic behaviors of the gear system, as verified by the above analysis.Then, DTE of these models are analyzed, which is one of the important indexes that evaluate the transmission performance of the gear system.As shown in Figure 7a, the amplitude of DTE of the dynamic nonlinear model is larger than those of other two models.Similarly, in addition to integer multiple of the meshing frequency, the circular vibration at 345 Hz plays an important role in the DTE response of the dynamic nonlinear model, as shown in Figure 7b.Again, this phenomenon cannot be revealed by other two models.On the other hand, it can be seen in Figure 8a,b (red coil represent phase trajectories, and blue dots represent poincare mapping points) that the DET response of the dynamic nonlinear model is the quasi-periodic motion while the DET response of the dynamic linear model is the periodic-1 motion.Obviously, the former is more complicated than the latter.
Through above analysis and comparisons, it can be concluded that the dynamic nonlinear model can reveal various important nonlinear characteristics of the gear system.Therefore, it is more realistic than other two dynamic models.

Model Validation
Although the above comparison has proven that the new dynamic model is more effective than the previous model, it is implemented based on the equations of motion derived by ourselves.In view of this, a series of virtual prototyping simulations and physical experiments are conducted to verify the correctness of the new dynamic model.
Firstly, the well-known ADAMS ® software is used to conduct virtual prototyping simulations.In this example, both the initial gear backlash and the bearing clearance are 50 μm; the rotational speed of the pinion is 240 rpm; other parameters are listed in Table 1.
Secondly, an experimental setup is developed to further verify the correctness of the proposed dynamic model, as shown in Figure 9.The pinion is connected to the motor shaft and is driven by the servomotor.The collision shaft is fixed on the frame by bolts.The gear is supported by a thrust bearing mounted on the collision shaft.The diameter of the center hole of the gear is designed to be larger than the diameter of the collision shaft, and the resultant clearance is used to simulate the bearing clearance.The accelerometer is attached to the end of the collision shaft and used to measure the real-time collision acceleration signal of the gear along the LOA, which is affected by the DTE of the gear system.
The acceleration signal is collected by the data acquisition card and processed by the upper computer.The parameters of the experimental setup are the same as those of the virtual prototyping model.As shown in Figures 10a and 11a, the average DTE responses of the proposed dynamic model and the virtual prototyping model are 62.62 μm and 61.57μm, respectively.They are very close and their error is only 1.677%.However, the average DTE response of the previous dynamic model is 84.71 μm, which is larger than those of both the proposed dynamic model and the virtual prototyping model.Besides, the DTE responses of both the proposed dynamic model and the virtual prototyping model in the frequency-domain mainly converge at integer multiple of the meshing frequency (i.e., 80.01 Hz and 160.00 Hz in the proposed dynamic model; 79.92 Hz and 160.07 Hz in the virtual prototyping model) as well as a special circular vibration frequency (i.e., 136.5 Hz in the proposed dynamic model; 127.91 Hz in the virtual prototyping model), as shown in Figures 10b and  11b.Their errors in the one-time meshing frequency, the two-times meshing frequency, and the radial vibration frequency are small, i.e., 0.112%, −0.044% and 6.716%, respectively.However, this special circular vibration frequency cannot be found in the previous dynamic model.It should be noted that, since the contact stiffness/damping cannot be time-varying in the virtual prototyping simulations, its simulation results are slightly different from those of the proposed dynamic model.Based on these comparisons between the proposed dynamic model and the virtual prototyping model, as listed in Table 2, it can be seen that the proposed dynamic model is more realistic than the previous dynamic model.
On the other hand, the DTE responses of both the proposed dynamic model and the experimental setup in the frequency-domain mainly converge at the integer multiple of the meshing frequency (i.e., 80.01 Hz and 160.00 Hz in the proposed dynamic model; 79.91 Hz and 160.03Hz in the virtual prototyping model) as well as a special circular vibration frequency (i.e., 136.5 Hz in the proposed dynamic model; 144.00 Hz in the virtual prototyping model), as shown in Figures 10b and 12b.Their errors in the one-time meshing frequency, the two-times meshing frequency and the radial vibration frequency are small,

Upper computer
Adapter Gear Accelerometer

Data acquisition card Motor driver Motor Pinion
i.e., 0.125%, −0.019% and −5.208%, respectively.Since the DTE cannot be directly measured in the experimental study and the dynamic response is measured in the form of acceleration signal using an accelerometer, there is a slight difference between the experimental results and the theoretical results.Furthermore, the influences of different bearing clearances on the vibration amplitude of the gear system are studied, as shown in Figure 13.The numerical simulations exhibit the same variation tendency as the physical experiments, i.e., the vibration amplitude of the gear system becomes large with the increase of the bearing clearance.Based on these comparisons of the numerical simulations and the experimental results, as listed in Table 2, it can be seen that the proposed dynamic model can provide more authentic results than the previous dynamic model.

Stability Analysis
Due to the existence of gear backlash, bearing clearance and time-varying contact stiffness/damping, the gear system is a complex dynamic system with strong nonlinearity.To thoroughly understand the complicated dynamic characteristics of the gear system, the rotational speed, the gear backlash and the bearing clearance are chosen as bifurcation parameters to study their effects on system responses.The equations of motion of the gear system are solved by the fourth-fifth order Runge-Kutta method.Its nonlinear dynamic characteristics are analyzed via bifurcation diagrams, time-domain waveforms, FFT spectra, Poincaré maps, and phase diagrams.

Effects of Rotational Speed of Pinion
The rotational speed n1 of the pinion is one of the key parameters that affect the dynamic behaviors of the gear system.In this section, the global characteristics of the new and previous dynamic models are compared using bifurcation diagrams.Then the nonlinear behaviors of the gear system are further investigated in detail.
Let the initial backlash b0 be 50 μm and the initial bearing clearance cr be 20 μm.Other parameters are listed in Table 1. Figure 14a shows the bifurcation diagram of DTE of the new dynamic model with respect to the rotational speed n1 of the pinion.When n1 ∈ [100, 970] rpm, the gear system is in quasi-periodic motion state.When n1 ∈ [980, 1030] rpm, the gear system changes from quasi-periodic motion to period-1 motion through hopping.In the rotational speed range of n1 ∈ [1040, 1640] rpm, the gear system returns to the quasiperiodic motion.When n1 ∈ [1650, 1820] rpm, the gear system falls into period-1 motion through hopping and then becomes period-nT motion.Finally, when n1 is larger than 1830 rpm, the gear system enters into chaotic motion.
Figure 14b shows the bifurcation diagram of DTE of the previous dynamic model with respect to the rotational speed n1 of the pinion.When n1 is less than 1020 rpm, the gear system is in period-1 motion state.When n1 ∈ [1040, 2120] rpm, the gear system undergoes a series of transitions between period-1 and period-nT motion through hopping.Finally, when n1 is larger than 2140 rpm, the gear system enters into chaotic motion.From Figure 14, it can be seen that the new dynamic model exhibits more complicated behaviors than the previous dynamic model.When the rotational speed is low, the previous dynamic model is in period-1 motion, while the new dynamic model undergoes quasiperiodic motion.With the increase of the rotational speed, there exist hopping phenomena in each bifurcation diagram, causing load fluctuation between teeth.However, more hopping phenomena can be found in the new dynamic model.Furthermore, although both dynamic models exhibit evolution paths to chaotic motion, the new dynamic model enters into chaotic motion earlier than the previous dynamic model.Figures 17 and 18 show the dynamic response of the gear system at the rotational speed n1 = 1700 rpm and 1780 rpm, respectively.As shown in Figures 17a and 18a, compared with Figure 15a, their amplitudes of DTE are larger than that in low rotational speed.From Figures 17b and 18b, it can be seen that their vibration responses are mainly distributed at the meshing frequencies e nf .Compared with Figure 17b, Figure 18b shows that its amplitude is higher at map, indicating that the gear system is in period-1 motion state.However, compared with Figure 17c, the thickness of phase diagram becomes larger and the Poincaré map shows several twisted-closed curves in Figure 18c, indicating that the gear system is in period-nT motion state.In summary, 1700 rpm and 1780 rpm are two representative rotational speeds in the range of [1650, 1820] rpm.As shown in Figure 14a, this range represents the evolution path of the gear system from periodic motion to chaotic motion.
Figure 19 shows the dynamic response of the gear system at the rotational speed n1 = 2500 rpm.As shown in Figure 19a, its amplitude is large and irregular.From Figure 19b and Figure 19c, it can be seen that its frequency spectrum is continuous and its phase diagram is confusion.These phenomena indicates that the gear system enters into chaotic motion state and its motion becomes unpredictable.

Coupling Effects of Initial Gear Backlash and Bearing Clearance
In this section, the coupling effects of the initial gear backlash and the bearing clearance are researched using bifurcation diagrams with respect to the initial gear backlash b0 under different bearing clearance cr = 0, 10, 20, 30 μm when the pinion's rotational speed n1 = 2500 r/min.Other parameters are listed in Table 1.
As shown in Figure 20a, if there is no bearing clearance, i.e., cr = 0, the gear system is in stable period-nT motion state.But such working condition is rare in practical applications.If the bearing clearance is small, as shown in Figure 20b,c, with the increase of the initial gear backlash, the gear system will evolve from stable period-nT motion to unstable chaotic motion.Therefore, large gear backlash may decrease the stability of the gear system.If the bearing clearance is large, i.e., cr = 30 μm, the gear system is always in unstable chaotic motion state, as shown in Figure 20d.Therefore, it is not recommended to adopt excessive bearing clearance in practical applications.
Next, two typical working conditions, i.e., cr = 10 and 20 μm, are analyzed.As shown in Figure 20b, the stable zone is large if the bearing clearance is 10 μm.The gear system is in stable motion state when the initial backlash falls in (0, 61) μm.When the bearing clearance is 20 μm, however, the gear system will enter into the chaotic motion state if the initial backlash is larger than 44 μm, as shown in Figure 20c.Therefore, the increase of bearing clearance will decrease the stability zone of the gear system.
On the other hand, due to the coupling effects of the gear backlash and the bearing clearance, the gear system exhibits complicated nonlinear dynamic characteristics.As shown in Figure 20b, if the bearing clearance is 10 μm, there exists chaotic motion in the stable zone when the initial gear backlash is 28 μm and 37 μm.Similarly, as shown in Figure 20c, if the bearing clearance is 20 μm, there exists a mixture of chaotic motion and period-nT motion in the stable zone when the initial backlash is in (0, 27) μm.The system motion is stable only in (28, 43) μm.
To sum up, when analyzing the stability problem of the gear system, the coupling effects of the gear backlash and the bearing clearance should be considered.Generally speaking, large gear backlash and large bearing clearance may decrease system stability.Therefore, small gear backlash and small bearing clearance are advisable.Besides, in some working conditions, small gear backlash will also make the gear system enter an unstable motion state.Therefore, it should be chosen according to specific working conditions.

Conclusions
In this paper, we research dynamic modeling and stability analysis for the spur gear system considering gear backlash and bearing clearance.A nonlinear bearing collision model with time-varying contact stiffness/damping is suggested for describing the bearing collision behaviors, which can take into account the coupling effects that have been ignored by the linear collision model.The proposed dynamic model of the gear system considers not only the coupling effect of gear backlash and bearing clearance, but also the effect of dynamic working pressure angle and nonlinear bearing collision.Therefore, it can depict the dynamic behaviors of the gear system more realistically than the previous dynamic model does.Several important parameters, i.e., rotational speed, gear backlash, and bearing clearance, are chosen as bifurcation parameters to investigate their influences on system stability via bifurcation diagrams, time-domain waveforms, FFT spectra, Poincaré maps, and phase diagrams.Various complicated nonlinear behaviors, such as hopping, multiple periodic motion, quasi-periodic motion, and chaotic motion are revealed.This study can provide useful reference for the multi-clearance coupling research of the gear system in complicated working environments. a

O
and b O are the centers of the sleeve and the shaft, respectively; a e and b e are the vibration displacement vectors of the sleeve and the shaft, respectively; ab b a = -e e e is the eccentricity vector between the sleeve and the shaft; a R and b R are the radius of the sleeve and the shaft, respectively;θ is the angle between the eccentricity vector and the X axis.

Figure 2 .
Figure 2. Geometrical relationship between a shaft and a sleeve.
c R R is the initial bearing clearance, δ denotes the embedding amount.If ab r | |<c e , the shaft and the sleeve are in separation state, and the collision force is zero.If ab r | | c ≥ e

2 s
are the tooth thickness on the pitch circle of the pinion and the gear, respectively; ' α is the working pressure angle on the pitch circle; 0 α is the pressure angle on the standard pitch circle; and b0 is the initial backlash.According to the meshing principle of the involute gear pair, ' 1 s and ' 2 s are derived as (

1 θ and 2 θα 1 ω
are the rotational angles of the pinion and the gear, respectively; ( ) is the radial vibration displacement along the LOA in the X direction; ( ) is the radial vibration displacement along the LOA in the Y direction.( ) ( ) m cos n e t e t ω = is static transmission error related to manufacturing and installation accuracy of gears, m e is the error fluctuation amplitude, is the rotational frequency of the pinion, and 1 z is the teeth number of the pinion.

Figure 5 .
Figure 5. Trajectory of center of driven gear.However, the dynamic responses of the static nonlinear model, the dynamic linear model, and the dynamic nonlinear model are compared via the time-domain waveforms, FFT spectra, and phase diagrams when the rotating speed n1 = 400 rpm and the load torque T = 10 Nm.At first, the vibration responses of these models are analyzed.Since the vibration behaviors in the X and Y directions are similar and the latter is larger than the former, the vibration response in the Y direction is used as an example to compare them.As shown in Figure6a, the vibration amplitude of the dynamic nonlinear model is larger than those

Figure 6 .
Figure 6.Vibration response in time-domain and frequency-domain.

Figure 9 .
Figure 9. Experimental setup of gear system.

Figure 14 .
Figure 14.Bifurcation diagram of DTE with respect to rotational speed n1.

Figure 16a exhibits theFigure 16 .
Figure 16a exhibits the dynamic behaviors of the gear system at n1 = 1000 rpm.The vibration responses are mainly distributed at the meshing frequencies e nf and the circu-

e 4 f , e 5 f , e 6 f and e 7 f
. In addition, the phase diagram of the gear system in Figure17cis a closed ring and there is only one point in the Poincaré e

Figure 20 .
Figure 20.Bifurcation diagrams of DTE response with respect to initial gear backlash under different bearing clearance.

Table 1 .
Main parameters of the gear pair.

Table 2 .
Comparison results of numerical simulation, virtual prototyping simulation and physical experiment.