Effects of Tooth Surface Crack Propagation on Meshing Stiffness and Vibration Characteristic of Spur Gear System

: Tooth surface cracks are considered as the early stage of the development of tooth surface spalling failure. Understanding the excitation mechanism of surface cracks has a great signiﬁcance in the early diagnosis of spalling faults. However, there are few studies on the dynamic modelling of surface cracks, and the inﬂuence mechanism of surface cracking on the dynamic characteristics of a gear system is also not yet clear during its propagation process. Thus, an analytical calculation model of the meshing stiffness of gear with tooth surface crack is developed. Then, a dynamic model of a spur gear system with six degrees of freedom (DOF) is established based on the proposed surface crack calculation model. The effects of surface crack propagation on the meshing stiffness and dynamic characteristics of gear system are investigated. The results show that the side frequencies of dynamic transmission error (DTE) are more sensitive than those of the acceleration responses during the surface crack propagation, which is more favorable to the surface crack fault diagnosis. Compared to the traditional spalling fault model, the proposed model can accurately characterize the dynamic characteristics of a gear system with the early spalling defect.


Introduction
Gearboxes are widely applied in the power transmission field [1][2][3]. However, as an vital component of a gearbox, gears often operate in stressful conditions, such as excessive service load, poor lubrication, or high temperature [4], and thus some faults related to gears can be easily observed. Spalling is one of the typical gear tooth faults, which may lead to abnormal dynamic performance of the gearbox due to the reduction in gear mesh stiffness [5][6][7]. Therefore, the mesh stiffness computation models of gears with or without a spalling fault has garnered the increased attention of researchers in the field of fault diagnosis of gear systems [8][9][10].
Compared to the finite element model, many researchers preferred the analytical calculation model of mesh stiffness due to its faster calculation [11,12]. Weber [13] calculated the deformation of a spur gear under load using the theory of displacement superposition. Cornell [14] improved the Weber's model by considering the influence of fillet foundation on the gear. Kasuba and Evans [15] calculated the meshing stiffness of a healthy gear using the digital method. Yang and Lin [16] then developed the potential energy method, which was diffusely applied in the computation of gear meshing stiffness. Tian [17] further optimized the model in [16], considering shear stiffness. Based on the gear tooth fillet foundation deflection correction formula put forward by Sainsot [18], Chen and Shao [19] improved the original model by considering the deformation of tooth fillet foundation in the calculation model of the gear's comprehensive mesh stiffness.
Aiming to study the dynamical performance of the gear system with spalling failures, some spalling fault modelling methods have been proposed [20,21]. For instance, Chaari Appl. Sci. 2021, 11,1968 2 of 18 et al. [22] proposed an analytical model for computing the meshing stiffness of a spalled gear, and then discussed the influences of spalling parameters on the meshing stiffness. Ma et al. [23,24] presented a dynamical model of a spalled gear pair and investigated the effect of spalling failure on the mesh stiffness and dynamic performance of a gear system. Based on the slicing method, Han [25] studied the influence of a spalling defect on meshing stiffness, and Shao et al. [26] put forward a meshing stiffness calculation method for a gear with a spalled fault. Fakhfakh et al. [27] investigated the vibration characteristics of a gear system with a spalled defect by establishing the dynamical model of the gear system, and then verified the accuracy of the simulation results through experiments. Ma et al. [28] presented an improved mesh stiffness computation method for gears with a spalling fault by considering the influence of extended tooth contact. Saxena et al. [29] proposed a meshing stiffness calculation model of a spalling gear, and studied the influence of several spalled shapes, positions, and sliding friction on mesh stiffness. Jiang et al. [30] investigated the vibration performance of a helical gear system with spalled failure, considering the influence of sliding friction. Yu et al. [21] analyzed the influence of a nonlinear elliptical tooth surface contact pattern of a spalled gear on the dynamical response of a gear system, and validated the accuracy of the developed model through experiments. Luo et al. [31] put forward a new dynamical simulation model of a spur gear and verified the effectiveness of the model through various experimental tests, and then analyzed the impact of a spalled defect on the dynamical response of a gear system.
However, spalling faults are usually caused by the metal shedding of the tooth surface due to the gradual propagation of the surface crack under the action of surface stress [32,33]. As an early fault of the spalling failure, the surface crack fault is different from the spalling fault in its influence on meshing stiffness and the dynamic characteristic of the gear system. There exist few studies on the dynamic modelling of surface cracks, and the influence mechanism of a surface crack fault on the dynamic characteristics of a gear system is also not yet clear during its propagation process. Therefore, an analytical calculation model of meshing stiffness of the gear with tooth surface crack fault is proposed in the paper. Then, the effects of surface crack propagation on the meshing stiffness and dynamic characteristic of gear system are investigated, which is conducive to the diagnosis and monitoring of the early spalling failure of a gear system.
The following sections of this paper are arranged as follows. The proposed calculation model of meshing stiffness of a gear with a tooth surface crack is introduced in Section 2. Then, the six-DOF dynamical model of spur gear system is established in Section 3. Furthermore, the influence of surface crack propagation on the meshing stiffness and vibration characteristics of a spur gear system is analyzed and discussed in Section 4. Finally, some conclusions are obtained in Section 5.

