# Mesh Stiffness and Dynamic Modeling and Analysis of Modified Straight Bevel Gears
## 2 Time-Varying Meshing Stiffness Calculation Model
### Lead Crown Relief  
Due to the minor modification, this article considers the modification amount as a tooth shape error and calculates the meshing stiffness with its initial value given before. The LCR is illustrated in Figure 5, with the specific formula presented as follows [32]:  
\[E_{x}=B^{2}/8R_{i} \tag{16}\]  
Figure 5: Straight level gear.  
LCR varies along the tooth width direction. The non-coupled slice method is used to calculate LCR meshing stiffness. The deformation is set as \(\delta_{j}\) at the initial engagement. Based on the number of slices \(j\). By comparing the engagement deformation and modification amount for each slice: if \(\delta>E_{x}\), the gears are considered meshing; whereas if \(\delta<E_{x}\), the gears are considered non-meshing. \(E_{xj}\) represents the modification amount for the tooth pair after slicing [33]. The quantification for each slice can be expressed as follows:  
\[\delta k=\frac{1}{\frac{1}{\delta k_{b}}+\frac{1}{\delta k_{a}}+\frac{1}{\delta k _{s}}} \tag{17}\]  
The tooth stiffness, tooth bending stiffness, tooth shear stiffness, and tooth axial compression stiffness of each tooth after slicing are \(\delta k_{b}\), \(\delta k_{s}\) and \(\delta k_{a}\), respectively.  
\[\begin{array}{l}\mathbf{K}_{P}=diag(\delta k_{1}^{p},\delta k_{2}^{p},\dots,\delta k_{n}^{p}),\ k_{hj}=k_{h}/n\\ \mathbf{E}_{x}=diag(E_{x1},E_{x2},\dots,E_{xn}),\ \mathbf{\delta}^{m}=diag( \delta_{1},\delta_{2},\dots,\delta_{n})\\ \mathbf{F}_{d}=diag(F_{1},F_{2},\dots,F_{n}),\ F_{m}=sum(F_{1},F_{2},\dots,F_ {n})\end{array} \tag{18}\]  
\(\mathbf{K}_{P}\), \(\mathbf{\delta}^{m}\), \(\mathbf{\delta}^{p}\), \(\mathbf{E}_{x}\) and \(\mathbf{F}_{d}\) are the stiffness matrix of the pinion, total deformation matrix, deformation matrix of the pinion, total modification matrix, and total meshing force matrix. \(F_{m}\) represents the resultant force of the sliced meshing forces.  
Firstly, calculate the initial value of \(\mathbf{F}_{d}\) (when the gear is unmeshed, consider \(\delta-E_{xj}=0\)):  
\[\begin{array}{l}\lambda=\frac{d\delta_{j}^{p}k_{hj}}{d\delta_{j}^{p}k_{hj}+d \delta_{j}^{p}k_{hj}+d\delta_{j}^{p}k_{hj}},\\ \delta^{p}=\lambda\delta^{m}\\ \mathbf{K}_{P}(\mathbf{\delta}^{p}-\mathbf{E}_{xj})=\mathbf{F}_{d}\end{array} \tag{19}\]  
Iterating on \(\delta^{m}\) and \(\mathbf{F}_{d}\), set a meshing force error of \(F_{e}\), and the load on the gear is \(F\). When \(F_{m}-F>F_{e}\):  
\[\mathbf{\delta}^{m}=\mathbf{\delta}^{m}\left(1+\frac{F-F_{m}}{\xi F_{m}}\right) \tag{20}\]  
\[\mathbf{\delta}^{p}=\lambda\mathbf{\delta}^{m} \tag{21}\]  
where \(\xi\) is the cycle coefficient, the right choice can speed up the calculation, choose \(\xi=8\)[6].  
Bring \(\mathbf{\delta}^{p}\) into the calculation of \(\mathbf{F}_{d}\) again.  
\[\mathbf{K}_{P}(\mathbf{\delta}^{p}-\mathbf{E}_{x})=\mathbf{F}_{d} \tag{22}\]  
If the following equation is met, the iteration is completed.  
\[|F_{m}-F|\leq F_{e} \tag{23}\]  
After the iteration is completed, LCR meshing stiffness based on the slice method is:  
\[K_{e}=\frac{1}{\frac{1}{k_{f}^{p}}+\frac{1}{k_{f}^{3}}+\frac{F}{\max(\mathbf{ \delta}^{m})}} \tag{24}\]  
\[K=\sum_{i=1}^{2}K_{ei} \tag{25}\]