Contact stiffness and damping of spiral bevel gears under transient mixed lubrication conditions

Existing studies primarily focus on stiffness and damping under full-film lubrication or dry contact conditions. However, most lubricated transmission components operate in the mixed lubrication region, indicating that both the asperity contact and film lubrication exist on the rubbing surfaces. Herein, a novel method is proposed to evaluate the time-varying contact stiffness and damping of spiral bevel gears under transient mixed lubrication conditions. This method is sufficiently robust for addressing any mixed lubrication state regardless of the severity of the asperity contact. Based on this method, the transient mixed contact stiffness and damping of spiral bevel gears are investigated systematically. The results show a significant difference between the transient mixed contact stiffness and damping and the results from Hertz (dry) contact. In addition, the roughness significantly changes the contact stiffness and damping, indicating the importance of film lubrication and asperity contact. The transient mixed contact stiffness and damping change significantly along the meshing path from an engaging-in to an engaging-out point, and both of them are affected by the applied torque and rotational speed. In addition, the middle contact path is recommended because of its comprehensive high stiffness and damping, which maintained the stability of spiral bevel gear transmission.


Introduction
As crucial mechanical transmission parts, spiral bevel gears are widely used in aviation, automobile, and mining. They typically operate under heavy loads and high speeds owing to their favorable performances in terms of load carrying capacity and transmission stability. Meanwhile, lubricants are typically introduced to satisfy the demands of low-wear and long-life transmission of spiral bevel gears. In fact, direct surface contacts in some areas are inevitable because the height of the surface asperity is typically equal to or greater than the order of magnitude of the film thickness [1], particularly under heavy-load and low-speed conditions, termed as "mixed lubrication". In most gear transmission systems, the squeezing action caused by the dynamic variations of contact parameters significantly affects lubrication behaviors. This is known as "transient mixed lubrication". Lubricated contacts have been the focus of research because of their ubiquitous and formative effect on tribology, vibration, noise, and dynamic response in many driving parts, such as gears, bearings, and cam followers [2,3]. As two vital factors in tribo-dynamics, the stiffness and damping of transient mixed lubrication should be measured accurately. However, unlike full-film stiffness and damping or dry contact situation, investigations regarding transient mixed contact stiffness and damping are rare. Hence, we herein present a novel method to accurately predict the stiffness and damping in transient mixed lubrication. This method is sufficiently robust for managing any transient elastohydrodynamic lubrication (EHL) contact situations regardless of the severity of the asperity contacts, as well as smooth surface contacts in a unified stiffness and damping model based on transient mixed lubrication results.
Researches regarding EHL theories and experiments are progressing. However, a more realistic numerical simulation of EHL to further consider the twodimensional machined or stochastic roughness has been presented by Venner [4] in 1976, Kweh et al. [5] in 1992, and others. Subsequently, digitized threedimensional (3D) machined rough surfaces were adopted by Xu and Sadeghi [6], Zhu and Ai [7], and others. Following the aforementioned studies, Zhu and Hu [8,9] established the first unified mixed lubrication numerical simulation model coupled with 3D machined rough surfaces in point contact. In this unified mixed lubrication model, the asperity contact was regarded as a consequence of film thickness decreasing toward zero continuously, and the transition between dry contact and EHL contact is successive. Therefore, a unified Reynolds solution strategy is sufficient for addressing both asperity and EHL contacts [10]. Several techniques have been proposed or adopted to improve this model, such as the differential scheme [11], surface elastic deformation [12,13], surface plastic deformation [14,15], entrainment angle [16,17], and progressive mesh densification (PMD) method [18]. Moreover, the transient squeezing effect of oil film caused by dynamic changes in load, contact geometry, and speed should be considered, such as during startup, shutdown [19,20], and impact [21,22], as well as in various gears [23,24] and cam followers [25,26]. For spiral bevel gears, Pu et al. [27] presented a mixed lubrication model considering both the 3D machined roughness and velocity angle between two contact surfaces; based on the numerical simulation results, the friction and flash temperature were further predicted. Subsequently, Cao et al. [28] systematically analyzed the mixed lubrication, friction, flash temperature, and fatigue performances of spiral bevel gears in heel, middle, and toe contacts. Zhang et al. [29] and Gan et al. [30] coupled a mixed lubrication model with thermal analysis commercial software to investigate the temperature characteristics of spiral bevel gears. However, the aforementioned and most available studies regarding spiral bevel gears were conducted under the assumption of steady-state EHL. Recently, Wang et al. [31] established a transient mixed lubrication model using transient PMD method for spiral bevel gears. The analysis results revealed that the steady EHL considerably overestimated the lubrication performances of spiral bevel gears relative to unsteady EHL. Subsequently, Wang et al. [32] further investigated the transient mixed lubrication, friction, flash temperature, and fatigue behaviors under three types of contact trajectories.
Lubricated contact interfaces are crucial to tribodynamics analysis in various driving systems because of their nonlinear and time-varying contact stiffness and damping; furthermore, they are believed to be one of the most important sources of vibration and noise. Consequently, significant effort has been expended to predict the stiffness and damping of lubricated contact interfaces. Qin et al. [33] proposed a stiffness model of the EHL line contact, which is composed of oil film stiffness and contact body stiffness. Similarly, Zhang et al. [34,35] compared several EHL contact stiffness models, and the effects of load, speed, regular waviness roughness, and starved oil on stiffness and damping were investigated. Zhou et al. [36] published an oil film stiffness model both in the tangential and normal directions with steady-state non-Newtonian EHL, and then used this model to analyze the stiffness of a modified spur gear. Subsequently, Zhou and Xiao [37] considered normal and tangential damping for gear drive research. Xiao et al. [38] investigated the synthetic stiffness, scuffing, and meshing efficiency of a modified gear pair under extreme conditions based on non-Newtonian transient thermal EHL. However, the studies mentioned above focused on the EHL line contact. Xiao and www.Springer.com/journal/40544 | Friction Shi [24] presented a stiffness and damping model in point contact transient non-Newtonian thermal EHL for crowned herringbone gears, and the effects of the operating conditions and contact parameters on the stiffness and damping of oil films were revealed. Xiao et al. [39] introduced oil film stiffness and damping into a dynamic study of spur gears and compared the results obtained with those of the traditional dynamic model; it highlights that the oil film stiffness contributed effectively to alleviation, whereas oil film damping benefitted vibration and frictional heat reduction. Cao et al. [40] published a tribo-dynamic model for spiral bevel gears using the Hertz contact stiffness without involving oil film stiffness.
In retrospect, studies investigating transient mixed contact stiffness and damping are scarce, although the asperity contact may significantly affect the transient EHL performance and subsurface elastic deformation, in which case disregarding the roughness is impractical. Therefore, a novel mixed contact stiffness and damping model is proposed herein based on a recently published transient mixed lubrication model [31,32]. Moreover, the lack of accurate predictions of transient mixed contact stiffness and damping hinders the comprehensive understanding of the dynamic response of spiral bevel gears. Hence, systematic research was performed to reveal the effects of roughness, transmission torque, rotational speed, and contact trajectory on the transient mixed contact stiffness and damping.