Proposed Meshing Stiffness Calculation Model with Tooth Surface Crack
As demonstrated in Figure 1, the tooth surface crack fault is usually located at the middle of the tooth surface and is symmetrical along the center plane of gear tooth width. The boundary of tooth surface cracks is assumed as consisting of straight lines for the convenience of modeling and calculation. L 1 represents the crack length, L 2 stands for the crack width, and D refers to the crack depth.
The potential energy method is diffusely adopted in the gear mesh stiffness calculation due to its fast calculation and simple modeling [16,17], which is also applied in this paper. The gear tooth is assumed to be a cantilever beam with a variable cross section when establishing the meshing stiffness calculation model of a gear tooth with surface cracking, which is displayed in Figure 2. The bending energy U b , shear energy U s , and axial compressive energy U a deposited in a tooth can be expressed as [16,17], where k b , k s , and k a signify the bending stiffness, shear stiffness, and axial compressive stiffness, respectively. F refers to the meshing force acting on the gear tooth. E and G are the Young's modulus and shear modulus, respectively, and it is assumed that the material of the gear is uniform and the Young's modulus and shear modulus of its surface and internal parts are the same in this paper. d denotes the distance from the meshing point to the dedendum circle, and h is the distance from the meshing point to the center line of the gear tooth. F a and F b represent the horizontal and vertical components of the meshing force F, respectively, which are expressed as, The potential energy method is diffusely adopted in the gear mesh stiffness calculation due to its fast calculation and simple modeling [16,17], which is also applied in this paper. The gear tooth is assumed to be a cantilever beam with a variable cross section when establishing the meshing stiffness calculation model of a gear tooth with surface cracking, which is displayed in Figure 2. The bending energy Ub, shear energy Us, and axial compressive energy Ua deposited in a tooth can be expressed as [16,17], Figure 2. The cantilever beam model of gear with tooth surface crack fault. The potential energy method is diffusely adopted in the gear mesh stiffness calculation due to its fast calculation and simple modeling [16,17], which is also applied in this paper. The gear tooth is assumed to be a cantilever beam with a variable cross section when establishing the meshing stiffness calculation model of a gear tooth with surface cracking, which is displayed in Figure 2. The bending energy Ub, shear energy Us, and axial compressive energy Ua deposited in a tooth can be expressed as [16,17], where kb, ks, and ka signify the bending stiffness, shear stiffness, and axial compressive stiffness, respectively. F refers to the meshing force acting on the gear tooth. E and G are the Young's modulus and shear modulus, respectively, and it is assumed that the material According to Equations (1) and (2), the stiffness K b , Ks, and K a can be deduced as, According to the characteristic of involute tooth profile of healthy gear tooth, the expressions of parameters d, x, h, h x , A x , and I x are as follows, where h x is the cross-sectional height with x distance from the dedendum. R b and R f are the base circle radius and dedendum circle radius, respectively. α represents the gear rotation angle, α 2 and α 3 refer to half of tooth angle of the base circle and dedendum circle, respectively. A x is the area of the effective section and I x represents the area moment of inertia. The meshing positions of the gear with tooth surface cracking are presented in Figure 3. Point G represents the position of dedendum, point P denotes the position at the beginning of the crack, point S refers to the position at the end of the crack, and point T stands for the position of the addendum. As displayed in Figure 3a, the parameters A x and I x remain unchanged when the contact line is located between G and S. Thus, the mesh stiffness of gear tooth pair in this area is the same as that of healthy gear teeth. As illustrated in Figure 3b, the crack part can still bear shear, axial compressive, and bending force when the contact line is located between S and P, and the mesh stiffness of gear teeth in this area is also the same as that of normal gear. When the contact line passes through the end point P of the crack, the crack turns to an open state due to the presence of meshing force, and the parameters A x and I x of the cracked tooth zone will decrease, which is presented in Figure 3c.  Through the above analysis, the cross section area Ax' and the moment of inertia Ix' of the gear with tooth surface crack are calculated as follows, Through the above analysis, the cross section area A x ' and the moment of inertia I x ' of the gear with tooth surface crack are calculated as follows, where D' represents the reduced height of the cross section due to the tooth surface crack, it can be denoted as, where, h x1 is the cross-sectional height with x 1 distance from the dedendum. According to the research conclusions of Yang and Sun, the magnitude of the Hertz contact stiffness k h of healthy meshing tooth pairs is a constant throughout the meshing line. For the meshing gear teeth with surface crack, the crack will not affect the effective contact width of the gear teeth during the whole meshing process. Therefore, the Hertz contact stiffness of surface cracked teeth is consistent with the healthy case. Its expression is [34], where W and ν represent tooth width and Poisson's ratio, respectively. Under the action of meshing force, deformations of the tooth fillet-foundation also affect the gear meshing stiffness. The formula widely adopted to compute the filletfoundation deformation of gear tooth was developed by Sainsot et al. in Ref. [18]. It is indicated as, where k f denotes the tooth fillet-foundation stiffness. The parameters µ f and S f are displayed in Figure 2 in [18], and the parameters L, Q, M, and P are functions of θ f and h f , which can be presented as [18], where the meaning of symbols h f , θ f , A i , B i , C i , D i , E i , and F i can be observed in [18]. Finally, the comprehensive meshing stiffness of one tooth pair can be calculated as, where the subscripts 1 and 2 mean the pinion (driving gear) and gear (driven gear), respectively. The comprehensive meshing stiffness of double tooth pairs is deduced as, damping of the bearing in the x and y direction, respectively. The subscript i = p, g refers to the pinion and gear, respectively.  The displacements of the pinion and gear are projected onto the meshing line, then the relative displacement δ between the pinion and gear can be obtained, where θp and θg denote the angular displacement of the pinion and gear, respectively. xp and xg represent the lateral displacement of the pinion and gear along the x direction, respectively. yp and yg refer to the lateral displacement of the pinion and gear along the y direction, respectively.
The motion equation of the driven gear is indicated as follows, DTE denotes the error between the actual position and theoretical position of the gear under the action of dynamic mesh force. The characteristics of DTE have a direct impact on the vibration, noise and working stability of gear system. The calculation formula of DTE is expressed as,

