Nonlinear Dynamic Characteristic Analysis of a Coated Gear Transmission System

: Coatings can signiﬁcantly improve the load-carrying performance of a gear surface, but how they a ﬀ ect the vibration characteristic of the system is an urgent issue to be solved. Taking into account the nonlinear factors like the variable mesh sti ﬀ ness, friction, backlash, and transmission error, a six-degree-of-freedom spur gear transmission system with coatings is presented. Meanwhile, the ﬁnite element method is applied to acquire the time-varying mesh sti ﬀ ness of the coated gear pair in the engagement process. With the support of the time-history curve, phase curve, Poincare map, and fast Fourier transform spectrum, the dynamic characteristics and the e ﬀ ects of the coating elastic modulus on vibration behaviors of a gear transmission system are minutely dissected by using a numerical integration approach. Numerical cases illustrate that the dynamic characteristic of a gear transmission system tends toward a one-period state under the given operating condition. They also indicate that, compared with softer coatings, sti ﬀ er ones can properly enhance the transmission performance of the coated gear pair. Numerical results are also compared with previous studies, and can establish a theoretical basis for dynamic design and vibration control of the coated gear transmission system.


Introduction
Due to the characteristics of compact structure, stable transmission, and constant speed ratio, the gear mechanism has been widely applied in various industrial fields, such as automobiles, ships, and aircrafts. With the progression of these fields, the geared system is developing in the direction of high speed, heavy load, and lightweight. It brings a new research task for developing innovative materials, in which the surface-strengthening technology is a key. Currently, coat plating has been extensively applied in gear components to obtain higher power densities and fatigue lives [1][2][3]. On the other hand, with the increasingly strict restriction on transmission noise, the prediction and control of its vibration have attracted tremendous attention. Although the coating can improve the load-carrying performance of gear surfaces, how will it affect the vibration behavior of the system? The vibration mechanism of a coated gear is quite complex, which is intimately related to the gear dynamic behavior [1,[4][5][6]. Hence, it is essential to constitute an appropriate and accurate mathematical model to obtain the nonlinear dynamic characteristic for a coated gear transmission system, which can provide theoretical supports for its dynamic design and vibration control.

Nonlinear Dynamic Model
The engagement process of a spur gear pair is analyzed on a rigid system, which is shown in Figure 1. The tooth of the driving wheel (pinion) contacts the mating tooth of the driven one (gear) along the action line N p N g , which is the inner common tangent of two base circles. Due to the meshing characteristics of involute gear, the actual effective action line is B p B g , limited by each addendum circle [25]. D p and D g are the critical points of the meshing area of single and double teeth. P is the pitch point where the relative sliding speed between the tooth surfaces reverses, resulting in an inversive direction of the friction force. Based on the Coulomb friction law, the friction forces on pinion (F p f ) and gear (F g f ) in the meshing tooth pairs changing in time can be described as: where i is equal to 1 and 2, representing single-tooth and double-tooth engagement areas, respectively. Theoretically, it is incorporated with the effects of both dry-contact and lubricant friction. However, lubricant friction has very minor impact on gear-pair torsional behaviors and the predicted motions are not prominent though different friction expressions [8,26]. Thus, a constant µ is still used in the current paper. λ(t) is the direction coefficient of tooth friction, expressed by: where "sgn" and "mod" are the sign and modulus functions, as defined in [8,27,28]. φ P , φ B g , and φ D p are individually the involute outspread angles at points P, B g , and D p .

Nonlinear Dynamic Model
The engagement process of a spur gear pair is analyzed on a rigid system, which is shown in Figure 1. The tooth of the driving wheel (pinion) contacts the mating tooth of the driven one (gear) along the action line NpNg, which is the inner common tangent of two base circles. Due to the meshing characteristics of involute gear, the actual effective action line is BpBg, limited by each addendum circle [25]. Dp and Dg are the critical points of the meshing area of single and double teeth. P is the pitch point where the relative sliding speed between the tooth surfaces reverses, resulting in an inversive direction of the friction force. Based on the Coulomb friction law, the friction forces on pinion ( pf F ) and gear ( gf F ) in the meshing tooth pairs changing in time can be described as: where i is equal to 1 and 2, representing single-tooth and double-tooth engagement areas, respectively. Theoretically, it is incorporated with the effects of both dry-contact and lubricant friction. However, lubricant friction has very minor impact on gear-pair torsional behaviors and the predicted motions are not prominent though different friction expressions [8,26]. Thus, a constant μ is still used in the current paper. ( ) t λ is the direction coefficient of tooth friction, expressed by: where "sgn" and "mod" are the sign and modulus functions, as defined in [8,27,28].
are individually the involute outspread angles at points P, Bg, and Dp. The friction arms of pinion and gear with ith tooth pair can be represented by: where g D φ are the involute outspread angles at point Dg and α is the pressure angle at the pitch point. By simplifying the realistic gear system into a mass-spring damping system, a lumped parameter dynamic model having six-DOF (one rotational and two translational for each gear disc) is built [29]. The model schematic diagram with the time-varying mesh stiffness, transmission error, backlash, and friction force is shown in Figure 2 The friction arms of pinion and gear with ith tooth pair can be represented by: where φ D g are the involute outspread angles at point D g and α is the pressure angle at the pitch point. By simplifying the realistic gear system into a mass-spring damping system, a lumped parameter dynamic model having six-DOF (one rotational and two translational for each gear disc) is built [29]. The model schematic diagram with the time-varying mesh stiffness, transmission error, backlash, and friction force is shown in Figure 2, where the xand y-axes are individually parallel to the off-line of action (OLOA) and line of action (LOA) directions of the gear pair. It is composed of a spur gear pair with masses (m p ,m g ), moments of inertia (I p ,I g ), and base circle radii (r p ,r g ), and the two gears are supported on roller bearings. The dynamic equations mentioned above contain time-varying mesh stiffness, friction, backlash, and static transmission error, which causes the strong nonlinear problem with parametric excitation in the system. A fourth-to fifth-order Runge-Kutta algorithm with a fixed time-step ( ω ) is performed to numerically integrate these coupled nonlinear differential formulas [31].
The total integration time is supposed to be ω . Thus, the transient motion can be eliminated and the steady-state forced response can be gained.

Calculation of Time Varying Mesh Stiffness
To find the connection between the rectilinear and the torsional mesh stiffness for gears, a method from [15,18] by using FE method is introduced. The current numerical calculations are performed via FE analysis code ABAQUS under a quasistatic condition. Table 1 lists parameters of a pair of Forschungsstelle für Zahnräder und Getriebebau (FZG) standard gear. The bore diameters of the pinion and gear are 45 and 65 mm, respectively. Since the contact ratio is 1.33, the three meshing teeth pairs of the FE model are applied to calculate the mesh stiffness. In order to decrease the element number in the model and reduce the calculation workload, it is necessary to refine the mesh locally. As seen in Figure 3, the gear-bore surface is fixed, and a torque T = ⋅ is applied on the body of the pinion. There are 302,075 linear hexahedral elements of type C3D8R in the model. The total number of nodes is 380,284, and the average aspect ratio of the elements is 5.59. The gear and pinion surfaces are guaranteed to touch each other at any time. The thickness of the gear coating is 55 m [32], and the coating is perfectly bonded on the gear substrate. The coefficient of friction (COF) between the coated gear surfaces is μ . The value of COF has been experimentally estimated and described afterwards.  Both the torsional and transversal motions in the system are considered in the proposed nonlinear dynamic model. In light of the second Newtonian law, the equations of motion in the x, y, and θ directions for the pinion and gear can be individually formulated as follows: where an overdot means differentiation with respect to time t.
The damping mechanisms in the gear mesh and bearings are assumed to have a linear feature, so the meshing forces are described in the following form: Coatings 2020, 10, 39 where k m (t) is the time-varying mesh stiffness, which will be discussed in Section 3. c m is the mesh damping which is proportional to the integral average stiffness k mm and can be evaluated approximately using the following equation [23,30]: where ξ is the damping ratio (0.03-0.17).
The nonlinear displacement function f (δ) can be expressed as the following: where b is half of the total gear backlash. δ is the relative displacement along the LOA, which is also known as the dynamic transmission error (DTE), defined by: where e(t) denotes the static transmission error (STE). It can commonly be approximated as a harmonic function: where ω e is the fundamental frequency of STE, ω e = ω p z p . The dynamic equations mentioned above contain time-varying mesh stiffness, friction, backlash, and static transmission error, which causes the strong nonlinear problem with parametric excitation in the system. A fourth-to fifth-order Runge-Kutta algorithm with a fixed time-step (π/(150ω p )) is performed to numerically integrate these coupled nonlinear differential formulas [31]. The total integration time is supposed to be 100π/ω p . Thus, the transient motion can be eliminated and the steady-state forced response can be gained.

Calculation of Time Varying Mesh Stiffness
To find the connection between the rectilinear and the torsional mesh stiffness for gears, a method from [15,18] by using FE method is introduced. The current numerical calculations are performed via FE analysis code ABAQUS under a quasistatic condition. Table 1 lists parameters of a pair of Forschungsstelle für Zahnräder und Getriebebau (FZG) standard gear. The bore diameters of the pinion and gear are 45 and 65 mm, respectively. Since the contact ratio is 1.33, the three meshing teeth pairs of the FE model are applied to calculate the mesh stiffness. In order to decrease the element number in the model and reduce the calculation workload, it is necessary to refine the mesh locally. As seen in Figure 3, the gear-bore surface is fixed, and a torque T p = 60 N · m is applied on the body of the pinion. There are 302,075 linear hexahedral elements of type C3D8R in the model. The total number of nodes is 380,284, and the average aspect ratio of the elements is 5.59. The gear and pinion surfaces are guaranteed to touch each other at any time. The thickness of the gear coating is 55 µm [32], and the coating is perfectly bonded on the gear substrate. The coefficient of friction (COF) between the coated gear surfaces is µ. The value of COF has been experimentally estimated and described afterwards.   Table 2. For a DLC coating, the multiple variety of its structures as well as its diamond bond (sp 3 ) and graphitic bond (sp 2 ) content lead to a diversity of mechanical properties. Its elastic modulus fluctuates from 60 to 650 GPa [33]. To get the COF, the test is performed on a four-ball tribometer [34] using GL-5 heavy-duty extreme pressure gear oil (85w-90) as the lubricant. The friction curve of the DLC-coated sample is shown in Figure 4. It can be found that the COF during the stabilization period is individually about 0.067, which is imported into the above numerical calculation as a constant. The angular deflections of the measurement points on the end surface circle of the pinion bore are obtained. The angular displacement of the pinion body is denoted by the minimal angular deflection p ϑ . The rectilinear mesh stiffness can be evaluated as follows: The time-varying mesh stiffness is calculated by selecting 20 meshing points (evenly distributed) from the tooth root to top. To verify the correctness of the current FE method, its result is compared with that of the commercial software KISSsoft (Gleason Company in Bubikon,  Table 2. For a DLC coating, the multiple variety of its structures as well as its diamond bond (sp 3 ) and graphitic bond (sp 2 ) content lead to a diversity of mechanical properties. Its elastic modulus fluctuates from 60 to 650 GPa [33]. To get the COF, the test is performed on a four-ball tribometer [34] using GL-5 heavy-duty extreme pressure gear oil (85w-90) as the lubricant. The friction curve of the DLC-coated sample is shown in Figure 4. It can be found that the COF during the stabilization period is individually about 0.067, which is imported into the above numerical calculation as a constant.   The curve-fitted plot of mesh stiffness variation along with the rotation angle is evinced in Figure 6. The meshing process of the gear pair inevitably appears to be alternating single and double teeth meshing. It is observed that the gear-mesh stiffness rises with the increase of coating elastic modulus. The main cause is that the deformation of the stiff coating is smaller than that of the soft one under the same load. Moreover, due to the influence of the substrate, when the coating elastic modulus is doubled, the mesh stiffness is not doubled, but the increase in amplitude is approximately the same.   The angular deflections of the measurement points on the end surface circle of the pinion bore are obtained. The angular displacement of the pinion body is denoted by the minimal angular deflection ϑ p . The rectilinear mesh stiffness can be evaluated as follows:

Analysis of the Nonlinear Dynamic Responses
The time-varying mesh stiffness is calculated by selecting 20 meshing points (evenly distributed) from the tooth root to top. To verify the correctness of the current FE method, its result is compared with that of the commercial software KISSsoft (Gleason Company in Bubikon, Switzerland) [10]. It can be seen from Figure 5 that the two groups of results have good accuracy and coincide with each other.   The curve-fitted plot of mesh stiffness variation along with the rotation angle is evinced in Figure 6. The meshing process of the gear pair inevitably appears to be alternating single and double teeth meshing. It is observed that the gear-mesh stiffness rises with the increase of coating elastic modulus. The main cause is that the deformation of the stiff coating is smaller than that of the soft one under the same load. Moreover, due to the influence of the substrate, when the coating elastic modulus is doubled, the mesh stiffness is not doubled, but the increase in amplitude is approximately the same.   The curve-fitted plot of mesh stiffness variation along with the rotation angle is evinced in Figure 6. The meshing process of the gear pair inevitably appears to be alternating single and double teeth meshing. It is observed that the gear-mesh stiffness rises with the increase of coating elastic modulus. The main cause is that the deformation of the stiff coating is smaller than that of the soft one under the same load. Moreover, due to the influence of the substrate, when the coating elastic modulus is doubled, the mesh stiffness is not doubled, but the increase in amplitude is approximately the same. Figure 6. The meshing process of the gear pair inevitably appears to be alternating single and double teeth meshing. It is observed that the gear-mesh stiffness rises with the increase of coating elastic modulus. The main cause is that the deformation of the stiff coating is smaller than that of the soft one under the same load. Moreover, due to the influence of the substrate, when the coating elastic modulus is doubled, the mesh stiffness is not doubled, but the increase in amplitude is approximately the same.

Analysis of the Nonlinear Dynamic Responses
The main parameters in the current dynamic model are cited from [27], which are shown in Table 3. The present section describes the influence of coating on the nonlinear dynamic responses. It is a fact that the gear transmission system may exhibit multiperiod, quasiperiod, and chaos motion. To demonstrate the dynamic behaviors more visibly, different diagrams, such as the time-history curve, phase curve, Poincare map, and FFT (fast Fourier transform) spectrum, are shown as follows. For the case of an uncoated gear pair (E c = E s ), it can be found from Figure 7a that the vibration displacement of the driving wheel along OLOA is symmetrical about x p = 0, with quasiperiod characteristic, and difference between each cycle is small. The phase diagram (Figure 7b) is a two-piece space zone composed of several closed curves, some narrow and some wide, and approximately symmetrical about the center point. On the Poincare map (Figure 7c), the phase points gather in a small area, and the system tends toward a single periodic response. As shown in Figure 7d, the spectral lines of the spectrum map are distributed discretely on the points of the combined frequency. There are many superharmonic components in the system, and the amplitude components of the first three harmonic frequencies are higher. From Figure 8, we can find that the vibration amplitude along OLOA decreases with the growth of the coating elastic modulus.
two-piece space zone composed of several closed curves, some narrow and some wide, and approximately symmetrical about the center point. On the Poincare map (Figure 7c), the phase points gather in a small area, and the system tends toward a single periodic response. As shown in Figure 7d, the spectral lines of the spectrum map are distributed discretely on the points of the combined frequency. There are many superharmonic components in the system, and the amplitude components of the first three harmonic frequencies are higher. From Figure 8, we can find that the vibration amplitude along OLOA decreases with the growth of the coating elastic modulus.   Figure 9 shows the time-history curve, phase curve, Poincare map, and FFT spectrum of the driving gear along LOA. It is manifest from the time-history curve that the vibration displacement along LOA fluctuates around 4.8 μm. Comparing Figure 9a with Figure 7a, the peak vibration amplitude along LOA is larger than that along the OLOA. The main reason is that the meshing force along LOA is greater than the friction along OLOA. The phase curve is a banded area which is neither repeated nor closed in the space, and it is clustered and distributed in a circular shape on both sides. It can be seen from the Poincare map that the phase points are radially distributed, but concentrated in a very narrow area, which shows that the response of the driving wheel along LOA is a one-period motion. Figure 9c demonstrates that the system also has many superharmonic components. Similar to Figure 8, Figure 10 exhibits that the vibration amplitude along LOA decreases with the increase in the coating elastic modulus.  Figure 9 shows the time-history curve, phase curve, Poincare map, and FFT spectrum of the driving gear along LOA. It is manifest from the time-history curve that the vibration displacement along LOA fluctuates around 4.8 µm. Comparing Figure 9a with Figure 7a, the peak vibration amplitude along LOA is larger than that along the OLOA. The main reason is that the meshing force along LOA is greater than the friction along OLOA. The phase curve is a banded area which is neither repeated nor closed in the space, and it is clustered and distributed in a circular shape on both sides. It can be seen from the Poincare map that the phase points are radially distributed, but concentrated in a very narrow area, which shows that the response of the driving wheel along LOA is a one-period motion. Figure 9c demonstrates that the system also has many superharmonic components. Similar to Figure 8, Figure 10 exhibits that the vibration amplitude along LOA decreases with the increase in the coating elastic modulus.
neither repeated nor closed in the space, and it is clustered and distributed in a circular shape on both sides. It can be seen from the Poincare map that the phase points are radially distributed, but concentrated in a very narrow area, which shows that the response of the driving wheel along LOA is a one-period motion. Figure 9c demonstrates that the system also has many superharmonic components. Similar to Figure 8, Figure 10 exhibits that the vibration amplitude along LOA decreases with the increase in the coating elastic modulus.  As shown in Figure 11, the DTE fluctuates between 4.5 and 9 μm. The phase diagram shows two concentrated areas, continuous but not repeated. According to the Poincare map, the system is still regarded as one-period motion. As seen from the FFT spectrum, the amplitudes of the first and second frequencies are larger, while the amplitudes of the remaining components are smaller. Figure  12 shows the DTE curve of coated gear models with different elastic moduli. When the coating elastic modulus increases, the DTE of the gear pair decreases as a whole. As shown in Figure 11, the DTE fluctuates between 4.5 and 9 µm. The phase diagram shows two concentrated areas, continuous but not repeated. According to the Poincare map, the system is still regarded as one-period motion. As seen from the FFT spectrum, the amplitudes of the first and second frequencies are larger, while the amplitudes of the remaining components are smaller. Figure 12 shows the DTE curve of coated gear models with different elastic moduli. When the coating elastic modulus increases, the DTE of the gear pair decreases as a whole.
As shown in Figure 11, the DTE fluctuates between 4.5 and 9 μm. The phase diagram shows two concentrated areas, continuous but not repeated. According to the Poincare map, the system is still regarded as one-period motion. As seen from the FFT spectrum, the amplitudes of the first and second frequencies are larger, while the amplitudes of the remaining components are smaller. Figure  12 shows the DTE curve of coated gear models with different elastic moduli. When the coating elastic modulus increases, the DTE of the gear pair decreases as a whole.  The dynamic meshing forces of coated gear pairs are influenced by the DTE and mesh stiffness based on Equation (5). Meshing forces under different coating elastic moduli are shown in Figure 13. The dynamic meshing force fluctuates around 1750 N. Its amplitude decreases with the increase of coating elastic modulus. The dynamic meshing forces of coated gear pairs are influenced by the DTE and mesh stiffness based on Equation (5). Meshing forces under different coating elastic moduli are shown in Figure 13. The dynamic meshing force fluctuates around 1750 N. Its amplitude decreases with the increase of coating elastic modulus. Meanwhile, the dynamic analysis method of the cracked gear can be applied to the coated gear, because both the crack and coating affect the time-varying mesh stiffness. Under the same external load, coatings with different elastic moduli have different time-varying mesh stiffness ( Figure 6). The mesh stiffness decreases with the increased length of the tooth root crack [25,31]. It is found from [17] that the amplitude of the vibration displacement declines with the increase of mesh stiffness, which agrees well with the result of the current study (Figures 8 and 10). From [35], the amplitude of DTE increases greatly when the impact of a crack on a gear body is involved. It verifies the law of the effect of the coating elastic modulus on DTE in the present work ( Figure 12). Chen et al. [36] demonstrated that the tooth root crack could bring serious oscillations to the gear dynamic meshing force. This means that the fluctuation situation of the dynamic meshing force becomes more severe with the decrease of mesh stiffness, which is consistent with the conclusion of the research ( Figure 13).

Conclusions
The nonlinear dynamic characteristic of a DLC-coated gear transmission system is analyzed by using a numerical integration approach. The effect of the coating elastic modulus on the dynamic response is investigated. Several illustrative numerical cases are introduced. According to the results, the following conclusions are put forward. Meanwhile, the dynamic analysis method of the cracked gear can be applied to the coated gear, because both the crack and coating affect the time-varying mesh stiffness. Under the same external load, coatings with different elastic moduli have different time-varying mesh stiffness ( Figure 6). The mesh stiffness decreases with the increased length of the tooth root crack [25,31]. It is found from [17] that the amplitude of the vibration displacement declines with the increase of mesh stiffness, which agrees well with the result of the current study (Figures 8 and 10). From [35], the amplitude of DTE increases greatly when the impact of a crack on a gear body is involved. It verifies the law of the effect of the coating elastic modulus on DTE in the present work ( Figure 12). Chen et al. [36] demonstrated that the tooth root crack could bring serious oscillations to the gear dynamic meshing force. This means that the fluctuation situation of the dynamic meshing force becomes more severe with the decrease of mesh stiffness, which is consistent with the conclusion of the research (Figure 13).

Conclusions
The nonlinear dynamic characteristic of a DLC-coated gear transmission system is analyzed by using a numerical integration approach. The effect of the coating elastic modulus on the dynamic response is investigated. Several illustrative numerical cases are introduced. According to the results, the following conclusions are put forward.
Under the given operating condition, the dynamic characteristic of a gear transmission system tends toward a one-period state, and there are many superharmonic components in the system. With the increase of coating elastic modulus, the gear-mesh stiffness and damping increase. The amplitudes of vibration displacements and velocities along OLOA and LOA decrease, respectively. Also, the oscillating curves of DTE move down as a whole. Moreover, increasing the coating elastic modulus can effectively reduce the amplitude of dynamic meshing force. It indicates that compared with softer coatings, stiffer ones can properly improve the transmission performance of the coated gear system.