Transient mixed lubrication theory
In this study, a transient mixed lubrication model [31,32] considering the time-varying angle of entraining velocity, as shown in Fig. 1, was used to numerically simulate the transient lubrication characteristics of spiral bevel gears. This transient mixed lubrication model should accurately accomplish two objectives: one is the solution of the basic Reynolds equations to obtain the film thickness and pressure distribution, and the other is the ongoing transmission of the transient squeezing effect between adjacent time steps. For the solutions of film thickness and pressure distribution in the elliptical contact area of spiral bevel gears, a Reynolds equation involving both wedge terms in the x-and y-axis directions and a squeeze term was utilized, expressed as follows: where the right side of the equation includes both the wedge and squeeze terms, in which ue(t) and θe(t) are the variant entraining speed and entraining speed angle, respectively. The film thickness h(x, y, t) is related to the contact position and time, and it can be obtained as the sum of the time-dependent contact geometry, surface elastic deformation, and roughness amplitude, expressed as follows: where h0(t) indicates the transient normal approach of contact surfaces before elastic deformation; x 2 /2Rx(t) and y 2 /2Ry(t) represent the original contact geometries; 1(x, y, t) and 2(x, y, t) are the roughness amplitudes of the real machined surface. In addition, the Boussinesq integral method is employed to calculate elastic deformation v(x, y, t), which is expressed as 548 Friction 10(4): 545-559 (2022) The typically used Barus viscosity-pressure equation is employed to obtain the lubricant viscosity under different pressures, and it is expressed as The Dowson-Higginson density-pressure equation is utilized to determine the effect of pressure on density, and it is expressed as 9 0 9 0.6 10 1 1 1.7 10 The applied load, which should be equal to the pressure over the entire solution domain Ω, is expressed as For the ongoing transmission of the transient squeezing effect between adjacent time steps, the first time step was regarded as a quasi-steady-state as no transient squeezing effect caused by the dynamic alteration in contact parameters had occurred prior to this step. Furthermore, the transient squeezing effect is related to the results of the steady-state solution at the first transient step or the results of the last transient step at the subsequent transient steps. Moreover, to accelerate the numerical calculation and improve the convergence accuracy, the transient PMD method [31] was adopted. Details for solving this transient mixed lubrication model are available in the authors' previous publications [31,32].

