Dynamic Analysis of a Fault Planetary Gear System under Nonlinear Parameter Excitation

In order to facilitate lubrication and avoid the gear stuck due to thermal expansion, there needs to be a gap between the tooth profiles. As a strong nonlinear factor, the backlash will affect the motion state of the planetary gear system. When the gear failures occur, the motion state of the system will accordingly change. In this study, the meshing stiffness of the gear pair with tooth tip chipping fault is calculated by combining the analytic geometry method and the potential energy method. Then, a new nonlinear dynamic model including tooth backlash, time-varying mesh stiffness, and manufacturing error is established to study the dynamic response of the system. The equations of motion are derived by the Lagrangian method and solved by the numerical integration method. Taking the excitation frequency and tooth backlash as the variation parameters, respectively, the dynamic characteristics of the system are analyzed by comparing the global bifurcation diagrams between the health system and the fault system, and the path of the system into chaos is revealed. At the same time, the local characteristics of the system are revealed through the phase diagrams and Poincaré maps. The results show that with the variation of excitation frequency and tooth backlash, the fault system presents a more complex motion state. This study can provide the theoretical support for dynamic design and fault diagnosis of planetary gear transmission systems under the environment of gear fault-prone.


Introduction
Planetary gear systems are widely used in aerospace, agricultural machinery, construction machinery, and other fields because of their high transmission ratio and transmission efficiency [1]. For the purpose that the gears can be fully lubricated and avoid jamming, tooth backlash between the engaged teeth is indispensable. However, the backlash is a strong nonlinear factor for the dynamic characteristics of the gear system. Due to the nonlinear factors, such as backlash, gear fault, and so on, the system exhibits different motion states with the variational excitation frequency. Hence, for better dynamic design and fault diagnosis for planetary gear systems, it is necessary to analyze their dynamic characteristics.
In recent decades, many scholars have conducted modeling analyses on the gear transmission system. In 1994, Kahraman [2] proposed a dynamic model of the planetary gear system and studied the load-sharing characteristics. Later, the vibration modes of the planetary gear system [3][4][5] and the suppression law of the meshing phase to the planetary modal response are analyzed [6]. en, the load distribution coefficient is deduced based on the model [7][8][9]. Wang et al. [10] revealed the nonlinear phenomena and evolution mechanisms of bifurcation and chaos via a threedegree-of-freedom torsional vibration model. en, Shen et al. [11] studied the dynamic characteristics of spur gear systems using the incremental harmonic balance method. In the planetary gear system, gear errors are inevitable, such as installation errors [12], manufacturing errors [13], and geometric errors [14], and the influence on the system is also analyzed [15].
Backlash and meshing stiffness are nonlinear factors for the gear system. To analyze the dynamic characteristics of the planetary gear system, the nonlinear dynamic model containing the tooth backlash [16][17][18] and time-varying meshing stiffness was established [19][20][21][22]. en, Huang et al. [23] and Pan and Vicuña [24] introduced the fractal backlash to the model. However, under high-speed and light-load conditions, gears may appear as tooth backside contact.
erefore, a multistate dynamic model was established by Liu et al. [25], and the variation law of the meshing force under different meshing conditions was analyzed. e influence of bearing clearance [26] and gear surface modification [27] on the dynamic response of the planetary gear system was revealed. Later, a nonlinear dynamic model of the multistage gear transmission was established by Zhao and Ji [28] and Xiang et al. [29], and the nonlinear dynamic characteristics of the system are analyzed at the same time.
Due to poor lubrication, impact load, and stress concentration, gear failures often occur and the dynamic response of the system will change accordingly. For early fault detection, Shi et al. [30] analyzed the fault characteristics under variable load. Pan et al. [31] analyzed the frequency components of the vibration signal in the fault and the healthy state. Considering the flexible ring gear and bearing fault, Liu et al. [32] presented a rigid-flexible coupling planetary gear dynamic model. By establishing the meshing stiffness model under spalling fault, Luo et al. [33] and Xiang et al. [34] analyzed the dynamic characteristics of the planetary gear system. Later, Shen et al. [35] proposed a purely torsional model to analyze the dynamic characteristics of planetary gear under wear fault. Yang et al. [36] proposed a nonlinear dynamic model with tooth backlash and bearing clearance and analyzed the vibration response under crack fault [37][38][39]. Luo et al. [40] established the meshing stiffness model for spalling and pitting failures and compared the dynamic response under different failure types.
In the previous literature, the dynamic response of the planetary gear system under crack failure, pitting failure, and spalling was researched. However, there are limited research studies on the nonlinear dynamic characteristics of planetary gear systems with tooth tip chipping fault under nonlinear parameters excitation. erefore, in order to reveal the dynamic characteristics of the fault system, this study established a nonlinear dynamic model containing tooth backlash, time-varying meshing stiffness, and static transmission error and analyzed the chaos and bifurcation characteristics of the system via choosing the excitation frequency and tooth backlash as variable parameters. e rest of the study is organized as follows. Section 2 establishes the planetary gear dynamic model. In Section 3, the dynamic characteristics of the system are analyzed by taking the tooth backlash and rotation speed as control variables. Finally, some conclusions are given in Section 4.

