In this work, I propose a comprehensive analytical method to calculate the time-varying meshing stiffness of herringbone gears accounting for tooth profile modification and load deformation. Herringbone gears are widely used in high-power transmission systems due to their high contact ratio, low axial force, and excellent load-carrying capacity. The time-varying meshing stiffness is a critical excitation source for vibration and noise in gear systems. To accurately predict the dynamic behavior, it is essential to incorporate the effects of profile modifications, which alter the tooth geometry and contact conditions. My approach extends the potential energy method to herringbone gears, considering the unique features of double helical gears, including the presence of a run-out groove (undercut) and the complex contact line distribution.
1. Introduction
Herringbone gears, also known as double helical gears, combine two helical gear sections with opposite helix angles on a single gear blank. This design cancels out the axial thrust forces, making them ideal for heavy-duty applications. However, the meshing process is inherently dynamic, and the time-varying mesh stiffness (TVMS) is a primary internal excitation. Tooth profile modification, such as tip relief and root relief, is commonly applied to improve load distribution and reduce transmission errors. This modification inevitably changes the tooth geometry, contact ratio, and effective stiffness. Therefore, developing an accurate method to compute the TVMS of modified herringbone gears is of great theoretical and practical importance.
Existing research has focused on spur and helical gears, with some studies on herringbone gears using finite element methods (FEM) or simplified analytical models. However, few works have simultaneously considered the geometric features of the modified tooth profile, the influence of the run-out groove, and the variation in contact ratio under load. My work aims to fill this gap by deriving a tooth surface equation for herringbone gears with profile modification and using a slice-based potential energy method to compute the TVMS. I also analyze how key parameters such as modification amount, modification length, modification curve order, run-out groove width, and input torque affect the meshing characteristics.
2. Tooth Surface Equation with Profile Modification
To model the modified tooth profile, I start with a hypothetical rack cutter whose transverse profile consists of four segments: an involute line, a tip relief curve, a root relief curve, and a transition arc, as illustrated conceptually. The tip and root relief curves are described by a power function \( y = a x^n \), where \( a \) is the coefficient determined by the modification parameters, and \( n \) is the order of the modification curve. The transverse profile equations in the cutter coordinate system \( S_0 \) are expressed as follows for each segment \( i \):
For the involute segment (segment 1, CD):
$$ \mathbf{R}_{S0,1} = \begin{bmatrix} x_D – l_1 \cos \alpha_t \\ \pm (y_D + l_1 \sin \alpha_t) \\ 0 \\ 1 \end{bmatrix} $$
For the tip relief segment (segment 2, DE):
$$ \mathbf{R}_{S0,2}(\Delta_1, h_1, n) = \begin{bmatrix} x_D + u_1 \cos \alpha_t – a_1 u_1^n \sin \alpha_t \\ \pm (y_D – u_1 \sin \alpha_t – a_1 u_1^n \cos \alpha_t) \\ 0 \\ 1 \end{bmatrix} $$
For the root relief segment (segment 3, BC):
$$ \mathbf{R}_{S0,3}(\Delta_2, h_2, n) = \begin{bmatrix} x_C + u_2 \cos \alpha_t + (-1)^{n+1} a_2 u_2^n \sin \alpha_t \\ \pm (y_C – u_2 \sin \alpha_t + (-1)^{n+1} a_2 u_2^n \cos \alpha_t) \\ 0 \\ 1 \end{bmatrix} $$
For the transition arc (segment 4, AB):
$$ \mathbf{R}_{S0,4} = \begin{bmatrix} x_{O_e} – d_1 \sin \theta \\ \pm (y_{O_e} – c_1 \cos \theta) \\ 0 \\ 1 \end{bmatrix} $$
In these equations, \( \alpha_t \) is the transverse pressure angle, \( \Delta_1 \) and \( \Delta_2 \) are the modification amounts for the tip and root, \( h_1 \) and \( h_2 \) are the modification lengths, and \( u_1, u_2 \) are parameters along the local coordinate axes. The coefficients \( a_1 \) and \( a_2 \) are defined as \( a_1 = \Delta_1 \cos \alpha_t / u_a^n \) with \( u_a = h_1 / \cos \alpha_t + \Delta_1 \sin \alpha_t \), and \( a_2 = \Delta_2 \cos \alpha_t / u_b^n \) with \( u_b = h_2 / \cos \alpha_t – \Delta_2 \sin \alpha_t \). The sign “±” distinguishes the left and right tooth flanks.
By applying coordinate transformations from the cutter to the gear, I obtain the tooth surface equation for the herringbone gear in the gear coordinate system \( S_2 \). The final tooth surface position vector for a given tooth is:
$$ \mathbf{R}_{S2,i} = \begin{bmatrix} R_{S1,ix} \cos \phi_i – R_{S1,iy} \sin \phi_i + r (\cos \phi_i + \phi_i \sin \phi_i) \\ R_{S1,ix} \sin \phi_i + R_{S1,iy} \cos \phi_i + r (\sin \phi_i – \phi_i \cos \phi_i) \\ R_{S1,iz} \end{bmatrix} $$
where \( R_{S1,ix}, R_{S1,iy}, R_{S1,iz} \) are the components of the cutter surface in the intermediate coordinate system, \( \phi_i \) is the rotation angle of the gear during generation, and \( r \) is the pitch circle radius. This equation forms the basis for computing the contact point coordinates \( x(\alpha_z) \) and \( y(\alpha_z) \) and the pressure angle \( \alpha_z \) for any axial position \( z \) and time \( t \).
3. Time-Varying Meshing Stiffness Calculation Principle
3.1 Overall Stiffness Composition
The meshing stiffness of herringbone gears is decomposed into transverse plane stiffness (in the plane perpendicular to the gear axis) and axial stiffness (along the gear axis). Under the action of a meshing force \( F \), the elastic potential energy stored in the gear pair includes contributions from Hertzian contact, tooth bending, shear, radial compression, axial bending, axial torsion, and the foundation (gear body) stiffness. The total stiffness for a single tooth pair in mesh is given by:
$$ k = \frac{1}{ \sum_{j=1}^2 \left( \frac{1}{k_{tbj}} + \frac{1}{k_{tsj}} + \frac{1}{k_{taj}} + \frac{1}{k_{tfj}} + \frac{1}{k_{abj}} + \frac{1}{k_{atj}} + \frac{1}{k_{afj}} \right) + \frac{1}{k_h} } $$
where the subscript \( j = 1, 2 \) denotes the driving and driven gears, respectively. The terms are defined as: \( k_{tb} \) = transverse bending stiffness, \( k_{ts} \) = transverse shear stiffness, \( k_{ta} \) = radial compressive stiffness, \( k_{tf} \) = gear body radial stiffness, \( k_{ab} \) = axial bending stiffness, \( k_{at} \) = axial torsional stiffness, \( k_{af} \) = gear body axial stiffness, and \( k_h \) = Hertzian contact stiffness.
3.2 Transverse Stiffness Calculation
To compute the transverse stiffness components, I use a slice method. The entire tooth width is divided into \( N \) thin slices of thickness \( \Delta z \). For each slice, the tooth is modeled as a cantilever beam. The elastic potential energies for transverse bending, transverse shear, and radial compression are:
$$ dU_{tb} = \frac{F^2}{2 d k_{tb}} = \int_{x_M}^{x(\alpha_z)} \frac{3 M_t^2}{4 E \Delta z y^3} dx $$
$$ dU_{ts} = \frac{F^2}{2 d k_{ts}} = \int_{x_M}^{x(\alpha_z)} \frac{1.2 (F \cos\beta \cos\alpha_z)^2}{4 G \Delta z y} dx $$
$$ dU_{ta} = \frac{F^2}{2 d k_{ta}} = \int_{x_M}^{x(\alpha_z)} \frac{(F \cos\beta \sin\alpha_z)^2}{4 E \Delta z y} dx $$
Here, \( M_t = F \cos\beta [ \cos\alpha_z (x(\alpha_z) – x) – \sin\alpha_z y(\alpha_z) ] \) is the bending moment in the transverse plane, \( \beta \) is the helix angle, \( E \) is Young’s modulus, \( G \) is the shear modulus, and \( \Delta z = (2B + B_t) / N \), where \( B \) is the width of each single helical gear and \( B_t \) is the width of the run-out groove.
Summing the contributions from all slices yields the total transverse stiffness components:
$$ k_{tb} = \sum_{i=1}^N \frac{1}{ \int_{x_M}^{x(\alpha_z)} \frac{3 \cos^2\beta [\cos\alpha_z (x(\alpha_z)-x) – \sin\alpha_z y(\alpha_z)]^2}{2E \Delta z y^3} dx } $$
$$ k_{ts} = \sum_{i=1}^N \frac{1}{ \int_{x_M}^{x(\alpha_z)} \frac{0.6 \cos^2\beta \cos^2\alpha_z}{G \Delta z y} dx } $$
$$ k_{ta} = \sum_{i=1}^N \frac{1}{ \int_{x_M}^{x(\alpha_z)} \frac{\cos^2\beta \sin^2\alpha_z}{2E \Delta z y} dx } $$
The Hertzian contact stiffness for each contact line segment is:
$$ k_h = \frac{\pi E L}{4(1-\mu^2)} $$
where \( \mu \) is Poisson’s ratio and \( L \) is the total contact line length. The gear body radial stiffness is computed using the formula from Sainsot et al.:
$$ k_{tf} = \sum_{i=1}^N \frac{\Delta z}{ \cos^2\beta \cos^2\alpha_z \left\{ E \left[ L^* \left( \frac{u_f}{S_f} \right)^2 + M^* \frac{u_f}{S_f} + P^* (1 + Q^* \tan^2\alpha_z) \right] \right\} } $$
where \( L^*, M^*, P^*, Q^* \) are dimension coefficients based on the gear geometry, and \( u_f, S_f \) are defined according to the cantilever beam model.
3.3 Axial Stiffness Calculation
The axial meshing force component induces bending and torsional deformations of the tooth. The axial bending moment \( M_a \) and torsional moment \( M_T \) about the X and Y axes are:
$$ M_a = F \sin\beta \left( \frac{1}{N_s} \sum_{i=1}^N x(\alpha_z) – x \right) $$
$$ M_T = F \sin\beta \left( \frac{1}{N_s} \sum_{i=1}^N y(\alpha_z) \right) $$
Here, \( N_s \) is the number of engaged tooth slices. The axial bending and torsional stiffnesses are derived from the potential energy:
$$ k_{ab} = \frac{1}{ \int_{x_M}^{x(\alpha)} \frac{3 \sin^2\beta \left( \frac{1}{N_s} \sum_{i=1}^N x(\alpha_z) – x \right)^2}{E B^3 y} dx } $$
$$ k_{at} = \frac{1}{ \int_{x_M}^{x(\alpha)} \frac{3 \sin^2\beta \left( \frac{1}{N_s} \sum_{i=1}^N y(\alpha_z) \right)^2}{G(4B y^3 + B^3 y / \cos^2\beta)} dx } $$
The axial foundation stiffness \( k_{af} \) accounts for the bending deformation of the gear body, including the run-out groove region. Its potential energy is:
$$ U_{af} = \frac{F^2}{2 k_{af}} = \int_0^{r_f} \frac{M_{af}^2}{2E I_{af}} dx_1 $$
where \( M_{af} = F \sin\beta \left( \frac{1}{N_s} \sum_{i=1}^N x(\alpha_z) – x_1 \right) \) and \( I_{af} \) is the area moment of inertia of the gear body at distance \( x_1 \) from the center. The expression for \( I_{af} \) depends on the position relative to the hole radius \( r_0 \), the run-out groove radius \( r_t \), and the root radius \( r_f \).
4. Meshing Characteristics Under Load
4.1 Determination of Actual Meshing Start and End Points
To evaluate the effect of profile modification and load deformation on the contact state, I determine the actual start and end points of meshing. The total deformation \( \delta(\alpha_z) \) along the line of action for the engaged tooth pairs is:
$$ \delta(\alpha_z) = \frac{F \cos\beta}{\sum_{m=1}^{N_c} k_m(\alpha_z)} $$
where \( N_c \) is the number of simultaneously engaged tooth pairs. The geometric deviation \( \Delta(\alpha_z) \) caused by the profile modification is computed from the modified tooth surface. The actual meshing points are those where the sum of the deformation and the geometric deviation is positive, i.e., the tooth flanks are in contact. I implemented an iterative process that discretizes the tooth profile from root to tip and checks for contact conditions. This process yields the actual meshing radii \( r’_{s1} \) and \( r’_{e1} \) for the driving gear.
To validate my method, I compared the results with a finite element analysis (FEA) using ANSYS. The FE model was built with the derived modified tooth surface and a torque applied to the driving gear. The meshing start and end radii from the analytical method were 37.882 mm and 41.420 mm, while the FEA results were 38.004 mm and 41.226 mm. The errors were only 0.32% and 0.47%, respectively, confirming the accuracy of my approach.
4.2 Contact Ratio
Profile modification reduces the effective tooth thickness in the tip and root regions, leading to a decrease in the transverse contact ratio \( \varepsilon_\alpha \). The actual transverse contact ratio for modified herringbone gears is calculated as:
$$ \varepsilon’_\alpha = \frac{ \sqrt{ r’^2_{e1} – r^2_{b1} } + \sqrt{ r’^2_{s2} – r^2_{b2} } – \sqrt{ \frac{a^2}{(1+u)^2} – r^2_{b1} } – \sqrt{ \frac{u^2 a^2}{(1+u)^2} – r^2_{b2} } }{ p_{bt} } $$
where \( a \) is the center distance, \( u \) is the speed ratio, \( p_{bt} \) is the base pitch, and \( r_{b1}, r_{b2} \) are the base radii. The overlap ratio (longitudinal contact ratio) \( \varepsilon_\beta \) depends only on the helix angle and tooth width and is not affected by profile modification:
$$ \varepsilon_\beta = \frac{B \sin\beta}{\pi m_n} $$
The total contact ratio is the sum: \( \varepsilon = \varepsilon’_\alpha + \varepsilon_\beta \).
4.3 Time-Varying Contact Line Length
The contact line length \( L \) varies with the meshing position. For modified herringbone gears, the actual contact line length depends on the transverse and overlap ratios. The expression for \( L \) is piecewise:
If \( \varepsilon’_\alpha \ge \varepsilon_\beta \):
$$ L(t) = \begin{cases} \frac{2vt}{\sin\beta_b} & 0 \le vt \le \varepsilon_\beta p_{bt} \\ \frac{2B}{\cos\beta_b} & \varepsilon_\beta p_{bt} \le vt \le \varepsilon’_\alpha p_{bt} \\ \frac{2B}{\cos\beta_b} – \frac{2vt – 2\varepsilon’_\alpha p_{bt}}{\sin\beta_b} & \varepsilon’_\alpha p_{bt} \le vt \le (\varepsilon_\beta + \varepsilon’_\alpha) p_{bt} \end{cases} $$
If \( \varepsilon’_\alpha < \varepsilon_\beta \):
$$ L(t) = \begin{cases} \frac{2vt}{\sin\beta_b} & 0 \le vt \le \varepsilon’_\alpha p_{bt} \\ \frac{2\varepsilon’_\alpha p_{bt}}{\sin\beta_b} & \varepsilon’_\alpha p_{bt} \le vt \le \varepsilon_\beta p_{bt} \\ \frac{2\varepsilon’_\alpha p_{bt}}{\sin\beta_b} – \frac{2vt – 2\varepsilon_\beta p_{bt}}{\sin\beta_b} & \varepsilon_\beta p_{bt} \le vt \le (\varepsilon_\beta + \varepsilon’_\alpha) p_{bt} \end{cases} $$
Here, \( v = 2\pi n_j r_{bj} \) is the velocity along the line of action, and \( \beta_b \) is the base helix angle.
5. Results and Discussion: Parametric Influence Analysis
I used the proposed method to analyze the effects of various design parameters on the meshing stiffness and contact ratio of a specific herringbone gear pair. The basic parameters of the gear pair are summarized in Table 1.
| Parameter | Value | Parameter | Value |
|---|---|---|---|
| Number of teeth \( Z_1 \) (driving) | 34 | Number of teeth \( Z_2 \) (driven) | 31 |
| Normal module \( m_n \) (mm) | 2 | Normal pressure angle \( \alpha_n \) (°) | 22.5 |
| Helix angle \( \beta \) (°) | 30 | Overall gear width (mm) | 58 |
| Center distance \( a \) (mm) | 75.5 | Input torque (N·m) | 500 |
| Run-out groove width \( B_t \) (mm) | 10 | Young’s modulus \( E \) (Pa) | 2.1×10¹¹ |
| Poisson’s ratio \( \mu \) | 0.3 | Input speed \( n \) (rpm) | 7500 |
5.1 Effect of Run-Out Groove Width
I varied the run-out groove width \( B_t \) from 0 to 40 mm while keeping the single helical gear width constant. The results are shown in Table 2. As \( B_t \) increases, the gear body becomes stiffer, leading to an increase in the mean mesh stiffness. The contact ratio remained unchanged because the groove width does not affect the tooth geometry itself. The fluctuation amplitude of the mesh stiffness was also relatively insensitive to the groove width.
| Groove Width \( B_t \) (mm) | 0 | 10 | 20 | 30 | 40 |
|---|---|---|---|---|---|
| Transverse contact ratio \( \varepsilon’_\alpha \) | 1.101 | 1.101 | 1.101 | 1.101 | 1.101 |
| Total contact ratio \( \varepsilon \) | 3.011 | 3.011 | 3.011 | 3.011 | 3.011 |
| Mean mesh stiffness (N/(mm·μm)) | 17.617 | 17.650 | 17.866 | 18.294 | 18.772 |
| Stiffness fluctuation (N/(mm·μm)) | 1.777 | 1.781 | 1.782 | 1.773 | 1.777 |
5.2 Effect of Modification Amount
I studied the influence of the modification amount \( \Delta_g \) (applied equally to tip and root) on the meshing characteristics. As shown in Table 3, increasing the modification amount reduces the tooth thickness in the contact regions, leading to a decreased contact ratio and lower mean mesh stiffness. The stiffness fluctuation peaked when the total contact ratio was close to an integer value (e.g., near 3.0). For example, with a modification amount of 10 μm, the total contact ratio was 3.011, resulting in a relatively high fluctuation. When the modification amount was too large (15 μm), the transverse contact ratio dropped below 1.0, which could cause impact and reduce smoothness.
| Modification Amount \( \Delta_g \) (μm) | 0 | 5 | 10 | 12 | 15 |
|---|---|---|---|---|---|
| Transverse contact ratio \( \varepsilon’_\alpha \) | 1.260 | 1.249 | 1.101 | 1.035 | 0.982 |
| Total contact ratio \( \varepsilon \) | 3.170 | 3.159 | 3.011 | 2.945 | 2.892 |
| Mean mesh stiffness (N/(mm·μm)) | 18.324 | 18.066 | 16.923 | 16.141 | 15.421 |
| Stiffness fluctuation (N/(mm·μm)) | 1.660 | 1.941 | 2.039 | 1.893 | 1.570 |
5.3 Effect of Modification Length
The modification length \( h_g \) was varied from 0.6 to 1.2 mm. The results in Table 4 indicate that a longer modification zone reduces the effective meshing area, lowering the contact ratio and the mean stiffness. The stiffness fluctuation was highest when the total contact ratio was close to an integer, e.g., at \( h_g = 1.0 \) mm where \( \varepsilon = 3.011 \).
| Modification Length \( h_g \) (mm) | 0.6 | 0.8 | 1.0 | 1.2 |
|---|---|---|---|---|
| Transverse contact ratio \( \varepsilon’_\alpha \) | 1.228 | 1.153 | 1.101 | 1.013 |
| Total contact ratio \( \varepsilon \) | 3.138 | 3.063 | 3.011 | 2.923 |
| Mean mesh stiffness (N/(mm·μm)) | 17.806 | 17.226 | 16.923 | 16.126 |
| Stiffness fluctuation (N/(mm·μm)) | 1.489 | 1.942 | 2.039 | 1.849 |
5.4 Effect of Modification Curve Order
The modification curve order \( n \) was set to 2, 3, 4, and 6. As shown in Table 5, a higher order modification curve results in a larger curvature near the tooth tip, which means the tooth thickness is reduced less at the same modification length and amount. Consequently, the contact ratio and mean mesh stiffness increase with the order. However, for orders higher than 4, the change becomes negligible because the modification amounts are very small (on the order of micrometers), and the geometric difference is minimal. The stiffness fluctuation also tends to increase slightly with the order, as the contact becomes smoother.
| Modification Curve Order \( n \) | 2 | 3 | 4 | 6 |
|---|---|---|---|---|
| Transverse contact ratio \( \varepsilon’_\alpha \) | 1.196 | 1.239 | 1.249 | 1.251 |
| Total contact ratio \( \varepsilon \) | 3.106 | 3.149 | 3.159 | 3.161 |
| Mean mesh stiffness (N/(mm·μm)) | 17.885 | 18.035 | 18.069 | 18.082 |
| Stiffness fluctuation (N/(mm·μm)) | 1.705 | 1.865 | 1.938 | 1.968 |
5.5 Effect of Input Torque
I varied the input torque from 250 to 875 N·m. The results are summarized in Table 6. As the torque increases, the elastic deformation of the teeth becomes larger, allowing more of the modified tip and root regions to come into contact. This increases the effective contact ratio and mean mesh stiffness. At a torque of 500 N·m, the total contact ratio was nearly equal to 3, leading to a high stiffness fluctuation. For torque values above 750 N·m, the modification regions were fully engaged, and the contact ratio and stiffness reached a plateau, indicating that the influence of torque becomes insignificant. At very low torques (e.g., 250 N·m), the transverse contact ratio was less than 1.0, which could cause severe meshing impact.
| Input Torque (N·m) | 250 | 375 | 500 | 625 | 750 | 875 |
|---|---|---|---|---|---|---|
| Transverse contact ratio \( \varepsilon’_\alpha \) | 0.948 | 1.024 | 1.101 | 1.175 | 1.253 | 1.254 |
| Total contact ratio \( \varepsilon \) | 2.858 | 2.934 | 3.011 | 3.085 | 3.163 | 3.164 |
| Mean mesh stiffness (N/(mm·μm)) | 15.153 | 16.086 | 16.923 | 17.792 | 18.103 | 18.116 |
| Stiffness fluctuation (N/(mm·μm)) | 1.370 | 1.707 | 2.039 | 1.654 | 1.697 | 1.701 |
6. Validation of the Proposed Method
To further verify my analytical model, I compared the TVMS results with those from a three-dimensional finite element model. The FE model included the exact modified tooth geometry and the run-out groove. The comparison, shown conceptually, demonstrates excellent agreement in the trends and magnitudes of the mesh stiffness. The mean mesh stiffness from my analytical method was 17.885 N/(mm·μm), while the FE result was 18.570 N/(mm·μm), yielding a small error of 3.83%. This error is well within the acceptable range for engineering applications, confirming the robustness of the proposed calculation method for herringbone gears.
7. Conclusions
In this work, I developed a comprehensive analytical method for calculating the time-varying meshing stiffness of herringbone gears with tooth profile modification. The method accounts for the geometric features of the modified tooth profile, the effect of load-induced deformation on the contact state and contact ratio, and the presence of the run-out groove. The key findings are summarized as follows:
1. The proposed method, based on the potential energy approach and a slice model, accurately predicts the TVMS of modified herringbone gears. The results are verified by finite element analysis, with errors under 4%.
2. Profile modification, including tip and root relief, reduces the tooth thickness in the modified regions. Under load, only part of the modified region participates in meshing, leading to a decrease in the transverse contact ratio and the mean mesh stiffness compared to unmodified gears.
3. The influence of the run-out groove is significant: increasing the groove width increases the gear body stiffness and thus the mean mesh stiffness, but does not affect the contact ratio or the stiffness fluctuation amplitude.
4. The mean mesh stiffness of modified herringbone gears decreases with increasing modification amount and modification length, but increases with the order of the modification curve. The effect of the curve order becomes negligible beyond the 4th order due to the small magnitudes involved.
5. The input torque has a substantial impact on the TVMS. As torque increases, the contact ratio and mean stiffness rise due to increased tooth deformation, until the modified regions are fully engaged, after which the stiffness stabilizes. The stiffness fluctuation is most pronounced when the total contact ratio is close to an integer value.
These insights provide a valuable theoretical basis for the design and optimization of herringbone gears, helping to reduce vibration and improve the dynamic performance of transmission systems.