Research on the dynamic performance of double-helical planetary gear transmission under friction and wear

Abstract

In this paper, the main research is to investigate the effects of mechanical, friction, temperature, and tooth wear on the dynamic performance of a double-helical planetary gear transmission system. To achieve this, a new lumped parameters model is presented, which combines dynamics and lubrication. The model allows for the study of excitation effects and considers various parameters such as material, lubricating oil, and working conditions. By proposing a multi-source excitation calculation method, the paper provides a novel approach that has not been reported in the existing literature. Additionally, the paper establishes a tooth wear model and a lubrication model for the double-helical gear pair. It shows that under constant torque, the displacement and acceleration amplitudes of the system initially decrease and then increase with increasing input speed. Under constant speed, the displacement and acceleration amplitudes increase with increasing torque. Amplitudes decrease during running-in wear, increase in stable wear, and rapidly increase in severe wear. Boundary lubrication minimizes vibration, while elastohydrodynamic lubrication maximizes it.

Introduction

A double-helical planetary gear transmission (PGT) system has significant advantages, including small axial force, low noise, high load-bearing capacity, and high transmission power, which is widely used in heavy-load and high-power equipment such as marine, aerospace, aviation, and petroleum.^1–3 There is an impact load in the gear transmission process, which is due to the influence of transmission error, time-varying mesh stiffness, tooth surface friction, and tooth wear. As a result, the thickness, pressure, and temperature distribution of oil film in the lubricating contact area are also changed. In return, the oil film shear force and the friction force generated by the contact of the asperity peaks and the friction torque formed by their combined action change the dynamic characteristics.

With the development of mechanics theory and the advancement of experimental measurement techniques, the theoretical foundation of the dynamics of herringbone gear transmission systems has become increasingly refined. Khoozani et al. ⁴ proposed a lumped parameter model for herringbone PGT systems considering the gyroscopic effect, studying the inherent characteristics of these gear systems. Shuai et al. ⁵ established a dynamic model for herringbone PGT systems, investigating the influence of sun gear float characteristics on the system's dynamic behavior. Siyu et al., ⁶ Hu et al., ⁷ and Siyu and Jinyuan ⁸ proposed a dynamic model for herringbone gear transmission systems considering damping, eccentric mass, and time-varying mesh stiffness, calculating the natural frequencies and mode shapes of these systems, and determining the effects of tooth backlash, transmission error, and crossover angle on their dynamic behavior. Liu et al. ⁹ presented a dynamic analysis model for herringbone gear reducers, conducting theoretical and experimental research on the vibration characteristics of these gear systems. Kang and Kahraman^10,11 conducted theoretical and experimental investigations on the dynamic behavior of herringbone gear systems, studying the influence of key design and manufacturing parameters on their dynamic characteristics, concluding that deviations in left and right crossover angles and the axial gear support conditions have a significant impact on system vibration. Paouris et al. ¹² proposed a friction dynamics model for predicting the contact efficiency and dynamic response of hypoid gears in automotive applications. This model couples the dynamic model with an analytical friction model (viscous and boundary). Fatourehchi et al. ¹³ conducted an integrated analysis of friction dynamics and wheel center planetary systems to study the close correlation between transmission efficiency and noise, vibration and harshness (NVH) optimization. It was determined that there may be conflicting requirements between NVH optimization and transmission efficiency. Tatar et al. ¹⁴ developed a three-dimensional (3D) dynamic model for planetary gear rotor systems, studying modal sensitivity analysis of the number of planetary gears, planetary gear meshing errors, planetary gear masses, gear mesh stiffness, and planetary gear speeds. Dai et al. ¹⁵ proposed a new dynamic model for PGT, studying the behavioral variations in PGT caused by gear pin-hole position errors, and obtained the frequency modulation mechanism of gear vibration in PGT. Liu et al. ¹⁶ used angular displacement and translational displacement as generalized coordinates to establish lumped parameter dynamic sub-models for sun gears and planetary gears, studying the coupled dynamic analysis of flexible planetary gears under high-speed and variable speed conditions. Jiang et al. ¹⁷ comprehensively considered the characteristics of permanent magnet synchronous motors, internal gear structures, and load torques, and based on motor theory and Maxwell's equations, established a dynamic model for permanent magnet gear transmission systems. Statistical analysis of gear fatigue load spectra was conducted using the Rainflow counting method and Goodman theory. Nie et al. ¹⁸ proposed a new dynamic model for PGT systems by considering all vibration sources and their corresponding transmission path effects based on geometric parameters and material properties, studying the dynamic response of the system. Luo and Li ¹⁹ established a thermodynamic model for gear transmission systems, introduced gear meshing thermal stiffness, analyzed the system's equivalent viscous damping using the Coulomb friction damping model, and studied the thermal dynamic response of the system. Sanchez-Espiga et al. ²⁰ investigated the influence of meshing phase, tooth thickness error, and planetary carrier pin-hole position error on the load-sharing behavior of PGT systems.

Practices have shown that gear transmission systems with superior lubrication performance can extend the service life of equipment while reducing failures such as pitting, wear, and adhesive wear on the gear surfaces. In recent years, helical gear planetary drive systems have been operating at increasingly faster speeds, transmitting higher power, and under more severe operating conditions. Li and Kahraman ²¹ proposed a model for predicting the load-related power loss of gear pairs based on the elastohydrodynamic lubrication (EHL) model. The EHL model includes all key contact parameters of gear pairs that change with time, such as rolling and sliding velocities, curvature radii, and normal loads, enabling continuous analysis of gear contact from the root to the tip under any lubrication conditions. Simon ²² studied the influence of manufacturing parameters on mixed EHL conditions, including full film lubrication, mixed lubrication, and boundary lubrication. Based on this, suggestions were made for improving the efficiency of mixed EHL and face milling of hyperbolic gears. Zhao et al. ²³ established a coupled analysis model for the mixed lubrication and wear characteristics of involute gears considering the microtopography of the gear surfaces. The influence of lubricant composition on the distribution of lubrication performance parameters and wear characteristics was analyzed. Ouyang et al. ²⁴ proposed a high-speed spur gear lubrication dynamics model based on the frictional dynamics theory. The lubrication and dynamic performance of high-speed spur gears were studied by combining finite-element models, and the time-varying mesh stiffness and static transmission error were derived precisely in the form of a Fourier series.

Through the development of scientific research and the accumulation of production practice, researchers have deepened their understanding of the wear mechanism and proposed empirical formulas for predicting wear volume and physical models for describing wear. Kahraman et al. ²⁵ combined the finite-element gear contact mechanics model, sliding distance calculation algorithm, and Archard wear equation for predicting contact surface wear to propose a surface wear prediction model. The influence of tooth profile deviation on gear surface wear was studied. Shen et al. ²⁶ proposed a pure torsional dynamic model for planetary gears considering gear wear, treating gear wear as a tooth deviation, and analyzing the effect of gear wear on the vibration response of planetary gears. Liu ²⁷ proposed a dynamic analysis model for a grade planetary gear system considering the wear of the meshing gear surface and evaluated the gear wear status using vibration signals in the planetary gear system. Tian et al. ²⁸ calculated the wear depth based on the Archard wear theory, obtained the load distribution through 3D load tooth contact analysis, and considered the influence of angular misalignment by introducing the variation of meshing clearance between two meshing tooth surfaces. The effect of angular deviation on wear prediction of planetary gear sets with rotating carriers was studied. Sánchez et al. ²⁹ proposed an analytical model for calculating the mesh stiffness, load sharing, and transmission error of spur cylindrical gears, considering the effects of surface wear. The variation of transmission error with wear cycle number was investigated, and an optimal gear modification method to minimize dynamic loads caused by transmission errors was proposed. Xiao et al. ³⁰ obtained the basic parameters of the EHL and wear model for herringbone gears through the analysis of meshing characteristics and kinematics. The influence of small gear torque, speed, and surface roughness on maximum pressure, temperature rise, and minimum oil film thickness was discussed under full film lubrication conditions. Meanwhile, the rough load ratio, flash temperature, and tooth surface wear depth during the gear start–stop stage were studied under mixed lubrication point contact.