Transient mixed contact stiffness and damping model
In engineering practice, a lubricant is introduced to effectively separate the contact surfaces of mechanical transmission components. However, owing to the irregular asperities of surfaces, the direct contact of asperity is almost inevitable, namely mixed lubrication, particularly in situations involving low speeds and heavy loads. In a mixed lubrication area, the high-pressure peak attached to the direct contact of asperities may result in significant elastic deformation in the contact surfaces, thereby significantly affecting the stiffness and damping of the entire mixed lubrication contact area. In addition, the dynamic variations of the contact force, contact geometry, entraining, and sliding velocity may result in a significant oil film squeezing effect, which not only can change the shape of the oil film, but also effect the contact condition of asperities, and hence fluctuate the stiffness and damping of the mixed lubrication contact. Consequently, we propose a contact stiffness and damping model of transient mixed lubrication, in which the stiffness and damping of the oil film as well as the rough subsurface of the solid material are considered as one; this stiffness and damping model is expected to predict any transient mixed lubrication contact regardless of the severity of the asperity contacts.

Stiffness in transient mixed lubrication
The relationship between load and displacement in a non-conforming EHL contact can be characterized using two models [41], i.e., the nonlinear force model (Eq. (7)), and linear force model (Eq. (8)), as expressed below: where F is the applied load; f indicates the EHL lift force; k is the stiffness of the EHL contact; n is the exponent of the nonlinear force model; ho is the normal approach of the contact surfaces before elastic deformation, i.e., h o is a convergent regulator in the solution of the unified Reynolds equation system when using the Newton-Raphson iterative method, and it is obtained using Eq. (2) until the convergence accuracy is satisfied.
For the transient mixed lubrication of spiral bevel gears, to reveal the relation between the applied load F and normal approach displacement h o , a typical transient mixed lubrication case depicted by the load-displacement relationship is shown in Fig. 2(a). As shown, owing to the dynamic load variation from the engaging-in to www.Springer.com/journal/40544 | Friction the engaging-out point, the transient normal approach displacement exhibited an approximately linear response, and this linear relationship became more prominent as the limits to two extremely close transients were reached. Accordingly, we adopted the linear force model in this study to evaluate the transient mixed contact stiffness of spiral bevel gears, and the solution principle is illustrated in Fig. 2.
In transient mixed lubrication, the lubricant and solid material of the contact subsurface can effectively support the applied load through compression deformation in a mixed lubrication area, and the entire contact can be equivalent to a massless spring. The normal approach increments h0(t) at adjacent time steps (as described in Fig. 2(b)), which involve the direct contact of asperity and the squeezing effect of the lubricant, can be obtained from the numerical simulation results of transient mixed lubrication. After utilizing the linear force model with the corresponding applied load increments F(t), the transient mixed contact stiffness as a function of time is expressed as