Results and Discussions
According to the parameters of a spur gear system shown in Table 1, the meshing stiffness and dynamic responses of the gear system with surface crack is calculated. Then, The displacements of the pinion and gear are projected onto the meshing line, then the relative displacement δ between the pinion and gear can be obtained, where θ p and θ g denote the angular displacement of the pinion and gear, respectively. x p and x g represent the lateral displacement of the pinion and gear along the x direction, respectively. y p and y g refer to the lateral displacement of the pinion and gear along the y direction, respectively.
The motion equation of the pinion can be obtained as follows, The motion equation of the driven gear is indicated as follows, DTE denotes the error between the actual position and theoretical position of the gear under the action of dynamic mesh force. The characteristics of DTE have a direct impact on the vibration, noise and working stability of gear system. The calculation formula of DTE is expressed as,

Results and Discussions
According to the parameters of a spur gear system shown in Table 1, the meshing stiffness and dynamic responses of the gear system with surface crack is calculated. Then, the effects of surface crack propagation on the meshing stiffness and vibration characteristics are studied and discussed. As displayed in Figure 5, taking the surface crack parameters as L 1 = 1, 2, 3 mm, L 2 = 40 mm, D = 3 mm, and x 1 = 2 mm, the influence of the crack length (L 1 ) on the meshing stiffness of single-tooth pair and double-tooth pairs is studied. It can be obtained that the total meshing stiffness of the cracked gear is smaller than that of the healthy gear at the beginning of the meshing process. The stiffness drop zone is enlarged with the increase in the crack length, while the decrease in the magnitude of the meshing stiffness under different crack lengths remains the same. The maximum ratios of the stiffness reduction of single and double tooth pairs are 8.2% and 7.0%, respectively. the effects of surface crack propagation on the meshing stiffness and vibration characteristics are studied and discussed. As displayed in Figure 5, taking the surface crack parameters as L1 = 1, 2, 3 mm, L2 = 40 mm, D = 3 mm, and x1 = 2 mm, the influence of the crack length (L1) on the meshing stiffness of single-tooth pair and double-tooth pairs is studied. It can be obtained that the total meshing stiffness of the cracked gear is smaller than that of the healthy gear at the beginning of the meshing process. The stiffness drop zone is enlarged with the increase in the crack length, while the decrease in the magnitude of the meshing stiffness under different crack lengths remains the same. The maximum ratios of the stiffness reduction of single and double tooth pairs are 8.2% and 7.0%, respectively.  As displayed in Figure 6, the influence of the crack width (L2) on the meshing stiffness of single-tooth pair and double-tooth pair is investigated when the surface crack parameters are set as L1 = 2 mm, L2 = 15, 30, 45 mm, D = 3 mm, and x1 = 2 mm. It can be observed As displayed in Figure 6, the influence of the crack width (L 2 ) on the meshing stiffness of single-tooth pair and double-tooth pair is investigated when the surface crack parameters are set as L 1 = 2 mm, L 2 = 15, 30, 45 mm, D = 3 mm, and x 1 = 2 mm. It can be observed from the calculation results that the reduction in the magnitude of the meshing stiffness gradually goes up with the increase in the crack width, but the mesh period where the meshing stiffness decreases remain unchanged. The maximum ratios of the stiffness reduction of single and double tooth pairs are 13.0% and 8.7%, respectively. from the calculation results that the reduction in the magnitude of the meshing stiffness gradually goes up with the increase in the crack width, but the mesh period where the meshing stiffness decreases remain unchanged. The maximum ratios of the stiffness reduction of single and double tooth pairs are 13.0% and 8.7%, respectively. The effects of the crack depth (D) on the meshing stiffness of the single-tooth pair and double-tooth pair are analyzed when setting the surface crack parameters as L1=2 mm, L2 = 40 mm, D = 1, 2, 3 mm and x1 = 2 mm. We can find in Figure 7 that the reduction in the magnitude of the meshing stiffness is gradually enlarged when the crack depth increases, while the areas where the mesh stiffness decreases remain the same. The maximum ratios of the stiffness reduction of single and double tooth pairs are 10.8% and 7.2%, respectively.

