In mechanical transmission systems, spur gears are fundamental components widely used in applications such as new energy vehicles and aerospace machinery. The mesh stiffness, as a primary source of internal excitation, plays a crucial role in gear dynamics. Accurate evaluation of mesh stiffness is essential for improving transmission accuracy and optimizing gear design in these fields. Traditional methods for calculating mesh stiffness include experimental approaches, potential energy methods, finite element methods, and hybrid techniques. While potential energy methods are popular due to their computational efficiency, they often overlook the influence of driving speed and the accompanying centrifugal effects. In high-speed operations, the driving speed significantly impacts the dynamic behavior of spur gears, and centrifugal forces can alter the deformation response, thereby affecting mesh stiffness. This paper proposes an original computational algorithm based on Euler beam theory to calculate the dynamic mesh stiffness of spur gears, incorporating centrifugal effects into the velocity field. By deriving the equations of motion using Hamilton’s principle and solving them with the Newmark method, we investigate the nonlinear relationship between centrifugal effects and dynamic mesh stiffness. The results demonstrate that both the natural frequency and dynamic mesh stiffness increase with driving speed under centrifugal influence. Additionally, materials with high elastic moduli tend to suppress the impact of driving speed on dynamic mesh stiffness, while higher density materials amplify it. This study provides insights for further analysis of vibration and noise in spur gears under centrifugal effects.
The dynamic behavior of spur gears is complex due to the interaction between mechanical loads and rotational effects. In particular, centrifugal forces arising from high driving speeds can lead to significant changes in gear deformation and stiffness. Existing studies often focus on static conditions, neglecting the dynamic nature of gear meshing. For instance, potential energy methods and finite element analyses have been extended to account for profile modifications and driving speed, but they typically treat centrifugal effects quasi-statically. This limitation highlights the need for a dynamic approach that integrates centrifugal forces into the stiffness calculation. Our method addresses this by modeling spur gears as flexible Euler beams and incorporating centrifugal terms into the kinetic and potential energy expressions. The governing equations are derived using Hamilton’s principle, resulting in a system that includes mass, damping, and stiffness matrices influenced by driving speed. The dynamic mesh stiffness is then computed by solving for gear displacements under varying operating conditions, providing a more realistic representation of high-speed gear behavior.
The motion equations for a rotating flexible spur gear are derived based on Euler beam theory. The gear is discretized into beam elements, and the total displacement vector at any point accounts for both rotational and centrifugal-induced deformations. The position vector \(\mathbf{P}_T\) of a point on the gear is given by:
$$ \mathbf{P}_T = \begin{bmatrix} x – u_r + u \cos \theta – v \sin \theta \\ z_{i,x} + v_r + u \sin \theta + v \cos \theta \end{bmatrix} $$
where \(x\) is the axial coordinate, \(u\) and \(v\) are displacements due to centrifugal expansion, \(u_r\) and \(v_r\) are rotational displacements, \(\theta\) is the angular position, and \(z_{i,x}\) is the distance from the neutral axis. The velocity vector \(\dot{\mathbf{P}}_T\) is derived by differentiating \(\mathbf{P}_T\) with respect to time, incorporating the driving speed \(\dot{\theta}\). The kinetic energy \(T_i\) of the gear is expressed as:
$$ T_i = \frac{1}{2} \rho \int \dot{\mathbf{P}}_T \cdot \dot{\mathbf{P}}_T dV $$
where \(\rho\) is the density and \(dV\) is the volume element. The axial strain \(\varepsilon_{xx}\) considers geometric nonlinearities:
$$ \varepsilon_{xx} = u’ + \frac{1}{2} \left[ (u’)^2 + (v’)^2 \right] $$
where the prime denotes differentiation with respect to \(x\). The potential energy \(U\) is then:
$$ U = \frac{1}{2} \int E \varepsilon_{xx}^2 dV $$
with \(E\) being the elastic modulus. Applying Hamilton’s principle:
$$ \delta \int_{t_1}^{t_2} (U – T) dt = 0 $$
yields the equation of motion:
$$ (\mathbf{M}_t + \mathbf{M}_r) \ddot{\mathbf{X}} + (\mathbf{C}_r + \mathbf{C}_p) \dot{\mathbf{X}} + (\mathbf{K}_e + \mathbf{K}_v) \mathbf{X} = \mathbf{F} + \mathbf{F}_v $$
Here, \(\mathbf{M}_t\) and \(\mathbf{M}_r\) are the translational and rotational mass matrices, \(\mathbf{C}_r\) and \(\mathbf{C}_p\) are the Rayleigh and gyroscopic damping matrices, \(\mathbf{K}_e\) is the elastic stiffness matrix, \(\mathbf{K}_v\) is the centrifugal stiffness matrix, \(\mathbf{F}\) is the external force vector, and \(\mathbf{F}_v\) is the centrifugal force vector. The centrifugal stiffness matrix \(\mathbf{K}_v\) is symmetric and proportional to the square of the driving speed \(\dot{\theta}^2\), reflecting the centrifugal hardening effect. The matrices are defined for each beam element, with examples provided below.
The translational mass matrix \(\mathbf{M}_t\) for an element of length \(\Delta L\) and cross-sectional area \(A\) is:
$$ \mathbf{M}_t = \rho A \Delta L \begin{bmatrix}
\frac{1}{3} & 0 & 0 & \frac{1}{6} & 0 & 0 \\
0 & \frac{13}{35} & \frac{11\Delta L}{210} & 0 & \frac{9}{70} & -\frac{13\Delta L}{420} \\
0 & \frac{11\Delta L}{210} & \frac{\Delta L^2}{105} & 0 & \frac{13\Delta L}{420} & -\frac{\Delta L^2}{140} \\
\frac{1}{6} & 0 & 0 & \frac{1}{3} & 0 & 0 \\
0 & \frac{9}{70} & \frac{13\Delta L}{420} & 0 & \frac{13}{35} & -\frac{\Delta L^2}{105} \\
0 & -\frac{13\Delta L}{420} & -\frac{\Delta L^2}{140} & 0 & -\frac{\Delta L^2}{105} & \frac{\Delta L^2}{105}
\end{bmatrix} $$
The rotational mass matrix \(\mathbf{M}_r\), which accounts for rotary inertia, is:
$$ \mathbf{M}_r = \frac{\rho I}{\Delta L} \begin{bmatrix}
0 & 0 & 0 & 0 & 0 & 0 \\
0 & \frac{6}{5} & \frac{\Delta L}{10} & 0 & -\frac{6}{5} & \frac{\Delta L}{10} \\
0 & \frac{\Delta L}{10} & \frac{2\Delta L^2}{15} & 0 & -\frac{\Delta L}{10} & -\frac{\Delta L^2}{30} \\
0 & 0 & 0 & 0 & 0 & 0 \\
0 & -\frac{6}{5} & -\frac{\Delta L}{10} & 0 & \frac{6}{5} & -\frac{\Delta L}{10} \\
0 & \frac{\Delta L}{10} & -\frac{\Delta L^2}{30} & 0 & -\frac{\Delta L}{10} & \frac{2\Delta L^2}{15}
\end{bmatrix} $$
where \(I\) is the moment of inertia. The gyroscopic damping matrix \(\mathbf{C}_p\) is:
$$ \mathbf{C}_p = \rho A \Delta L \dot{\theta} \begin{bmatrix}
0 & -\frac{7}{10} & -\frac{\Delta L}{10} & 0 & -\frac{3}{10} & \frac{\Delta L}{15} \\
\frac{7}{10} & 0 & 0 & \frac{3}{10} & 0 & 0 \\
\frac{\Delta L}{10} & 0 & 0 & \frac{\Delta L}{15} & 0 & 0 \\
0 & -\frac{3}{10} & -\frac{\Delta L}{15} & 0 & -\frac{7}{10} & \frac{\Delta L}{10} \\
\frac{3}{10} & 0 & 0 & \frac{7}{10} & 0 & 0 \\
-\frac{\Delta L}{15} & 0 & 0 & \frac{\Delta L}{10} & 0 & 0
\end{bmatrix} $$
The elastic stiffness matrix \(\mathbf{K}_e\) is:
$$ \mathbf{K}_e = \frac{E}{\Delta L} \begin{bmatrix}
A & 0 & 0 & -A & 0 & 0 \\
0 & \frac{12I}{\Delta L^2} & \frac{6I}{\Delta L} & 0 & -\frac{12I}{\Delta L^2} & \frac{6I}{\Delta L} \\
0 & \frac{6I}{\Delta L} & 4I & 0 & -\frac{6I}{\Delta L} & 2I \\
-A & 0 & 0 & A & 0 & 0 \\
0 & -\frac{12I}{\Delta L^2} & -\frac{6I}{\Delta L} & 0 & \frac{12I}{\Delta L^2} & -\frac{6I}{\Delta L} \\
0 & \frac{6I}{\Delta L} & 2I & 0 & -\frac{6I}{\Delta L} & 4I
\end{bmatrix} $$
The centrifugal stiffness matrix \(\mathbf{K}_v\) is:
$$ \mathbf{K}_v = \frac{\rho A \Delta L \dot{\theta}^2}{210} \begin{bmatrix}
-70 & 0 & 0 & -35 & 0 & 0 \\
0 & K_{22}^v & K_{23}^v & 0 & K_{25}^v & K_{26}^v \\
0 & K_{23}^v & K_{33}^v & 0 & -K_{26}^v & K_{36}^v \\
-35 & 0 & 0 & -70 & 0 & 0 \\
0 & K_{25}^v & -K_{26}^v & 0 & K_{22}^v & -K_{23}^v \\
0 & K_{26}^v & K_{36}^v & 0 & -K_{23}^v & K_{36}^v
\end{bmatrix} $$
where the elements are defined as \(K_{22}^v = -156 + \frac{504I}{A \Delta L^2}\), \(K_{23}^v = -22\Delta L + \frac{42I}{A \Delta L}\), \(K_{25}^v = -54\Delta L – \frac{504I}{A \Delta L^2}\), \(K_{26}^v = 13\Delta L + \frac{42I}{A \Delta L}\), \(K_{33}^v = -4\Delta L^2 + \frac{56I}{A}\), and \(K_{36}^v = \frac{3A \Delta L^2 – 14I}{A}\). The centrifugal force vector \(\mathbf{F}_v\) is obtained by integrating over the element:
$$ \mathbf{F}_v = -\frac{\rho A \Delta L \dot{\theta}^2}{2} \left( L + \frac{\Delta L}{3}, 0, 0, L + \frac{2\Delta L}{3}, 0, 0 \right)^T $$
To compute the dynamic mesh stiffness, we employ the Newmark method to solve the equation of motion for gear displacements under varying driving speeds. The meshing process of spur gears is simulated using a single-tooth model, where the external force vector \(\mathbf{F}_i\) changes with rotation. For a meshing point \(i\), the force vector is:
$$ \mathbf{F}_i = \begin{bmatrix} 0 & 0 & 0 & \cdots & F_i \sin(\beta_i) & F_i \cos(\beta_i) & F_i \cos(\beta_i) z_{i,x} & \cdots & 0 & 0 & 0 \end{bmatrix}^T $$
where \(F_i\) is the meshing force, and \(\beta_i\) is the pressure angle at point \(i\), calculated as:
$$ \beta_i = \arccos \left( \frac{R_{bp}}{\sqrt{x_i^2 + z_{i,x}^2}} \right) – \arctan \left( \frac{z_{i,x}}{x_i} \right) $$
with \(R_{bp}\) being the base circle radius. The meshing velocity \(v_i\) at point \(i\) is:
$$ v_i = \dot{\theta} \sqrt{x_i^2 + z_{i,x}^2} $$
The time step \(\Delta t_i\) for moving between meshing points is:
$$ \Delta t_i = \frac{2 \sqrt{\Delta x_i^2 + \Delta y_i^2}}{v_i + v_{i+1}} $$
where \(\Delta x_i\) and \(\Delta y_i\) are elastic deflections. Initial conditions assume zero velocity and acceleration, with initial displacement computed using Hooke’s law: \(\mathbf{X}_1 = \mathbf{F}_1 / \mathbf{K}\). Iterating with the Newmark method yields the dynamic displacement vector \(\mathbf{X}_i\), from which the single-tooth dynamic stiffness \(k_{pi}\) for the pinion is:
$$ k_{pi} = \frac{F_i}{\Delta x_{i,x} \cos\left(\frac{\pi}{2} – \beta_i\right) + \Delta x_{i,y} \cos \beta_i} $$
Similarly, for the gear, \(k_{gi}\) is computed. The comprehensive dynamic mesh stiffness \(k_{ms}\) for a single tooth pair is:
$$ k_{ms} = \frac{k_{pi} k_{gi}}{k_{pi} + k_{gi}} $$
For double-tooth contact regions, the dynamic mesh stiffness is in series. The algorithm flowchart ensures accurate computation by iterating through meshing points while accounting for centrifugal effects.
To validate the proposed method, we compare results with Ansys simulations under quasi-static conditions. The parameters for the spur gear pair are listed in Table 1.
| Parameter | Pinion | Gear |
|---|---|---|
| Number of Teeth | 27 | 41 |
| Mass (kg) | 0.22 | 0.34 |
| Elastic Modulus (GPa) | 207 | 207 |
| Poisson’s Ratio | 0.3 | 0.3 |
| Module (mm) | 2.5 | 2.5 |
| Tooth Width (mm) | 10 | 10 |
| Pressure Angle (°) | 20 | 20 |
At a low driving speed of \(\dot{\theta}_p = 0.01\) rpm, centrifugal effects are negligible, and the single-tooth dynamic stiffness (STDS) matches the static stiffness from Ansys. However, at higher speeds like \(\dot{\theta}_p = 300\) rpm, the STDS shows greater fluctuations due to centrifugal forces, which Ansys does not fully capture. This discrepancy arises because Ansys neglects centrifugal effects in dynamic simulations and has inherent errors in Euler beam element integration. Figure 4 illustrates the comparison, where the STDS oscillates around the static stiffness, with amplitude increasing with driving speed due to dynamic excitation amplification by centrifugal effects.
Further analysis examines dynamic mesh stiffness under different driving speeds. Figure 5 compares static and dynamic mesh stiffness over a normalized meshing period \(t/t_c\). The dynamic mesh stiffness fluctuates around the static value, with fluctuations intensifying as driving speed increases. This is attributed to reduced time intervals between meshing points, preventing full recovery of elastic deformations and leading to cumulative effects. Initial meshing points exhibit larger fluctuations due to higher kinetic energy. Centrifugal effects amplify these fluctuations, especially in double-tooth contact regions where vibration energy is greater. Figure 6 shows single-tooth dynamic stiffness for various driving speeds, confirming increased volatility with speed. For instance, at 300 rpm, fluctuations are minor, but at 9,000 rpm, they become significant. Figure 7 demonstrates that dynamic mesh stiffness increases with driving speed due to centrifugal hardening, while static stiffness remains unaffected. This highlights the importance of centrifugal effects in high-speed spur gears applications.
The influence of material properties on dynamic mesh stiffness is critical for gear design. We analyze four common materials: hard aluminum alloy, cast iron, ceramics, and carbon fiber nylon, with parameters in Table 2.
| Material | Elastic Modulus (GPa) | Density (kg/m³) | Specific Modulus (m) |
|---|---|---|---|
| Hard Aluminum Alloy | 70 | 2.7 | 25.92 |
| Cast Iron | 207 | 7.89 | 26.24 |
| Carbon Fiber Nylon | 230 | 1.76 | 130.68 |
| Ceramics | 410 | 3.15 | 130.16 |
Figure 10 plots the growth rate and volatility rate of dynamic mesh stiffness against driving speed. The growth rate, defined as the percentage increase in stiffness relative to static conditions, and volatility rate, representing fluctuation amplitude, both rise with driving speed. Aluminum alloy shows the highest growth due to its low elastic modulus, which enhances centrifugal effects. Cast iron and carbon fiber nylon have similar elastic moduli, but cast iron’s higher density leads to greater volatility. Ceramics, with high elastic modulus, exhibit suppressed growth. Carbon fiber nylon, being lightweight and strong, offers better stability at high speeds. Table 4 compares the influence of mass matrices on growth rates, indicating that including rotational mass matrices improves accuracy, especially at high speeds.
| Mass Matrix | Material | Growth Rate at 2,000 rpm (%) | Growth Rate at 6,000 rpm (%) | Growth Rate at 10,000 rpm (%) |
|---|---|---|---|---|
| \(\mathbf{M}_t\) | Hard Aluminum Alloy | 14.31 | 50.24 | 70.10 |
| \(\mathbf{M}_t + \mathbf{M}_r\) | Hard Aluminum Alloy | 15.71 | 53.85 | 86.50 |
| \(\mathbf{M}_t\) | Cast Iron | 12.78 | 46.77 | 65.57 |
| \(\mathbf{M}_t + \mathbf{M}_r\) | Cast Iron | 13.45 | 49.93 | 69.02 |
| \(\mathbf{M}_t\) | Ceramics | 3.39 | 21.73 | 38.26 |
| \(\mathbf{M}_t + \mathbf{M}_r\) | Ceramics | 3.83 | 24.60 | 47.10 |
| \(\mathbf{M}_t\) | Carbon Fiber Nylon | 3.30 | 20.34 | 36.56 |
| \(\mathbf{M}_t + \mathbf{M}_r\) | Carbon Fiber Nylon | 3.30 | 21.74 | 37.71 |
Natural frequency analysis under centrifugal effects reveals additional insights. Table 2 compares natural frequencies from the proposed method and finite element analysis, showing errors below 5%, validating the approach. Campbell diagrams in Figures 8 and 9 illustrate frequency changes with driving speed. Without rotational mass matrices, certain modes like the first and third frequencies remain constant, while others increase. With rotational mass matrices, frequency splitting and mode coupling occur, such as a frequency veering at 12,790 rpm, indicating centrifugal hardening. This aligns with existing theories on centrifugal effects in rotating systems.
In conclusion, the dynamic mesh stiffness of spur gears is significantly influenced by driving speed, centrifugal effects, and material properties. The proposed method, based on Euler beam theory and Hamilton’s principle, provides accurate stiffness calculations by integrating centrifugal forces into the dynamic equations. Key findings include:
– Dynamic mesh stiffness fluctuates around static values, with amplitude increasing with driving speed due to centrifugal effects.
– Centrifugal hardening enhances stiffness and natural frequencies, particularly in double-tooth contact regions.
– Materials with high elastic moduli suppress driving speed impacts, while high density materials amplify fluctuations.
– Including rotational mass matrices improves computational accuracy at high speeds.
These results aid in optimizing spur gears for high-speed applications, reducing vibration and noise. Future work could extend this approach to helical gears or include thermal effects for comprehensive analysis.