Damping in transient mixed lubrication
The dynamic meshing motion of a gear pair may result in the transient compression process of mixed lubrication contact, and the compression speed of the mixed lubrication contact, also known as the normal approach speed, can be calculated using the amount of compression (normal approach displacement) and compression time. Moreover, assuming that the entire mixed lubrication contact is a massless system with both stiffness and damping, the dynamic equation in a transient mixed lubrication contact is expressed as where 0 ( ) h t  is the normal approach speed defined as For each transient t, after the stiffness is obtained using Eq. (9), the damping of the transient mixed lubrication can be calculated by employing the equation introduced above.

Results and discussion
Herein, a contact stiffness and damping model for the transient mixed lubrication of spiral bevel gears is presented, and the detailed calculation process is shown in Fig. 3. This investigation not only highlights the effect of mixed lubrication on the stiffness and damping of spiral bevel gears, but also reveals the effects of the operating conditions. Hence, this study can provide strong evidence for the optimization of the design of vibration and impact reduction for the reasonable selection of operating conditions as well as a more accurate stiffness and damping model for the kinematic investigation of spiral bevel gears.

Comparison with Hertz contact
Owing to the lack of a comprehensive analysis regarding the transient mixed lubrication of spiral bevel gears, the stiffness and damping are typically expressed by the simplified Hertz contact, which may result in considerable errors because the effects of direct contact of asperities, squeezing effect of lubricant, and the interaction between them are not considered. Hence, a comparison of the stiffness and damping between the transient mixed lubrication and Hertz contact was performed. Subsequently, a systematic investigation was performed to determine the effects of real machined roughness and various operating conditions, such as transmission torque, rotational speed, and contact trajectories.
To numerically investigate the stiffness and damping of the transient mixed lubrication of spiral bevel gears, the transient mixed lubrication model [32] was used. The space contact geometry of the spiral bevel gears, as shown in Fig. 4, and the loaded tooth contact analysis (LTCA) method [28] were employed first to obtain the contact parameters, including the contact load, contact geometry, and speed along the contact trajectory of the curved surface. Subsequently, the contact parameters were imported into the transient mixed lubrication model. The geometric parameters are listed in Table 1, whereas the machine parameters of the spiral bevel gears are described in Ref. [28]. The torque applied to the pinion was 100 N·m, and the rotational speed of the gear was 500 rpm. In addition, a conversion strategy that changed the rotating angle into the pertinent time was processed to facilitate the solution of the transient EHL. The time-varying contact parameters comprised the load, equivalent contact radii of the curvatures, and entraining speed in both the x-and y-axes, as shown in Fig. 5.
Variations in the normal displacement at each transient step are required to obtain the evolution of the mixed contact stiffness and damping along the contact path. Therefore, the transient mixed lubrication model of our previous study was modified to further include the output of the normal approach of contact bodies at each transient step. In this transient mixed lubrication model, a real 3D machined roughness surface was utilized to simulate the real surfaces of both the pinion and gear, as shown in Fig. 6, and both roughness values were adjusted to 0.707 μm, whereas the root mean square roughness was 1.0 μm. In addition, the numerical solution of the transient EHL was    Table 2.
Based on the results of the transient mixed lubrication simulation, the stiffness and damping model presented above was employed to compute where c ( ) v t is the normal approach of the Hertz contact at time t and is equivalent to the elastic deformation of the Hertz contact center，and c ( ) v t  is the compression speed defined as c . The comparisons of stiffness and damping between the transient mixed lubrication and Hertz contact, as well as the percentage of Hertz contact K (black line) to transient mixed lubrication contact R (red line) are plotted in Fig. 7. As shown in Fig. 7, except for the early period of stable engagement, the contact stiffness and damping of the Hertz contact were higher than those of the transient mixed lubrication up to 60% and 120%, respectively, which was primarily due to the considerable variation in film thickness resulting from the dynamic change of contact parameters and the squeezing action. In addition, the stiffness and damping of both the transient mixed lubrication and Hertz contact were subject to some fluctuations, particularly in the engaging-in and engaging-out periods; this is believed to be due to the uneven increase in load in equal time intervals, as shown in Fig. 8. Nevertheless, the complex contact geometry of the spiral bevel gears may ultimately account for this   phenomenon. In other words, the applied load increment is critical to the contact stiffness and damping. In the transient mixed lubrication of spiral bevel gears, the mixed contact stiffness in the period of stable engagement was much lower than that during engaging-in and engaging-out, owing to the slow increase in the applied load, as depicted in Fig. 5. However, the mixed damping during the period of stable engagement was higher owing to the compression speed or normal approach speed; however, it decelerated under a relatively stable high load. In addition, the applied load was affected by abrupt changes at the instantaneous alternate points of the gear pair engagement; consequently, the compression of the mixed lubrication contact accelerated. Therefore, abrupt changes in stiffness and damping appeared at these alternate points, as shown in Fig. 7.