Effects of a Single Crack Parameter on DTE
The influence of a single crack parameter on the DTE is investigated and presented in Figure 8. The variable parameters of surface cracks used in this study are the same as those in Section 4.1.1. The input parameters ωp and Tp are set to a fixed value equal to 900 rpm and 100 Nm, respectively. As shown in Figure 8a, the occurrence of cracks will increase the amplitude of DTE responses. The greater the crack length, the earlier the DTE value returns to the healthy tooth level. It can be discovered from Figure 8b that in the area affected by the surface crack failure, the value of DTE goes up with the increase in the crack width, and the DTE value is highest when L2 = 45 mm. It can also be seen in The effects of the crack depth (D) on the meshing stiffness of the single-tooth pair and double-tooth pair are analyzed when setting the surface crack parameters as L 1 = 2 mm, L 2 = 40 mm, D = 1, 2, 3 mm and x 1 = 2 mm. We can find in Figure 7 that the reduction in the magnitude of the meshing stiffness is gradually enlarged when the crack depth increases, while the areas where the mesh stiffness decreases remain the same. The maximum ratios of the stiffness reduction of single and double tooth pairs are 10.8% and 7.2%, respectively. from the calculation results that the reduction in the magnitude of the meshing stiffness gradually goes up with the increase in the crack width, but the mesh period where the meshing stiffness decreases remain unchanged. The maximum ratios of the stiffness reduction of single and double tooth pairs are 13.0% and 8.7%, respectively. The effects of the crack depth (D) on the meshing stiffness of the single-tooth pair and double-tooth pair are analyzed when setting the surface crack parameters as L1=2 mm, L2 = 40 mm, D = 1, 2, 3 mm and x1 = 2 mm. We can find in Figure 7 that the reduction in the magnitude of the meshing stiffness is gradually enlarged when the crack depth increases, while the areas where the mesh stiffness decreases remain the same. The maximum ratios of the stiffness reduction of single and double tooth pairs are 10.8% and 7.2%, respectively.