Dynamic Model of the Planetary Gear System
e planetary gear system has a complex structure that is different from fixed shaft gearboxes. In order to study the nonlinear dynamic behavior, a pure torsional dynamic model of the planetary gear system is proposed, which contains a sun gear s, a ring gear r, a carrier c, and N planet gears p i (i � 1, 2, . . ., N), as shown in Figure 1. In this model, all gears are standard spur gears. e carrier with planet gears fixed on is connected with the input shaft, and the sun gear is attached with the output shaft. e ring gear is fixed. In order to simplify the model, all components are assumed to be rigid. e potential energy method is a common method for calculating meshing stiffness. In the method, the gear tooth is regarded as a variable cross-section cantilever beam fixed on the tooth root, as shown in Figure 2. e total meshing stiffness can be composed of bending stiffness k b , axial compression stiffness k a , shear stiffness, Hertz contact stiffness k h , and fillet foundation stiffness k f . When in the gear occurs tooth tip chipping failure, a part of the material will fall off the tooth tip, which will affect the cross-sectional area and the moment of inertia of the tooth faulty part, as shown in Figure 3. ereby, to estimate the meshing stiffness of gear pairs with tooth tip chipping, an accurate stiffness model is established.
In order to simplify the model, the fracture surface is simplified as a plane, as shown in Figure 3. e intersection line between the plane and the tooth profile surface is the curve L 3 , and the intersection line with the gear end surface denotes the straight L 2 . In the three-dimensional coordinate system, the involute tooth profile equation Ω 1 can be expressed as where u and v represent the coordinate of the involute in the UOV coordinate system, which can be described as In the three-dimensional coordinate system, the fracture surface equation Ω 2 can be written as e equation of the curve L 2 can be deduced by simultaneous equations of (1) and (3). e equation of line L 3 in Figure 3(a) can be expressed as    Figure 3: Diagram of the tooth tip chipping part for the tooth.

Shock and Vibration
x � As tooth tip chipping failure occurs, a part of the material will fall off the tooth tip, so the cross-sectional area and moment of inertia of the faulty part will vary accordingly, which can be deduced by the following equations: ereby, the meshing stiffness can be obtained by the potential energy method [41]: where the parameters of L * , M * , P * , Q * , u f , and S f are given in the study by Ma et al. [42]. e total meshing stiffness of external-external mesh can be expressed as In this model, the internal-external mesh is assumed to be in a healthy state, and the total meshing stiffness can be deduced in the same way: 4 Shock and Vibration where n represents the number of tooth pairs engaged at the same time. e subscript 1 indicates the driving gear, and the subscript 2 represents the driven gear.

Manufacturing Error and
Damping. In this model, the carrier is connected to the input shaft, and the meshing frequency can be expressed as Assuming that the gear teeth are exactly the same, the static transmission error between the engaged gear teeth can be described as Fourier series, where the meshing frequency is the fundamental frequency. To simplify the model, merely, the fundamental frequency is taken into account, which can be expressed as [29] e spn � e aspn sin ω m t + ϕ spn , (10a) e mesh damping can be expressed as where k mspn and k mrpn are the arithmetic mean of k spn and k rpn , I j (j � s, ce, r, pn) is the moment of inertia, and m j (j � s, c, r, pn) is the equivalent mass that can be calculated by the following equations:

Gear Backlash.
In order to facilitate gear lubrication and prevent the gear from jamming during the operation, there is a gap between the engaged gear teeth. e backlash can be described by a piecewise function:

Dynamic Differential Equations.
According to the Lagrangian equation, the differential equations of torsional vibration under the excitations can be obtained: where T out and T in represent the output torque and the input torque, respectively. F spn and F rpn are adopted to describe the meshing force which can be expressed as [20] F spn � c spn _ δ spn + k spn f δ spn , (15a) In order to reduce the dimension of the differential equations, generalized coordinates are introduced as follows: Due to the large numerical difference between the various nonlinear parameters, the solving progress is difficult to conduct. In order to facilitate the solution and improve the calculation efficiency, the dimensionless processing is carried out. By introducing the time scale ω n and displacement scale b c , the dimensionless parameters can be obtained. e dimensionless time displacement, velocity, and acceleration can be expressed as Shock and Vibration 5 τ � ω n * t, e dimensionless mesh stiffness, static transmission error, and backlash are expressed as By substituting equations (15a)-(19b) into equation (14), the differential equation of dimensionless parameters can be obtained as follows:

Numerical Simulation and Results Analysis
Since the relative displacement between the engaged teeth has the same law, the dimensionless relative displacement of the engaged teeth between the sun gear and the planet gear is taken as the object to analyze, and the fourth-order Runge-Kutta method is used to solve the dynamics differential equations. e basic parameters of each component in the model are given in Table 1, and the meshing parameters are given in Table 2. e displacement scale b c is 1e − 5 m. e input torque is 100 Nm, and the transmission ratio is 4. e α a is assigned a value of 25°, y b is 0 mm, and z c is 20 mm. e meshing stiffness of external-external mesh and internalexternal mesh is shown in Figures 4 and 5. To avoid the influence of transient response, the first 500 response cycles of the system are omitted.