Effect of roughness
The asperity contact of the transient mixed lubrication may cause a local high-pressure, which may result in considerable surface elastic deformation. In addition, the change in the rough surface morphology affects the transient squeezing effect; hence, it is important to study the effect of roughness on the transient mixed contact stiffness and damping. In this study, a real 3D machined roughness surface measured using an optical profilometer was used to simulate the transient mixed lubrication, as shown in Fig. 6. This adjustable roughness data matrix was set to RMS 0 (smooth), 1.0, and 2.0 μm. The obtained results of transient mixed stiffness and damping of spiral bevel gears under RMS 0.0 (K), 1.0 (B), and 2.0 μm (R) as well as the percentage of RMS 0.0 to RMS 2.0 μm are shown in Fig. 9. As shown, roughness significantly affected the mixed contact stiffness and damping. For the engaging-in period of the spiral bevel gears, the mixed stiffness and damping increased with the roughness; however, a complex tendency emerged in the period of stable engagement. The increase in roughness provided higher stiffness and damping in the early stage. Contrarily, it provided lower stiffness and damping in the later stage, a phenomenon attributable to the complex interaction between roughness and the propagation of the transient squeezing effect along the direction of the lubricant flow. For the same reason, roughness may aggravate the fluctuation of stiffness and damping, as depicted in Fig. 9. However, roughness appeared to impose a limited effect on stiffness and damping during the engaging-out period. To reveal the effect www.Springer.com/journal/40544 | Friction of roughness on the mixed contact stiffness and damping, Fig. 10 shows the transient mixed lubrication normal approach of spiral bevel gears under three RMSs. Figure 11 shows the transient mixed lubrication asperity contact area and load ratio of spiral bevel gears under different RMSs, whereas Fig. 12 depicts the transient mixed lubrication oil film thickness contours of three RMSs at the most severe moment (tM) of the rough peak contact. It is noteworthy that the yellow areas surrounded by red imply the direct contact of the asperity. As shown in Fig. 10, with the increase in the RMS, the normal approach reduced over the entire engagement period owing to the increase in roughness, resulting in a more severe direct contact of the asperity, as shown in Figs. 11(a) and 12. Furthermore, more direct contacts of the solid material can effectively resist the normal approach by bearing a higher load with asperities, as shown in Fig. 11(b).

Effect of transmission torque
The applied load is critical to the stiffness and damping. In this study, three torques were applied to the pinion, i.e., 90, 100, and 110 N·m, as shown in Fig. 13    period with the increase in the transmission torque, particularly for the engaging-in and engaging-out periods, owing to the greater applied torque, which resulted in a greater load increment. With the increase in the applied torques, the mixed damping increased, particularly in the period of stable engagement, owing to the contact pressure. Furthermore, the squeeze effect of the lubricant was the most prominent, and the resistance to compression motion strengthened.

Effect of rotational speed
Three rotational gear speeds, i.e., 500, 1,000, and 1,500 rpm, were selected to reveal the effect of rotational speed on the transient mixed contact stiffness and damping under RMS 0.5 μm. The obtained transient mixed contact stiffness and damping of the spiral bevel gears at 500 (K), 1,000 (B), and 1,500 rpm (R) as well as the percentage of 500 to 1,500 rpm are shown in Fig. 15. As shown, the rotational speed had a limited effect on the transient mixed contact stiffness of the spiral bevel gears. However, the rotational speed was critical to the transient mixed contact damping. Furthermore, the greater the rotational speed, the smaller the mixed damping, owing to the compression speeds of the transient mixed lubrication increasing with the rotational speed [24].