Effects of a Single Crack Parameter on DTE
The influence of a single crack parameter on the DTE is investigated and presented in Figure 8. The variable parameters of surface cracks used in this study are the same as those in Section 4.1.1. The input parameters ωp and Tp are set to a fixed value equal to 900 rpm and 100 Nm, respectively. As shown in Figure 8a, the occurrence of cracks will increase the amplitude of DTE responses. The greater the crack length, the earlier the DTE value returns to the healthy tooth level. It can be discovered from Figure 8b that in the area affected by the surface crack failure, the value of DTE goes up with the increase in the crack width, and the DTE value is highest when L2 = 45 mm. It can also be seen in

Effects of a Single Crack Parameter on DTE
The influence of a single crack parameter on the DTE is investigated and presented in Figure 8. The variable parameters of surface cracks used in this study are the same as those in Section 4.1.1. The input parameters ω p and T p are set to a fixed value equal to 900 rpm and 100 Nm, respectively. As shown in Figure 8a, the occurrence of cracks will increase the amplitude of DTE responses. The greater the crack length, the earlier the DTE value returns to the healthy tooth level. It can be discovered from Figure 8b that in the area affected by the surface crack failure, the value of DTE goes up with the increase in the crack width, and the DTE value is highest when L 2 = 45 mm. It can also be seen in Figure 8c

Effects of a Single Crack Parameter on Acceleration Response
As revealed in Figure 9, the influence of a single crack parameter on the time domain response of acceleration is investigated. It can be observed that the maximum acceleration response decreases gradually with the increase in crack length, while the acceleration response is positively correlated with the length and depth of the surface crack. In an overall view, the influence of the single crack parameter on the time-domain responses of acceleration are not as obvious as the time-domain responses of DTE.

Effects of a Single Crack Parameter on Acceleration Response
As revealed in Figure 9, the influence of a single crack parameter on the time domain response of acceleration is investigated. It can be observed that the maximum acceleration response decreases gradually with the increase in crack length, while the acceleration response is positively correlated with the length and depth of the surface crack.