In the process of gear transmission, it involves temperature field, thermal deformation, friction, tooth wear, lubrication, and vibration. These contents are intrinsically linked and form an organic whole, few articles study it as a whole. This paper takes the high-speed and heavy-duty double-helical planetary gear pair as the research object and comprehensively analyzes the dynamic performance of the gear system considering factors such as thermal deformation, tooth wear, and lubrication. The research content of this paper is as follows. Firstly, a numerical analysis lubrication model for the high-speed and heavy-duty double-helical gear pairs is proposed, considering the surface roughness, the non-Newtonian properties of lubricant oil, the wear of the tooth surface, and the thermal effect. Based on the meshing characteristics of double-helical gears, the basic parameters of the lubrication model are calculated, and the control equations of the non-Newtonian fluid lubrication model considering wear, tooth surface roughness, and thermal effects are established. Secondly, considering the coupling relationship between dynamic and tribological characteristics, a multi-degree-of-freedom gap-nonlinear dynamic model of a double-helical planetary transmission system is put forward by taking into consideration the influence factors, such as time-varying meshing stiffness, meshing impact, tooth surface wear, tooth surface friction, and tooth surface temperature rise effect. Based on the theory of tooth surface wear and point-contact thermal EHL, a numerical iterative solution combination Runge–Kutta method and multi-grid method are used to carry out the coupled dynamics equations. Furthermore, the influence of working conditions, tooth wear, and lubrication on the dynamic characteristics of the gear system is obtained.

Dynamic modeling methodology for double-helical planetary transmission system

Excitation analysis of double-helical planetary transmission system

Figure 1 is a schematic diagram of the double-helical gear planetary transmission system. The sun gear and the input shaft are connected by splines, the power is transmitted from the sun gear to the planetary gear. Then, the power is divided by the planetary gear and output by the carrier. The sun gear and planetary gear are double-helical gears with undercuts, the inner gear is composed of two helical gears with different rotation directions but the same other geometric dimensions. The dotted lines and arrows in the figure represent the power transmission path. Table 1 shows the structural parameters of the double-helical planetary transmission system.

Figure 1.

Diagram of double-helical planetary transmission system: (a) structure diagram and (b) three-dimensional (3D) solid model.

Table 1.

Structural parameters of double-helical planetary transmission system.

Item	Tooth numbers	Normal propagation module (mm)	Speed (r/min)	Power (kW)	Pressure angle (°)	Helix angle (°)
Sun	34	3	5000	2000	20	30
Planet	31
Ring	96

With operation numbers increasing for the gear transmission system, the roughness peak of the contact surface of the gear pair becomes smaller and smoother and even the thickness of the gear becomes thinner, so that the geometric parameters of the contact surface will change. However, the surface roughness of the contact area has a significant effect on the lubrication state. Therefore, the lubrication characteristics model should consider the effect of tooth wear.

Figure 2 is a diagram of tooth wear. Assuming that the surface of the contact area of the gear pair is smooth, l is the half-width of the wear plane, and its calculation formula is:

$$l=a\cdot k^l$$

(1)

where k^l represents the wear coefficient and a represents the contact radius, respectively.

Figure 2.

Schematic diagram of tooth wear.

The wear depth can be obtained from Figure 2.

$$h_w=R_y-\sqrt{R_y^2-l^2}$$

(2)

where R_x and R_y represent the curvature radius of the meshing point in the x and y directions, respectively. The contact surface radius is set to $R_x^{\mathrm{\prime }}$ :

$$R_x^{\mathrm{\prime }}=R_x-h_w$$

(3)

Figure 3 is a model of adhesive wear. Using the wear model proposed by Ludema, ³¹ the calculation method of single point cumulative wear on the contact surface is as formula (4).

$$\mathrm{\Delta }{h}_{1,2}=2ap{k}^{\iota }\left|{\displaystyle \frac{V_1-V_2}{V_{1,2}}}\right|$$

(4)

Figure 3.

Model of adhesive wear: (a) bond point formation and (b) bond point failure.

If the accumulated wear on the tooth surface is small, the geometric morphology and pressure distribution of the tooth surface can be approximately assumed to remain unchanged. When the number of meshing is sufficient and the accumulated wear reaches the threshold, it is necessary to update the contact surface and pressure distribution. Use the superscript k to indicate the update of the kth contact surface and pressure, and use the symbol c^k to indicate the number of cycles of the driving wheel during the kth update step (i.e. from the kth contact surface and pressure distribution update to the next update). The accumulated wear amount of point A on the tooth profile during the kth update step can be expressed as

$$\mathrm{\Delta }{h}_{1,2}^k=2a^k{p}^k{k}^{\iota }\left|{\displaystyle \frac{V_1^k-V_2^k}{V_{1,2}^k}}\right|$$

(5)

The cumulative wear amount of point A in the kth update step is

$$\mathrm{\Delta }{h}_{1,2}^{{\xi }^k}=c^k\mathrm{\Delta }{h}_{1,2}^k$$

(6)

Assuming that after k updates, the accumulated wear on the tooth surface reaches the allowable value max. Then the total cumulative wear and cycle times at point A are

$$h_{1,2}=\sum _{k=0}^K\mathrm{\Delta }{h}_{1,2}^{{\xi }^k}$$

(7)

$$c^{\mathrm{tot}}=\sum _{k=0}^K{c}^k$$

(8)

The increased temperature of the tooth surface Δf, which will cause the thermal expansion of the gear to deviate from the tooth profile involution. The degree of deviation can be expressed by the involution error. As shown in Figure 4b, the dotted line is the tooth profile shape before the gear is heated and deformed, and the solid line is the tooth profile shape after deformation, that is, the position where any point K on the tooth profile changes to K′. At this time, the tooth profile error caused by the temperature rise of the tooth surface can be expressed as ³² :