Effects of different contact trajectories
The contact trajectory significantly affected the mixed lubrication action, friction, and contact fatigue [28,32], signifying its great potential effect on the tribo-dynamics response of spiral bevel gears. Therefore, the transient mixed stiffness and  damping under different contact trajectories should be investigated.
In this study, three contact trajectories, i.e., the heel, middle, and toe contacts, were adopted as references [28,42] to illustrate this effect. The meshing center points of the heel, middle, and toe contacts (represented as 1, 0, and -1, respectively) along the pitch cone were located at 1/4, 2/4, and 3/4 of the width of the tooth flank Bw, respectively, as shown in Fig. 16; additionally, the three contact patterns along the action line are shown in this figure. The adjustment parameters of the three contact patterns are listed in Table 3, and the adjustment directions are shown in Fig. 4. In all of the three contact patterns, the torque applied to the pinion was 100 N·m, the rotational speed of the gear was 500 rpm, and 200 data points were extracted from each pattern. The dynamic variations of the load, equivalent entraining velocity, and radii of curvature along the x-and y-axes along the action track of the spiral bevel gears resulting from the LTCA are shown in Fig. 17.
After obtaining the transient mixed lubrication results of three contact patterns under an RMS of 0.5 μm, the contact stiffness and damping of the spiral bevel gears can be obtained immediately using the transient mixed contact stiffness and damping model, as shown in Fig. 18. As shown, the transient mixed contact stiffness and damping of the three contact trajectories were closely associated with the applied loads, as shown in Fig.  17. The toe contact provided the highest and least enduring stiffness and damping during the period of stable engagement, whereas the lowest and most enduring was provided by the heel contact. The engaging-in and engaging-out periods were the longest for the toe contact, whereas they were the shortest for the heel contact. Benefiting from the reduction in impact and vibration for high stiffness and damping, the middle contact is more efficient for maintaining the stability of spiral bevel gear transmission owing to its comprehensive high stiffness and high damping.

Conclusions
Herein, a novel contact stiffness and damping model of spiral bevel gears in transient mixed lubrication is proposed. The model considers the transient mixed lubrication contact as one. The linear force model was adopted after argumentation, and it is expected to precisely predict the timevarying contact stiffness and damping regardless of the severity of the asperity contacts as well as smooth the surface contact in a unified stiffness and damping model. Additionally, to determine factors that affect the transient mixed contact stiffness and damping, the roughness, applied torque, rotational speed, and contact trajectory were further investigated. The results obtained for the transient mixed contact stiffness and damping of the spiral bevel gears were smaller than those of the conventional Hertz contact, except at the early period of stable engagement. The applied load increment was critical to the transient mixed stiffness and damping, and the uneven increase in load in equal time intervals might cause fluctuations in the mixed stiffness and damping. Roughness significantly affected the mixed contact stiffness and damping; consequently, the fluctuation in stiffness and damping may be aggravated owing to the mutual effects of www.Springer.com/journal/40544 | Friction roughness and the propagation of the transient squeezing effect along the direction of the lubricant flow. The mixed contact stiffness and damping increased as the transmission torque increased during the entire engagement period. The rotational speed had a limited effect on the transient mixed stiffness of the spiral bevel gears. However, the rotational speed was critical to the transient mixed contact damping, i.e., the greater the rotational speed, the smaller the mixed damping. As for the effects of different contact trajectories, the toe contact provided the highest and least enduring stiffness and damping in the period of stable engagement, whereas the lowest and most enduring was provided by the heel contact. Meanwhile, the engaging-in and engaging-out periods were the longest for the toe contact, whereas they were the shortest for the heel contact. Consequently, the middle contact was more efficient for maintaining the stability of spiral bevel gear transmission owing to its comprehensive high stiffness and damping.
This transient mixed contact stiffness and damping model is expected to be applied in various transmission components with lubricated contacts. Moreover, it demonstrated potential in reducing vibration and impact, as well as more accurate predictions of the tribo-dynamic response for spiral bevel gears.