Effects of Surface Crack Propagation on Meshing Stiffness and Vibration Characteristics
The influence of a single crack parameter on meshing stiffness and vibration characteristics of the gear system is discussed in Section 4.1. In actuality, the depth, width, and length of the surface crack will change during the crack propagation process. Therefore, the influence of the surface crack propagation on meshing stiffness and dynamic characteristics of a gear system is also studied and discussed in the following chapters. In our study, it is assumed that the width of the surface crack remains unchanged when the crack grows, and the crack propagates only in the length and depth direction. As presented in Figure 10, the red dotted line means the surface crack propagation path, point P is the crack propagation initial position, and point Si (i = 1, 2, 3) denotes the crack vertex position in the process of crack propagation. In this study, the spalling failure is assumed to occur when the surface crack propagates to the position S3. Specific parameters of the surface crack propagation case are displayed in Table 2.

Effects of Surface Crack Propagation on Meshing Stiffness and Vibration Characteristics
The influence of a single crack parameter on meshing stiffness and vibration characteristics of the gear system is discussed in Section 4.1. In actuality, the depth, width, and length of the surface crack will change during the crack propagation process. Therefore, the influence of the surface crack propagation on meshing stiffness and dynamic characteristics of a gear system is also studied and discussed in the following chapters. In our study, it is assumed that the width of the surface crack remains unchanged when the crack grows, and the crack propagates only in the length and depth direction. As presented in Figure 10, the red dotted line means the surface crack propagation path, point P is the crack propagation initial position, and point S i (i = 1, 2, 3) denotes the crack vertex position in the process of crack propagation. In this study, the spalling failure is assumed to occur when the surface crack propagates to the position S 3 . Specific parameters of the surface crack propagation case are displayed in Table 2. Appl. Sci. 2021, 11, x FOR PEER REVIEW 12 of 19 The crack propagation path Figure 10. Schematic of surface crack propagation. The effects of surface crack propagation on the mesh stiffness are shown in Figure 11. It can be found that the meshing stiffness of the gear system reduces gradually with the propagation of the surface crack, and the decrease in mesh stiffness is greater when the crack evolves into a spalling fault. It can be concluded that the spalling fault has a larger affect area and amplitude on the meshing stiffness compared with the surface crack fault.

Effects of the Surface Crack Propagation Progress on DTE
As presented in Figure 12, the DTE value of the spalling fault is obviously greater than that of the crack fault, which is due to the meshing stiffness of the spalling fault on the gear tooth being smaller than that of the crack fault. Meanwhile, the effects of crack fault cases on the time domain of DTE are not obvious before the crack evolves into a spalling fault.

Effects of Surface Crack Propagation Progress on Meshing Stiffness
The effects of surface crack propagation on the mesh stiffness are shown in Figure 11. It can be found that the meshing stiffness of the gear system reduces gradually with the propagation of the surface crack, and the decrease in mesh stiffness is greater when the crack evolves into a spalling fault. It can be concluded that the spalling fault has a larger affect area and amplitude on the meshing stiffness compared with the surface crack fault.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 12 of 19 The crack propagation path Figure 10. Schematic of surface crack propagation. The effects of surface crack propagation on the mesh stiffness are shown in Figure 11. It can be found that the meshing stiffness of the gear system reduces gradually with the propagation of the surface crack, and the decrease in mesh stiffness is greater when the crack evolves into a spalling fault. It can be concluded that the spalling fault has a larger affect area and amplitude on the meshing stiffness compared with the surface crack fault.

Effects of the Surface Crack Propagation Progress on DTE
As presented in Figure 12, the DTE value of the spalling fault is obviously greater than that of the crack fault, which is due to the meshing stiffness of the spalling fault on the gear tooth being smaller than that of the crack fault. Meanwhile, the effects of crack fault cases on the time domain of DTE are not obvious before the crack evolves into a spalling fault.