$$G_{{\mathrm{T}}_i}(t)={\displaystyle \frac{-\mathrm{\Delta }f(t)\lambda r_{\mathrm{b}i}(r_{\mathrm{b}i}+u_{\mathrm{b}i})}{2\times [r_{\mathrm{b}i}+u_{\mathrm{b}i}\cos {\alpha }_{ki}+\mathrm{\Delta }f(t)\lambda r_{\mathrm{b}i}(1-\cos {\alpha }_{ki}]}\times \left[{\displaystyle \frac{s_{\mathrm{g}i}}{r_{\mathrm{b}i}}-2(inv{\alpha }_{ki}-inv\alpha )}\right]}$$

(9)

$${\alpha }_{ki}=\arccos {\displaystyle \frac{r_{\mathrm{b}i}+u_{\mathrm{b}i}}{r_{\mathrm{c}i}(t)}}$$

(10)

where the negative sign represents the theoretical involution profile within the actual gear profile, $\mathrm{\Delta }f$ represents the temperature rise of the gear, s_g represents tooth thickness, $r_{\mathrm{b}i}$ (i = 1, 2), respectively, represents the base circle radius of the two gears, $r_{\mathrm{c}i}$ (i = 1, 2), respectively, represent the curvature radius of the two gears. λ represents the linear expansion coefficient of the gear material, ${\alpha }_{ki}$ represents the tooth profile pressure angle of the contact point after thermal expansion. $u_{\mathrm{b}i}$ (i = 1, 2), respectively, represent the thermal expansion deformation of the two gears. ³²

$$u_{\mathrm{b}i}=\lambda r_{\mathrm{b}i}T{r}_{\mathrm{o}i}+{\displaystyle \frac{1+\kappa }{1-\kappa }\cdot {\displaystyle \frac{\lambda r_{\mathrm{b}i}[r_{\mathrm{b}i}^2(1-2\kappa )-r_{\mathrm{o}i}^2]}{r_{\mathrm{b}i}^2-r_{\mathrm{o}i}^2}\cdot (T{r}_{\mathrm{b}i}-T{r}_{\mathrm{o}i})}}$$

(11)

where $T{r}_{\mathrm{o}i}$ represents the temperature of the driving shaft and the driven shaft when they are running stably, $T{r}_{\mathrm{b}i}$ represents the temperature on the base circle of the two gears, $r_{\mathrm{o}i}$ represents the shaft diameter where the driving shaft and the driven shaft fit with the gear.

Figure 4.

Model of tooth profile error: (a) tooth surface wear and tooth profile deformation and (b) thermal expansion of the tooth profile.

Therefore, the tooth profile deviation considering wear and temperature rise is:

$$\begin{array}{r}(G_{\mathrm{T}+\mathrm{\Delta }\mathrm{T},h}^K{)}^{1,2}=(G_{\mathrm{T}}^0{)}^{1,2}-\sum _{i=1}^K{(h_i^k)}^{1,2}+(G_{\mathrm{T}+\mathrm{\Delta }\mathrm{T}}{)}^{1,2}\end{array}$$

(12)

where $(G_{\mathrm{T}}^0{)}^{1,2}$ represents the initial tooth profile deviation of driving gear 1 or driven gear 2 caused by manufacturing errors; ${\sum }_{i=1}^K{(h_i^k)}^{1,2}$ represents the wear level of driving gear 1 or driven gear 2 after k cycles. $(G_{\mathrm{T}+\mathrm{\Delta }\mathrm{T}}{)}^{1,2}$ represents the tooth profile deviation of the driving gear 1 or the driven gear 2 due to the temperature rise of the tooth surface.

Generalized displacement of meshing pair

According to the relative position relationship of the gear pairs, the elastic deformation coordination conditions of the inner and outer meshing pairs for the double-helical PGT system were established (Figures 5 and 6). The generalized relative displacement expression for the meshing pairs can be derived.

$${\psi }_{\mathrm{sp}i}=\{\begin{array}{cc}{\varphi }_{\mathrm{sp}}-{\alpha }_{\mathrm{p}i}& T_s:\mathrm{Counterclockwise}\\ {\varphi }_{\mathrm{sp}}-{\alpha }_{\mathrm{p}i}& T_s:\mathrm{Clockwise}\end{array}$$

(13)

where T_s is the external torque acting on the sun gear, α_pi represents one side (either left or right side) of the sun gear (subscript s) meshing with the same side of planet-i (subscript pi) located at an arbitrary position, and ϕ_sp is the transverse pressure angle, respectively.

Figure 5.

Relative displacement of inner and outer meshing gear pairs.

Figure 6.

Relative deformation of planetary gear and carrier.

Relative displacement on the left sun–planet gear meshing line is as follows:

$$\begin{array}{rl}{\delta }_{\mathrm{spiL}}& =[(x_{\mathrm{sL}}-x_{\mathrm{piL}})\sin {\psi }_{\mathrm{spiL}}+(y_{\mathrm{sL}}-y_{\mathrm{piL}})\cos {\psi }_{\mathrm{spiL}}+({\theta }_{\mathrm{szL}}{r}_{\mathrm{bsL}}+{\theta }_{\mathrm{pizL}}{r}_{\mathrm{bpiL}})]\cos \beta \\ & +[(r_{\mathrm{bsL}}{\theta }_{\mathrm{sxL}}+r_{\mathrm{bpiL}}{\theta }_{\mathrm{pixL}})\sin {\psi }_{\mathrm{spiL}}+(r_{\mathrm{bsL}}{\theta }_{\mathrm{syL}}+r_{\mathrm{bpiL}}{\theta }_{\mathrm{piyL}})\cos {\psi }_{\mathrm{spiL}}+(z_{\mathrm{sL}}-z_{\mathrm{piL}})]\sin \beta -e_{\mathrm{spiL}}\end{array}$$

(14)

Relative displacement on the right sun–planet gear meshing line is as follows:

$$\begin{array}{rl}{\delta }_{\mathrm{spiR}}& =[(x_{\mathrm{sR}}-x_{\mathrm{piR}})\sin {\psi }_{\mathrm{spiR}}+(y_{\mathrm{sR}}-y_{\mathrm{piR}})\cos {\psi }_{\mathrm{spiR}}+({\theta }_{\mathrm{szR}}{r}_{\mathrm{bsR}}+{\theta }_{\mathrm{pizR}}{r}_{\mathrm{bpiR}})]\cos \beta \\ & +[(r_{\mathrm{bsR}}{\theta }_{\mathrm{sxR}}+r_{\mathrm{bpiR}}{\theta }_{\mathrm{pixR}})\sin {\psi }_{\mathrm{spiR}}+(r_{\mathrm{bsR}}{\theta }_{\mathrm{syR}}+r_{\mathrm{bpiR}}{\theta }_{\mathrm{piyR}})\cos {\psi }_{\mathrm{spiR}}-(z_{\mathrm{sR}}-z_{\mathrm{piR}})]\sin \beta -e_{\mathrm{spiR}}\end{array}$$

(15)

Relative displacement on the left ring-planet gear meshing line is as follows:

$$\begin{array}{rl}{\delta }_{\mathrm{pirL}}& =[(x_{\mathrm{piL}}-x_{\mathrm{rL}})\sin {\psi }_{\mathrm{pirL}}+(y_{\mathrm{rL}}-y_{\mathrm{piL}})\cos {\psi }_{\mathrm{pirL}}+({\theta }_{\mathrm{rzL}}{r}_{\mathrm{brL}}-{\theta }_{\mathrm{pizL}}{r}_{\mathrm{bpiL}})]\cos \beta \\ & +[(r_{\mathrm{bpiL}}{\theta }_{\mathrm{pixL}}+r_{\mathrm{brL}}{\theta }_{\mathrm{rxL}})\sin {\psi }_{\mathrm{pirL}}+(r_{\mathrm{bpiL}}{\theta }_{\mathrm{piyL}}+r_{\mathrm{brL}}{\theta }_{\mathrm{ryL}})\sin {\psi }_{\mathrm{pirL}}-(z_{\mathrm{pirL}}-z_{\mathrm{rL}})]\sin \beta -e_{\mathrm{pirL}}\end{array}$$

(16)

Relative displacement on the right ring-planet gear meshing line is as follows:

$$\begin{array}{rl}{\delta }_{\mathrm{pirR}}& =[(x_{\mathrm{piR}}-x_{\mathrm{rR}})\sin {\psi }_{\mathrm{pirR}}+(y_{\mathrm{rR}}-y_{\mathrm{piR}})\cos {\psi }_{\mathrm{pirR}}+({\theta }_{\mathrm{rzR}}{r}_{\mathrm{brR}}-{\theta }_{\mathrm{pizR}}{r}_{\mathrm{bpiR}})]\cos \beta \\ & +[(r_{\mathrm{bpiR}}{\theta }_{\mathrm{pixR}}+r_{\mathrm{brR}}{\theta }_{\mathrm{rxR}})\sin {\psi }_{\mathrm{pirR}}+(r_{\mathrm{bpiR}}{\theta }_{\mathrm{piyR}}+r_{\mathrm{brR}}{\theta }_{\mathrm{ryR}})\sin {\psi }_{\mathrm{pirR}}+(z_{\mathrm{pirR}}-z_{\mathrm{rR}})]\sin \beta -e_{\mathrm{pirR}}\end{array}$$

(17)

Relative deformation and displacement of planet carrier and planet gear are as follows:

$$\{\begin{array}{c}{\delta }_{\mathrm{xpic}}=(x_{\mathrm{c}}-x_{\mathrm{p}i})-r_{\mathrm{c}}{\theta }_{\mathrm{zc}}\sin {\alpha }_{\mathrm{p}i}\\ {\delta }_{\mathrm{ypic}}=(y_{\mathrm{c}}-y_{\mathrm{p}i})+r_{\mathrm{c}}{\theta }_{\mathrm{zc}}\sin {\alpha }_{\mathrm{p}i}\\ {\delta }_{\mathrm{zpic}}=(z_{\mathrm{c}}-z_{\mathrm{p}i})-r_{\mathrm{c}}{\theta }_{\mathrm{yc}}\cos {\alpha }_{\mathrm{p}i}+r_{\mathrm{c}}{\theta }_{\mathrm{xc}}\cos {\alpha }_{\mathrm{p}i}\end{array}$$

(18)

where $r_{\mathrm{bsL}}$ , $r_{\mathrm{bsR}}$ , $r_{\mathrm{bpiL}}$ , $r_{\mathrm{bpiR}}$ , $r_{\mathrm{brL}}$ , and $r_{\mathrm{brR}}$ represent the base circle radii on of the sun gear, planetary gear, and inner gear, respectively, $e_{\mathrm{spiL}}$ , $e_{\mathrm{spiR}}$ , $e_{\mathrm{pirL}}$ , and $e_{\mathrm{pirR}}$ represent the comprehensive equivalent error on both sides of the inner and outer meshing gear teeth, respectively, $\psi$ is the angle between the meshing plane and the y-axis.

Meshing force and frictional force of gear pair

The meshing force and friction force of the meshing pair are as follows:

$$\{\begin{array}{c}{F}_{\mathrm{spiL}}=P_{\mathrm{spiL}}+D_{\mathrm{spiL}}\\ F_{\mathrm{spiR}}=P_{\mathrm{spiR}}+D_{\mathrm{spiR}}\\ F_{\mathrm{pirR}}=P_{\mathrm{pirR}}+D_{\mathrm{pirR}}\\ F_{\mathrm{pirL}}=P_{\mathrm{pirL}}+D_{\mathrm{pirL}}\end{array}\{\begin{array}{c}{F}_{\mathrm{fspiL}}=u_{\mathrm{spi}}{F}_{\mathrm{spiL}}\\ F_{\mathrm{fspiR}}=u_{\mathrm{spi}}{F}_{\mathrm{spiR}}\\ F_{\mathrm{fpirR}}=u_{\mathrm{pir}}{F}_{\mathrm{pirR}}\\ F_{\mathrm{fpirL}}=u_{\mathrm{pir}}{F}_{\mathrm{pirL}}\end{array}\{\begin{array}{c}{P}_{\mathrm{spiL}}=k_{\mathrm{spiL}}f(\delta {)}_{\mathrm{spiL}}\\ P_{\mathrm{spiR}}=k_{\mathrm{spiR}}f(\delta {)}_{\mathrm{spiR}}\\ D_{\mathrm{spiL}}=c_{\mathrm{spiL}}{\stackrel{\cdot }{\delta }}_{\mathrm{spiL}}\\ D_{\mathrm{spiR}}=c_{\mathrm{spiR}}{\stackrel{\cdot }{\delta }}_{\mathrm{spiR}}\end{array}$$

(19)

where u is the friction coefficient for the meshing pair, k_spiL and k_spiR represent the left and right external meshing stiffness, respectively, c_spiL and c_spiR represent the left and right external meshing damping, respectively.

Differential equations of motion for each subsystem

When deriving the dynamic differential equation of the double-helical planetary transmission system, the following assumptions are made:

The contact of gear teeth and bearing support are connected by spring and damping.

The average meshing stiffness, supporting stiffness, mass, moment of inertia, and structural size of the planetary gears distributed evenly along the circumferential direction are equal.

The left and right sides of the double-helical teeth have the same geometric shape except for the opposite direction of rotation.

Figure 7 shows the dynamic model of the double-helical PGT system, with six degrees of freedom for each helical gear to translate along the x, y, and z directions and rotate about the x, y, and z axes.

Figure 7.

Dynamic model of double-helical planetary gear transmission (PGT) system.

Firstly, model derivation of three basic subsystems: (a) sun–planet i pair (left or right), (b) ring gear–planet i pair (left or right), and (c) planet carrier–planet i pair (left or right) or right). Then, the helical gears on the left and right sides are connected with the beam element formula, and finally, the dynamic model of the system is obtained. In this paper, only the dynamic equation of the left meshing pair of the herringbone gear is derived.

(a) Sun–planet meshing pair. The differential equation of motion of the left sun gear is as follows (Figure 8):

$$\{\begin{array}{c}{m}_{\mathrm{sL}}\left({\stackrel{\cdot \cdot }{x}}_{\mathrm{sL}}\text{-}2{\omega }_{\mathrm{c}}{\stackrel{\cdot }{\mathrm{y}}}_{\mathrm{sL}}\text{-}{\omega }_{\mathrm{c}}^2{\mathrm{x}}_{\mathrm{sL}}\right)+c_{\mathrm{sxL}}{\stackrel{\cdot }{\mathrm{x}}}_{\mathrm{sL}}+k_{\mathrm{sxL}}{f}_{\mathrm{sxL}}(x_{\mathrm{sL}})=-\sum _{i=1}^N{F}_{\mathrm{spiL}}\sin {\psi }_{\mathrm{spiL}}\cos \beta +\sum _{i=1}^N{F}_{\mathrm{fspiL}}\cos {\psi }_{\mathrm{spiL}}\\ m_{\mathrm{sL}}\left({\stackrel{\cdot \cdot }{y}}_{\mathrm{sL}}+2{\omega }_{\mathrm{c}}{\stackrel{\cdot }{x}}_{\mathrm{sL}}\text{-}{\omega }_{\mathrm{c}}^2{y}_{\mathrm{sL}}\right)+c_{\mathrm{syL}}{\stackrel{\cdot }{y}}_{\mathrm{sL}}+k_{\mathrm{syL}}{f}_{\mathrm{syL}}(y_{\mathrm{sL}})=-\sum _{i=1}^N{F}_{\mathrm{spiL}}\cos {\psi }_{\mathrm{spiL}}\cos \beta +\sum _{i=1}^N{F}_{\mathrm{fspiL}}\sin {\psi }_{\mathrm{spiL}}\\ m_{\mathrm{sL}}{\stackrel{\cdot \cdot }{z}}_{\mathrm{sL}}+c_{\mathrm{sz}}{\stackrel{\cdot }{z}}_{\mathrm{sL}}+k_{\mathrm{sz}}f(z_{\mathrm{sL}})=\sum _{i=1}^N{F}_{\mathrm{spiL}}\sin \beta \\ I_{\mathrm{sxL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{sxL}}+c_{\mathrm{sx}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{sxL}}+k_{\mathrm{sx}\theta \mathrm{L}}{\theta }_{\mathrm{sxL}}=-\sum _{i=1}^N{F}_{\mathrm{spiL}}\sin \beta \sin {\psi }_{\mathrm{spiL}}{r}_{\mathrm{bsL}}\\ I_{\mathrm{syL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{syL}}+c_{\mathrm{sy}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{syL}}+k_{\mathrm{sy}\theta \mathrm{L}}{\theta }_{\mathrm{syL}}=-\sum _{i=1}^N{F}_{\mathrm{spiL}}\sin \beta \cos {\psi }_{\mathrm{spiL}}{r}_{\mathrm{bsL}}\\ I_{\mathrm{szL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{szL}}=-\sum _{i=1}^N{F}_{\mathrm{spiL}}{r}_{\mathrm{bsL}}\cos \beta +T_s/2+\sum _{i=1}^N{F}_{\mathrm{fspiL}}{H}_{\mathrm{s}}\end{array}$$

(20)

Figure 8.

Dynamic model of sun gear–planet gear meshing pair.

The differential equation of motion of the left planet gear is as follows (Figure 9):

$$\{\begin{array}{c}{m}_{\mathrm{piL}}\left({\stackrel{\cdot \cdot }{x}}_{\mathrm{piL}}\text{-}2{\omega }_{\mathrm{c}}{\stackrel{\cdot }{\mathrm{y}}}_{\mathrm{piL}}\text{-}{\omega }_{\mathrm{c}}^2{\mathrm{x}}_{\mathrm{piL}}\right)+c_{\mathrm{pixL}}\stackrel{\cdot }{{\mathrm{x}}_{\mathrm{piL}}}+k_{\mathrm{pixL}}{f}_{\mathrm{pixL}}(x_{\mathrm{piL}})=F_{\mathrm{spiL}}\sin {\psi }_{\mathrm{spiL}}\cos \beta -F_{\mathrm{fspiL}}\cos {\psi }_{\mathrm{spiL}}\\ m_{\mathrm{piL}}\left({\stackrel{\cdot \cdot }{y}}_{\mathrm{piL}}+2{\omega }_{\mathrm{c}}{\stackrel{\cdot }{x}}_{\mathrm{piL}}\text{-}{\omega }_{\mathrm{c}}^2{y}_{\mathrm{piL}}\right)+c_{\mathrm{piyL}}{\stackrel{\cdot }{y}}_{\mathrm{piL}}+k_{\mathrm{piyL}}{f}_{\mathrm{piyL}}(y_{\mathrm{piL}})=F_{\mathrm{spiL}}\cos {\psi }_{\mathrm{spiL}}\cos \beta -F_{\mathrm{fspiL}}\sin {\psi }_{\mathrm{spiL}}\\ m_{\mathrm{piL}}{\stackrel{\cdot \cdot }{z}}_{\mathrm{piL}}+c_{\mathrm{pizL}}{\stackrel{\cdot }{z}}_{\mathrm{piL}}+k_{\mathrm{pizL}}f(z_{\mathrm{piL}})=-F_{\mathrm{spiL}}\sin \beta \\ I_{\mathrm{pixL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pixL}}+c_{\mathrm{pix}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{pixL}}+k_{\mathrm{pix}\theta \mathrm{L}}{\theta }_{\mathrm{pixL}}=-F_{\mathrm{spiL}}\sin \beta \sin {\psi }_{\mathrm{spiL}}{r}_{\mathrm{bpL}}\\ I_{\mathrm{piyL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{piyL}}+c_{\mathrm{piy}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{piyL}}+k_{\mathrm{piy}\theta \mathrm{L}}{\theta }_{\mathrm{piyL}}=-F_{\mathrm{spiL}}\sin \beta \cos {\psi }_{\mathrm{spiL}}{r}_{\mathrm{bpL}}\\ I_{\mathrm{pizL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pizL}}=-F_{\mathrm{spiL}}{r}_{\mathrm{bpL}}\cos \beta +F_{\mathrm{fspiL}}{H}_{\mathrm{p}i}\end{array}$$

(21)

Figure 9.

Dynamic model of ring gear–planetary gear meshing pair.

(b) Planet–ring meshing pair. The differential equation of motion of the left ring gear is as follows (Figure 10):

$$\{\begin{array}{c}{m}_{\mathrm{rL}}\left({\stackrel{\cdot \cdot }{x}}_{\mathrm{rL}}\text{-}2{\omega }_{\mathrm{c}}{\stackrel{\cdot }{y}}_{\mathrm{rL}}\text{-}{\omega }_{\mathrm{c}}^2{x}_{\mathrm{rL}}\right)+c_{\mathrm{rxL}}\stackrel{\cdot }{x_{\mathrm{rL}}}+k_{\mathrm{rxL}}{f}_{\mathrm{rxL}}(x_{\mathrm{rL}})=\sum _{i=1}^N{F}_{\mathrm{pirL}}\sin {\psi }_{\mathrm{pirL}}\cos \beta -\sum _{i=1}^N{F}_{\mathrm{fpirL}}\cos {\psi }_{\mathrm{pirL}}\cos \beta \\ m_{\mathrm{rL}}\left({\stackrel{\cdot \cdot }{y}}_{\mathrm{rL}}+2{\omega }_{\mathrm{c}}\underset{\mathrm{rL}}{\overset{\cdot }{x}}\text{-}{\omega }_{\mathrm{c}}^2{y}_{\mathrm{rL}}\right)+c_{\mathrm{ryL}}{\stackrel{\cdot }{y}}_{\mathrm{rL}}+k_{\mathrm{ryL}}{f}_{\mathrm{ryL}}(y_{\mathrm{rL}})=-\sum _{i=1}^N{F}_{\mathrm{pirL}}\cos {\psi }_{\mathrm{pirL}}\cos \beta +\sum _{i=1}^N{F}_{\mathrm{fpirL}}\sin {\psi }_{\mathrm{pirL}}\cos \beta \\ m_{\mathrm{rL}}{\stackrel{\cdot \cdot }{z}}_{\mathrm{rL}}+c_{\mathrm{rz}}{\stackrel{\cdot }{z}}_{\mathrm{rL}}+k_{\mathrm{rz}}f(z_{\mathrm{rL}})=-\sum _{i=1}^N{F}_{\mathrm{pirL}}\cos {\psi }_{\mathrm{pirL}}\sin \beta \\ I_{\mathrm{rxL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{rxL}}+c_{\mathrm{rx}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{rxL}}+k_{\mathrm{rx}\theta \mathrm{L}}{\theta }_{\mathrm{rxL}}=-\sum _{i=1}^N{F}_{\mathrm{pirL}}\sin \beta \sin {\psi }_{\mathrm{pirL}}{r}_{\mathrm{brL}}\\ I_{\mathrm{ryL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{ryL}}+c_{\mathrm{ry}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{\mathrm{ryL}}+k_{\mathrm{ry}\theta \mathrm{L}}{\theta }_{\mathrm{ryL}}=\sum _{i=1}^N{F}_{\mathrm{pirL}}\sin \beta \cos {\psi }_{\mathrm{pirL}}{r}_{\mathrm{brL}}\\ I_{\mathrm{rzL}}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{rzL}}=-\sum _{i=1}^N{F}_{\mathrm{pirL}}{r}_{\mathrm{brL}}\cos \beta +T_{\mathrm{r}}/2+\sum _{i=1}^N{F}_{\mathrm{fpirL}}\cos \beta \cdot H_{\mathrm{r}}\end{array}$$

(22)

Figure 10.

Planetary gear–carrier dynamic model.

The differential equation of motion of the left planet gear is

$$\{\begin{array}{c}{m}_{\mathrm{piL}}\left({\stackrel{\cdot \cdot }{x}}_{piL}\text{-}2{\omega }_c{\stackrel{\cdot }{y}}_{piL}\text{-}{\omega }_{\mathrm{c}}^2{x}_{piL}\right)+c_{pixL}{\stackrel{\cdot }{x}}_{piL}+k_{pixL}{f}_{pixL}(x_{piL})=-F_{pirL}\sin {\psi }_{pirL}\cos \beta +F_{fpirL}\cos {\psi }_{pirL}\cos \beta \\ m_{\mathrm{piL}}\left({\stackrel{\cdot \cdot }{y}}_{\mathrm{piL}}+2{\omega }_c{\stackrel{\cdot }{x}}_{piL}\text{-}{\omega }_{\mathrm{c}}^2{y}_{piL}\right)+c_{piyL}{\stackrel{\cdot }{y}}_{piL}+k_{piyL}{f}_{piyL}(y_{piL})=F_{pirL}\cos {\psi }_{pirL}\cos \beta -F_{fpirL}\sin {\psi }_{pirL}\cos \beta \\ m_{piL}{\stackrel{\cdot \cdot }{z}}_{piL}+c_{pizL}{\stackrel{\cdot }{z}}_{piL}+k_{pizL}f(z_{piL})=F_{pirL}\cos {\psi }_{pirL}\sin \beta \\ I_{pixL}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pixL}}+c_{\mathrm{pix}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{pixL}+k_{\mathrm{pix}\theta \mathrm{L}}{\theta }_{pixL}=F_{pirL}\sin \beta \sin {\psi }_{pirL}{r}_{brL}\\ I_{piyL}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{piyL}}+c_{\mathrm{piy}\theta \mathrm{L}}{\stackrel{\cdot }{\theta }}_{piyL}+k_{p\mathrm{iy}\theta \mathrm{L}}{\theta }_{piyL}=-F_{pirL}\sin \beta \cos {\psi }_{pirL}{r}_{brL}\\ I_{pizL}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pizL}}=F_{pirL}{r}_{brL}\cos \beta -F_{fpirL}\cos \beta \cdot H_{\mathrm{pi}1}\end{array}$$

(23)

(c) Carrier–planet gear kinematic differential equation. The differential equation of motion of the carrier is as follows:

$$\{\begin{array}{c}{m}_c{\stackrel{\cdot \cdot }{x}}_c=-k_{bx}[(x_c-x_{pi})-r_c{\theta }_{zc}\sin {\alpha }_{pi}]-c_{bx}\left[\left({\stackrel{\cdot }{x}}_c-{\stackrel{\cdot }{x}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\\ m_c{\stackrel{\cdot \cdot }{y}}_c=-k_{by}[(y_c-y_{pi})+r_c{\theta }_{zc}\sin {\alpha }_{pi}]-c_{by}\left[\left({\stackrel{\cdot }{y}}_c-{\stackrel{\cdot }{x}}_{pi}\right)+r_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\\ m_c{\stackrel{\cdot \cdot }{z}}_c=-k_{bz}[(z_c-z_{pi})-r_c{\theta }_{yc}\cos {\alpha }_{pi}+r_c{\theta }_{xc}\cos {\alpha }_{pi}]-c_{bz}\left[\left({\stackrel{\cdot }{z}}_c-{\stackrel{\cdot }{z}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{yc}\cos {\alpha }_{pi}+r_c{\stackrel{\cdot }{\theta }}_{xc}\cos {\alpha }_{pi}\right]\\ I_{cx}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{cx}}=-k_{b\theta x}({\theta }_{xc}-{\theta }_{xpi})-c_{b\theta x}({\stackrel{\cdot }{\theta }}_{xc}-{\stackrel{\cdot }{\theta }}_{xpi})-r_c\sin {\alpha }_{pi}{k}_{bz}[(z_c-z_{pi})-r_c{\theta }_{yc}\cos {\alpha }_{pi}+r_c{\theta }_{xc}\cos {\alpha }_{pi}]\\ -r_c\sin {\alpha }_{pi}{c}_{bz}\left[\left({\stackrel{\cdot }{z}}_c-{\stackrel{\cdot }{z}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{yc}\cos {\alpha }_{pi}+r_c{\stackrel{\cdot }{\theta }}_{xc}\cos {\alpha }_{pi}\right]\\ I_{cy}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{cy}}=-k_{b\theta y}({\theta }_{yc}-{\theta }_{ypi})-c_{b\theta y}({\stackrel{\cdot }{\theta }}_{yc}-{\stackrel{\cdot }{\theta }}_{ypi})+r_c\cos {\alpha }_{pi}{k}_{bz}[(z_c-z_{pi})-r_c{\theta }_{yc}\cos {\alpha }_{pi}+r_c{\theta }_{xc}\cos {\alpha }_{pi}]\\ +r_c\cos {\alpha }_{pi}{c}_{bz}\left[\left({\stackrel{\cdot }{z}}_c-{\stackrel{\cdot }{z}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{yc}\cos {\alpha }_{pi}+r_c{\stackrel{\cdot }{\theta }}_{xc}\cos {\alpha }_{pi}\right]\\ I_{cz}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{cz}}=r_c\sin {\alpha }_{pi}\left\{k_{bx}[(x_c-x_{pi})-r_c{\theta }_{zc}\sin {\alpha }_{pi}]+c_{bx}\left[\left({\stackrel{\cdot }{x}}_c-{\stackrel{\cdot }{x}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\right\}\\ -r_c\cos {\alpha }_{pi}\left\{k_{by}[(y_c-y_{pi})+r_c{\theta }_{zc}\sin {\alpha }_{pi}]+c_{by}\left[\left({\stackrel{\cdot }{y}}_c-{\stackrel{\cdot }{x}}_{pi}\right)+y_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\right\}+T_c\end{array}$$

(24)

Differential equation of motion of planet i

$$\{\begin{array}{c}{m}_p{\stackrel{\cdot \cdot }{x}}_{pi}=k_{bx}[(x_c-x_{pi})-r_c{\theta }_{zc}\sin {\alpha }_{pi}]+c_{bx}\left[\left({\stackrel{\cdot }{x}}_c-{\stackrel{\cdot }{x}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\\ m_p{\stackrel{\cdot \cdot }{y}}_{pi}=k_{by}[(y_c-y_{pi})+r_c{\theta }_{zc}\sin {\alpha }_{pi}]+c_{by}\left[\left({\stackrel{\cdot }{y}}_c-{\stackrel{\cdot }{x}}_{pi}\right)+y_c{\stackrel{\cdot }{\theta }}_{zc}\sin {\alpha }_{pi}\right]\\ m_p{\stackrel{\cdot \cdot }{z}}_{pi}=k_{bz}[(z_c-z_{pi})-r_c{\theta }_{yc}\cos {\alpha }_{pi}+r_c{\theta }_{xc}\cos {\alpha }_{pi}]+c_{bz}\left[\left({\stackrel{\cdot }{z}}_c-{\stackrel{\cdot }{z}}_{pi}\right)-r_c{\stackrel{\cdot }{\theta }}_{yc}\cos {\alpha }_{pi}+r_c{\stackrel{\cdot }{\theta }}_{xc}\cos {\alpha }_{pi}\right]\\ I_{px}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pxi}}=k_{b\theta x}({\theta }_{xc}-{\theta }_{xpi})+c_{b\theta x}({\stackrel{\cdot }{\theta }}_{xc}-{\stackrel{\cdot }{\theta }}_{xpi})\\ I_{py}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pyi}}=k_{b\theta y}({\theta }_{yc}-{\theta }_{ypi})+c_{b\theta y}({\stackrel{\cdot }{\theta }}_{yc}-{\stackrel{\cdot }{\theta }}_{ypi})\\ I_{pz}{\stackrel{\cdot \cdot }{\theta }}_{\mathrm{pzi}}=k_{b\theta z}({\theta }_{zc}-{\theta }_{zpi})+c_{\mathrm{b}\theta \mathrm{z}}({\stackrel{\cdot }{\theta }}_{zc}-{\stackrel{\cdot }{\theta }}_{zpi})-T_c\end{array}$$

(25)

The double-helical gear is connected by the left and right helical gears through the intermediate shaft segment, the Timoshenko beam element is equivalent to the intermediate shaft segment of the double-helical gear in the study. ³³

Lubrication model of double-helical gear pair

The coordinate system of the lubrication model

A herringbone gear transmission system is more widely used than spur and helical gears in high-power density conditions, such as in aircraft engines, large-scale warships, and gas turbines by its much higher load-carrying capacity and transmission stability. Due to the high loads and surface roughness that herringbone gears are subjected to, as well as the tooth profile modification commonly used in gear transmission, lubrication point contacts are formed on the gear surfaces. Although there is much literature on spur gear EHL and some on helical gear EHL performances, the previous research studies on helical or herringbone gear EHL were conducted without utilizing a comprehensive point contact EHL model and considering tooth surface characteristics. Meanwhile, the lubrication performance for heavy-loaded herringbone gears with profile modification in transient non-Newtonian thermal micro-EHL point contact is scarcely studied by researchers.^30,34

A herringbone gear drive with a clearance groove can be equivalent to a combination of two helical gear drives with opposite tooth directions if the machining and assembly errors are neglected. EHL analysis for crowned herringbone gears in point contact is also implemented according to the EHL solution to the equivalent crowned helical gears in point contacts (seen in Figure 11). Here, the coordinate origin is the midpoint of the contact line, the x-axis is vertical to the contact line and is in the tangent plane, the y-axis coincides with the contact line, and the z-axis is vertical to the contact line and is in the normal plane of the pinion. According to Figure 12, the contact domain of the gear drive is a narrow ellipse.

Figure 11.

Lubrication model of point contact.

Figure 12.

Contact domain of the gear drive.

The equations of the lubrication model

(a) Reynolds equation.

$${\displaystyle \frac{\mathrm{\partial }}{\mathrm{\partial }x}({\displaystyle \frac{\rho h^3}{\eta }{\displaystyle \frac{\mathrm{\partial }p}{\mathrm{\partial }x})+{\displaystyle \frac{\mathrm{\partial }}{\mathrm{\partial }y}({\displaystyle \frac{\rho h^3}{\eta }{\displaystyle \frac{\mathrm{\partial }p}{\mathrm{\partial }y})=-12u{\displaystyle \frac{\mathrm{\partial }(\rho h)}{\mathrm{\partial }x}}}}}}}}$$

(26)

(b) Film thickness equation with surface roughness. The surface of the gear meshing contact area will produce elastic deformation, so the thickness of the oil film needs to consider the initial geometric clearance, the elastic deformation, the surface roughness, and the tooth wear level. The initial geometric clearance of the double-helical gear pair is shown in Figure 13.

Figure 13.

Geometric clearance for double-helical gear pair.

(c) Film thickness equation.

$$h=h_0+{\displaystyle \frac{x^2}{2R_x}+{\displaystyle \frac{y^2}{2R_y}+\nu (x,y)+r(x,y)+h_w}}$$

(27)

where h₀ represents the initial clearance of the contact surface, R_x and R_y represent the curvature radius of the tooth profile along the x and y directions at the meshing point, respectively, $r(x,y)$ represents the surface roughness function, $\nu (x,y)$ represents elastic deformation, and $h_{\mathrm{w}}$ represents the tooth wear level.

$$r(x,y)=a_R\cdot \delta (x,y)$$

(28)

$$\nu (x,y)={\displaystyle \frac{2}{\pi E}\underset{\mathrm{\Omega }}{\int }{\displaystyle \frac{p(x,y)}{\sqrt{x^2+y^2}}}\mathrm{d}x\mathrm{d}y}$$

(29)

$$h_w=h_{i,1}+h_{i,2}$$

(30)

where $a_R$ represents the magnitude of the roughness, $\delta (x,y)$ represents a roughness function whose absolute value does not exceed 1, $h_{i,1}$ and $h_{i,2}$ represent the wear amount at the meshing point of the gear pair after the gear system operates for i cycles, respectively.

(d) Load balance equation. During the gear meshing process, the lubricating oil film was squeezed and the pressure was distributed over the entire contact area. Therefore, the load balance equation is:

$$\int \int p(x,y)\mathrm{d}x\mathrm{d}y=w$$

(31)

(e) Lubricant viscosity equation and density equation.

$$\eta ={\eta }_0\exp ((\ln {\eta }_0+9.67)(-1+(1+5.1\times {10}^{-9}p{)}^{z_0}({\displaystyle \frac{T-138}{T_0-138}{)}^{-s_0}))}$$

(32)

$$\begin{array}{r}\rho ={\rho }_0(1+{\displaystyle \frac{0.6\times {10}^{-9}p}{1+1.7\times {10}^{-9}p})(1-{\epsilon }_{\mathrm{T}}(T-T_0))}\end{array}$$

(33)

where ${\epsilon }_{\mathrm{T}}$ represents the coefficient of thermal expansion, z₀ represents the viscosity coefficient of lubricating oil, T₀ represents the coefficient of thermal expansion. s₀ represents the viscosity–temperature coefficient of lubricating oil.

(f) Energy equation.

$$\{\begin{array}{c}{c}_1{\rho }_1{v}_1{\displaystyle \frac{\mathrm{\partial }T}{\mathrm{\partial }x}=k_1{\displaystyle \frac{{\mathrm{\partial }}^{2}T}{\mathrm{\partial }{z}^2}}}\\ c_2{\rho }_2{v}_2{\displaystyle \frac{\mathrm{\partial }T}{\mathrm{\partial }x}=k_2{\displaystyle \frac{{\mathrm{\partial }}^{2}T}{\mathrm{\partial }{z}^2}}}\end{array}$$

(34)

where, c₁, c₂, ${\rho }_1$ , ${\rho }_2$ , k₁, and k₂ represent the specific heat, density, and thermal conductivity of gear and pinion, respectively.

(g) Viscosity equation for lubricating oil. Under high-speed or heavy-load conditions, lubricating oil is subjected to high pressures and generates instant local high temperatures within a small movement space. In such cases, the lubrication mechanism of the gear pair becomes highly complex. Due to the combined effects of thermal effects, contact surface roughness, and variations in contact zone load, the viscosity of the lubricant is not only affected by pressure and temperature but also by shear force and shear strain rate. Therefore, the rheological characteristics of the lubricating oil need to be considered under high-speed and heavy-load conditions. Engineering practice has shown that the Ree–Eyring fluid model can accurately simulate the rheological characteristics of lubricating oil under these conditions. The constitutive equation is given by equation (35) ³⁵ :

$$\stackrel{\cdot }{\gamma }={\displaystyle \frac{{\tau }_0}{{\eta }_0}\sinh ({\displaystyle \frac{\tau }{{\tau }_0})}}$$

(35)

The effective viscosity of the lubricating oil satisfies the following relationship:

$${\displaystyle \frac{1}{\eta *}={\displaystyle \frac{1}{\eta }{\displaystyle \frac{{\tau }_0}{\tau }\sinh ({\displaystyle \frac{\tau }{{\tau }_0})}}}}$$

(36)

In the above equations, ${\tau }_0$ represents the shear stress value when the shear strain rate and shear stress exhibit obvious nonlinearity, which depends on the type and molecular species of the fluid. $\tau$ represents the shear stress, and $\eta *$ represents the equivalent viscosity that takes into account the non-Newtonian characteristics of the lubricating oil.

Calculation process for tri-bo dynamic model of double-helical PGT

The calculation flow of the dynamic model for the herringbone planetary transmission system is shown in Figure 14, which is coupled with the structural vibration and tribological characteristics. Firstly, the basic parameters of the lubrication model of the herringbone gear pair are calculated including the entertainment speed, radius of curvature, and slip–roll ratio, based on geometrical parameters, the working condition parameters, and the surface roughness. The friction coefficient of the corresponding working condition is calculated by solving the lubrication model, and then the obtained friction force and friction torque are substituted into the dynamic equation of the herringbone planetary transmission system. The vibration displacement, vibration velocity, vibration acceleration, and dynamic meshing force F_p of the meshing pair of the herringbone planetary gear are obtained by the Runge–Kutta method. In return, using the obtained dynamic meshing force as the load input condition to re-solve the lubrication equations, obtaining a new friction coefficient, and substitute the updated value of the friction coefficient into the dynamic model to re-iteratively calculation until the dynamic load satisfies formula (30):

$${\displaystyle \frac{\left|\sum F_{\mathrm{p}i}(t)-\sum F_{\mathrm{p}i-1}(t)\right|}{\left|\sum F_{\mathrm{p}i}(t)\right|}\le {\epsilon }_{\mathrm{F}}}$$

(37)

Figure 14.

Calculation flowchart.

Result

Influence of different working conditions on the dynamic characteristics of the double-helical planetary transmission system

The operating conditions have a great influence on the vibration characteristics of the gear system, so it is necessary to study the effects of input power, rotational speed, and torque on the vibration characteristics of the gear system. Figure 15 shows the vibration characteristics of the external meshing pair for the double-helical planetary gear under different working conditions. Under the condition of the same input torque, the vibration displacement amplitude of the external meshing pair shows this change trend, firstly increasing rapidly, then decreasing slowly, and finally increasing rapidly with the increase of the rotational speed. The vibration acceleration amplitude of the external meshing pair shows an upward trend as a whole with the increase of the rotational speed. Under the condition of the same rotational speed, the amplitude of vibration displacement and vibration acceleration showed a slowly increasing trend with the increase of the input torque.

Figure 15.

Results of external meshing pair analysis model.

Figure 16 shows the analysis of the vibration characteristics of the internal meshing pair of the double-helical planetary gear under different working conditions. Under the condition of the same torque, the amplitude of vibration displacement first decreased and then increased with the increase of input speed, and the amplitude of vibration acceleration increased with the increase of speed. Under the condition of the same speed, the vibration displacement amplitude and vibration acceleration amplitude show an increasing trend with the increase of torque.

Figure 16.

Results of inner meshing pair dynamic model.

Influence of tooth wear on vibration characteristics for double-helical planetary gear meshing pair

Figure 17 shows the wear depth distribution of each gear in the planetary gear system along the meshing line. The results show the wear depth distribution along the meshing line while the sun gear, planetary gear, and ring gear have experienced 1 million, 2 million, and 3 million wear cycles, respectively. In the process of gear teeth from meshing in to meshing out, the wear depth first decreases and then increases. Since the sliding distance at the node is 0, the wear depth is the smallest. For planetary gears, the wear depth first decreases and then increases, and the wear amount at the tooth tip of the planetary gear is the largest, followed by the tooth root. The distribution law of the wear depth of the ring gear along the meshing line is the same as that of the planetary gear. The wear depth of the tooth profile increases with the number of cycles. Comparing the wear depth of sun gear, planetary gear, and ring gear, it can be seen that the maximum wear depth of sun gear is less than that of planetary gear but greater than that of ring gear. This is because the cycle number of the planetary gear is the most under the same operating time.

Figure 17.

Wear depth distribution: (a) distribution of wear depth of sun gear, (b) distribution of wear depth of planetary gear, and (c) wear depth distribution of ring.

Figure 18 shows the variation curves of the vibration displacement amplitude and vibration acceleration amplitude for double-helical planetary gear pair under different wear levels. When the meshing pair is in the running-in wear stage, the vibration displacement and acceleration amplitudes of the internal and external meshing pairs slowly decrease. Because the actual contact area of the new friction pair is small, but with the increase of the operation numbers of the gear system, the surface micro-roughness peaks are gradually ground away, the surface roughness value decreases, and the actual contact area increases. When entering the stable wear stage, the displacement and acceleration amplitude of the vibration of the internal and external meshing gear pairs increase very smoothly. The tooth wear at this stage is slow and steady. When the number of operations continues to increase and the gear system enters the severe wear stage, the vibration displacement and vibration acceleration amplitude of the external meshing pair and the internal meshing pair increase rapidly.

Figure 18.

Variation curve of the calculation results for gear system parameters with different wear levels: (a) vibration acceleration of external meshing pair and (b) vibration acceleration of inner meshing pair.

Influence of lubrication state on vibration characteristics for double-helical planetary gear meshing pair

Figure 19 shows the results of a dynamic model of a gear system with different lubrication conditions. It can be seen that in the state of EHL, the amplitude of vibration displacement, vibration velocity, and vibration acceleration for the inner and outer meshing pairs along the meshing line direction are the largest under different working conditions. Otherwise, in the boundary lubrication state, the amplitude of vibration displacement, vibration velocity, and vibration acceleration for the inner and outer meshing pairs along the meshing line direction are the largest under different working conditions. This is because the friction coefficient in the boundary lubrication state is larger than that in the mixed lubrication state, and the friction coefficient is the smallest in the EHL state. The results show that the larger the friction coefficient, the smaller the vibration for the gear system.

Figure 19.

Results of a dynamic model of gear system with different lubrication conditions.

Conclusion

Considering the coupling relationship between the dynamic characteristics and the tribological characteristics, a lumped translation–torsional coupled vibration model parameters model of double-helical PGT is established. Which comprehensively considers the factors, such as tooth wear, tooth surface friction excitation, and tooth surface temperature rise. Firstly, the Timoshenko beam element is used to couple the left and right sides of the double-helical gear. The dynamic differential equation of the double-helical planetary transmission system is derived. The eccentricity error of sun gear and planetary gear, time-varying mesh stiffness, tooth surface friction and wear, backlash, tooth surface roughness, and thermal effect of contact area on the dynamic behavior are considered in the model. Secondly, under the condition of the same input torque, the vibration displacement amplitude of the external meshing pair shows this change trend, firstly increasing rapidly, then decreasing slowly, and finally increasing rapidly with the increase of the rotational speed. The vibration acceleration amplitude of the external meshing pair shows an upward trend as a whole with the increase of the rotational speed. Under the condition of the same rotational speed, the amplitude of vibration displacement and vibration acceleration showed a slowly increasing trend with the increase of the input torque. When the meshing pair is in the running-in wear stage, the vibration displacement and acceleration amplitudes of the internal and external meshing pairs slowly decrease. Because the actual contact area of the new friction pair is small, but with the increase of the operation numbers of the gear system, the surface micro-roughness peaks are gradually ground away, the surface roughness value decreases, and the actual contact area increases. When entering the stable wear stage, the displacement and acceleration amplitude of the vibration of the internal and external meshing gear pairs increase very smoothly. The tooth wear at this stage is slow and steady. When the number of operations continues to increase and the gear system enters the severe wear stage, the vibration displacement and vibration acceleration amplitude of the external meshing pair and the internal meshing pair increase rapidly. Last but not least, comparing the effects of different lubrication states on the dynamic characteristics of the double-helical planetary transmission system, the larger the friction coefficient, the smaller the vibration for the gear system.