Bifurcation and Chaos of Fault and Health Systems with Excitation
Frequency. e motion state of the system will convert with the variation of the excitation frequency. erefore, the dimensionless excitation frequency Ω is chosen as the variable parameter, and the dimensionless backlash b is assigned the value of 2; the bifurcation diagram of the system in both health and fault conditions is shown in Figures 6 and 7. e largest Lyapunov exponent diagram of the health system is shown in Figure 8. As can be seen from figures, both the health and fault systems have rich chaotic and bifurcation characteristics. For the health system, when the excitation frequency Ω is between 0 and 0.69, the system presents a single-period motion and the corresponding largest Lyapunov exponent is less than zero. When it is between 0.69 and 1.28, the system undergoes a single-period motion to evolve into chaotic motion with small amplitude. When the excitation frequency changes in the range of 1.28-1.48, the system basically keeps cyclical movement, in which the largest Lyapunov exponent is less than zero. en, the system enters the chaotic motion eventually. For the fault system, it is similar to the state of motion experienced by the health system, but in the excitation frequency range corresponding to the periodic motion of the health system, the fault system shows obvious multiple periodic motions, as shown in Figure 9. Under different tooth backlashes, the bifurcation diagram of the system is displayed in Figure 10. It can be seen that the amplitude of chaotic motion augments with the increase of the initial backlash. e system enters the chaotic motion earlier as the speed increases if the backlash is designed large.
In order to further analyze the local characteristics of the health system and fault system, phase diagrams and Poincaré maps are used to describe the local characteristics of the system. Figures 11-15 are the phase diagram and Poincaré maps of the two systems at different excitation frequencies. It can be seen from Figure 11 that the health system exhibits a single-period motion at Ω � 0.5, while the fault system exhibits a multiple period motion. When Ω equals to 1, as shown in Figure 12, both of the two systems are in chaotic motion. When Ω is 1.34, as shown in Figure 13, the health system exhibits single-period motion, while the fault system presents multiple periodic motions. When Ω is 4, as shown in Figure 14, the health system exhibits quasiperiodic motion, while the faulty system appears as small amplitude chaotic motion. As Ω continues to increase, the health system exhibits the single-periodic motion at Ω � 1.47, while the faulty system occurs in chaotic motion, as shown in Figure 15. Hence, when the health system exhibits periodic motion, the fault system exhibits the multiperiod or pseudoperiodic motion at the corresponding excitation frequency.

Bifurcation and Chaos of Fault and Health Systems with
Different Initial Backlashes. In order to facilitate lubrication, there needs to be a gap between the engaged gear teeth and the value of the initial backlash that will affect the dynamic characteristics of the gear system. erefore, the dimensionless backlash b is chosen as the variable parameter to study the dynamic characteristics of the system with the dimensionless excitation frequency Ω to be 1.6. When the backlash varies between 0 and 3, the bifurcation diagrams of                 system motion state occurs jump variation at b � 1.095. en, the system enters a small amplitude chaotic motion at the interval of 1.095-1.29 corresponding to which the largest Lyapunov exponent is larger than zero. Afterwards, the system reenters the single-period motion when b varies in 1.29-1.53. As b continues to increase, and the system eventually enters chaotic motion. For the fault system, it is similar to the state of motion experienced by the health system, but in the backlash interval corresponding to the periodic motion state of the health system, the fault system      shows obvious small amplitude chaotic motion. Figure 19 shows the bifurcation characteristics of the two health systems with backlash. It can be seen from the figure that the fault system has more complex nonlinear characteristics and enters the chaotic motion earlier. After changing the speed, the bifurcation diagram of the fault system with the backlash is shown in Figure 20. When the speed increases, the amplitude of the chaotic motion will increase accordingly. In order to further analyze the local characteristics of the health system and fault system and reveal the path of the   system into chaos, phase diagrams and Poincaré maps are used to describe the local characteristics of the system. From Figures 21-24, it can be seen that the health system exhibits a single-period motion at b � 0.5, 1, and 1.4 and a small amplitude chaotic motion at 1.2. While, the fault system directly enters the chaotic motion.

Conclusions
is study proposed a meshing stiffness model of gear pair with tooth tip chipping fault by combining the analytic geometry method and potential energy method. en, a new nonlinear dynamic model of the planetary gear system under health and fault conditions is established considering the time-varying mesh stiffness, tooth backlash, and static transmission error. e model is derived by the Lagrangian method and solved by the numerical integration method. Taking the dynamic transmission error between the sun gear and planet gear as the research object, through the global bifurcation diagram, the variations of the two system motion states with the excitation frequency and tooth backlash are analyzed, and the local characteristics of the systems are analyzed via the phase diagrams and the Poincaré maps. Some conclusions can be obtained as follows: (1) With the variation of excitation frequency, both the health system and fault system show complex bifurcation and chaotic characteristics, and the chaotic motion is mixed with periodic windows; (2) e fault system has more complex nonlinear characteristics and enters the chaotic motion earlier. e fault system shows multiple periodic motions and a small amplitude chaotic motion in the interval