Effects of the Surface Crack Propagation Progress on DTE
As presented in Figure 12, the DTE value of the spalling fault is obviously greater than that of the crack fault, which is due to the meshing stiffness of the spalling fault on the gear tooth being smaller than that of the crack fault. Meanwhile, the effects of crack fault cases on the time domain of DTE are not obvious before the crack evolves into a spalling fault.  The effects of surface crack propagation on DTE spectra are also studied and presented in Figure 13. It is observed that the side frequency components appear near the mesh frequency fm and its harmonic frequencies (2fm, 3fm, …) in the spectrum diagrams of four fault cases, and the amplitudes of side frequency go up gradually with the surface crack propagation. In addition, the interval of the two adjacent side frequency components equals to the rotation frequency fn of the driving gear. It can be concluded that the surface crack propagation has a small influence on the time domain of DTE, but an obvious influence on the sideband components of the DTE spectrum, which is an important feature for the diagnosis of early spalling failure. The effects of surface crack propagation on DTE spectra are also studied and presented in Figure 13. It is observed that the side frequency components appear near the mesh frequency f m and its harmonic frequencies (2f m , 3f m , . . . ) in the spectrum diagrams of four fault cases, and the amplitudes of side frequency go up gradually with the surface crack propagation. In addition, the interval of the two adjacent side frequency components equals to the rotation frequency f n of the driving gear. It can be concluded that the surface crack propagation has a small influence on the time domain of DTE, but an obvious influence on the sideband components of the DTE spectrum, which is an important feature for the diagnosis of early spalling failure.  The effects of surface crack propagation on DTE spectra are also studied and presented in Figure 13. It is observed that the side frequency components appear near the mesh frequency fm and its harmonic frequencies (2fm, 3fm, …) in the spectrum diagrams of four fault cases, and the amplitudes of side frequency go up gradually with the surface crack propagation. In addition, the interval of the two adjacent side frequency components equals to the rotation frequency fn of the driving gear. It can be concluded that the surface crack propagation has a small influence on the time domain of DTE, but an obvious influence on the sideband components of the DTE spectrum, which is an important feature for the diagnosis of early spalling failure.

Effects of the Surface Crack Propagation Progress on Acceleration Response
The effects of surface crack propagation on the time domain of acceleration response are displayed in Figure 14. We can see that the spalling fault has an obvious impact compared with the crack fault, which results from the larger changes in meshing stiffness when the spalling gear just enters and exits the spalling area. It can also be seen that the influence of the tooth surface crack fault on the time domain of acceleration responses is also not obvious before the crack evolves into a spalling fault.  Figure 15 shows the influence of surface crack propagation on acceleration spectra. It can be found that the side frequency components also appear near the mesh frequency and its harmonic frequencies in the acceleration spectra of four fault cases. The side frequency bands change significantly when the crack developed into a spalling fault. Meanwhile, compared with the DTE spectra, the increased amplitudes of the side frequencies are not obvious before the surface crack propagates into the spalling fault.

Effects of the Surface Crack Propagation Progress on Acceleration Response
The effects of surface crack propagation on the time domain of acceleration response are displayed in Figure 14. We can see that the spalling fault has an obvious impact compared with the crack fault, which results from the larger changes in meshing stiffness when the spalling gear just enters and exits the spalling area. It can also be seen that the influence of the tooth surface crack fault on the time domain of acceleration responses is also not obvious before the crack evolves into a spalling fault.

Effects of the Surface Crack Propagation Progress on Acceleration Response
The effects of surface crack propagation on the time domain of acceleration response are displayed in Figure 14. We can see that the spalling fault has an obvious impact compared with the crack fault, which results from the larger changes in meshing stiffness when the spalling gear just enters and exits the spalling area. It can also be seen that the influence of the tooth surface crack fault on the time domain of acceleration responses is also not obvious before the crack evolves into a spalling fault.  Figure 15 shows the influence of surface crack propagation on acceleration spectra. It can be found that the side frequency components also appear near the mesh frequency and its harmonic frequencies in the acceleration spectra of four fault cases. The side frequency bands change significantly when the crack developed into a spalling fault. Meanwhile, compared with the DTE spectra, the increased amplitudes of the side frequencies are not obvious before the surface crack propagates into the spalling fault.  Figure 15 shows the influence of surface crack propagation on acceleration spectra. It can be found that the side frequency components also appear near the mesh frequency and its harmonic frequencies in the acceleration spectra of four fault cases. The side frequency bands change significantly when the crack developed into a spalling fault. Meanwhile, compared with the DTE spectra, the increased amplitudes of the side frequencies are not obvious before the surface crack propagates into the spalling fault. From the above results, we can find that sidebands of the spectra resulting from the crack propagation are more sensitive than the time domain responses. Thus, the sideband amplitude ratio (BAR) is introduced to further quantitatively analyze the effects of the tooth surface crack propagation on the spectra of DTE and acceleration. The calculation formula of the BAR is as follows [35], From the above results, we can find that sidebands of the spectra resulting from the crack propagation are more sensitive than the time domain responses. Thus, the sideband amplitude ratio (BAR) is introduced to further quantitatively analyze the effects of the tooth surface crack propagation on the spectra of DTE and acceleration. The calculation formula of the BAR is as follows [35], where A represents the amplitude of the ith side frequency on either side of the harmonic frequency f m . k denotes the number of side frequency. From Equation (21) we can learn that the higher the BAR is, the richer the sideband frequency components are. The BAR values of DTE and acceleration under different fault cases are computed and displayed in Table 3. It can be observed that the average BARs of DTE and acceleration tend to increase with the surface crack propagation, and reaches maximum after the spalling failure occurs. From Table 3, it can be calculated that the growth rates of the average DTE BARs under case #2, case #3, and case #4 are 1.00, 2.98, and 500.19, respectively, when compared to case #1, while the average BARs of acceleration are 0.26, 1.00, and 30.62, respectively. It can be concluded that the BAR of DTE increases faster and changes more obviously than that of the acceleration during the surface crack propagation, which is more conducive to the surface crack fault diagnosis.   Figure 16 shows the BAR changes of DTE and acceleration under different harmonic orders and fault cases. It can be found that BARs of DTE and acceleration tend to increase with the surface crack propagation at each harmonic frequency. The BAR value of DTE increases obviously at each harmonic frequency, while the BAR value of acceleration increases significantly only at the harmonic frequency of order 5-10. Therefore, compared with the acceleration response, the side band variations of DTE during the surface crack propagation are more propitious to the fault detection of a tooth surface crack.

Conclusions
The surface crack of a gear tooth is considered as an early fault that causes the spal defect. However, there has been few discussions on the dynamic modelling of surf cracking. Thus, an analytical calculation model of the mesh stiffness of gear with to surface crack is developed. Then, based on the presented model, the effects of the surf crack propagation on the meshing stiffness and dynamic characteristics of gear system investigated and discussed. Several conclusions are obtained: (1) The mesh stiffness of gear system will decrease gradually with the surface cra propagation, and it decreases significantly when the crack propagation evolves int spalling fault.
(2) Compared with the time domain responses of DTE and acceleration, the sideba frequencies in the frequency domain responses are more sensitive to the surface cra propagation.
(3) The side frequencies of DTE increase faster and change more obviously than th of the acceleration during the surface crack propagation, which is more conducive to surface crack fault diagnosis.

Conclusions
The surface crack of a gear tooth is considered as an early fault that causes the spalled defect. However, there has been few discussions on the dynamic modelling of surface cracking. Thus, an analytical calculation model of the mesh stiffness of gear with tooth surface crack is developed. Then, based on the presented model, the effects of the surface crack propagation on the meshing stiffness and dynamic characteristics of gear system are investigated and discussed. Several conclusions are obtained: (1) The mesh stiffness of gear system will decrease gradually with the surface crack propagation, and it decreases significantly when the crack propagation evolves into a spalling fault.
(2) Compared with the time domain responses of DTE and acceleration, the sideband frequencies in the frequency domain responses are more sensitive to the surface crack propagation.
(3) The side frequencies of DTE increase faster and change more obviously than those of the acceleration during the surface crack propagation, which is more conducive to the surface crack fault diagnosis.
It is worth mentioning that the calculation accuracy of the proposed model has not been verified in this paper. Therefore, the authors hope to use the finite element method, experimental method, or other effective methods to verify its computational accuracy in future research work, and finally apply it to the field of early spalling fault diagnosis.