Analytical and FE Determination of the Change of Single Stiffness for Cylindrical Gears with External Involute Teeth

The calculation of mesh stiffness with required accuracy is essential for determining the contact characteristics of gear pairs. The easiest approximation of the relative stiffness for the basic profile geometries is the so-called single stiffness. Standardized and analytical methods for the determination of the single and mesh stiffness of gears are used to achieve design goals considering the load capacity and the vibration excitation characteristics. Such methods involve the formulas of ISO 6336-1:2006 based on experimental relationships and the equations of Weber and Banaschek based on mechanical calculations. In this paper, guidelines are given to refine the analytical calculations. Our goal is to present the impact of the change of the applied pressure angle, module, load, rim thickness and tooth number on the maximal single stiffness. The profile geometry of the gears is generated with our program in  MATLAB. The profile of gears is calculated by the tool geometry and the kinematics of production. The geometry is imported into Abaqus. The sensitivity of the models to different parameters is examined and compared to those obtained by analytical calculations. The benchmarks for the single stiffness are the two most widely used analytical calculation methods in Europe such as ISO 6336-1:2006 formulas and Weber and Banaschek equations.


Introduction
Gears used in precision power drives must meet several requirements, covering both the load capacity and their vibration excitation characteristics. It is essential to take into due account the mesh stiffness of the gears in order to meet the expectations. The basic parameter of the calculations is the single stiffness. The characteristic of the load distribution between the teeth and over the face width of the teeth is of great importance in the design of the gear drives. To determine accurately the connection stresses, the stiffness of the mesh in the given engagement position is indispensable. The standardized basic parameter used to characterize it is the single stiffness of the contact. The single stiffness parameter represents the load required to produce 1 μm normal deformation over 1 mm face width of one or more pairs of deviation-free teeth in the contact. Estimating the level of vibration excitation is important issue that also requires a more precise definition of single tooth stiffness. The importance of these calculations has greatly increased due to the spread of electric drives in the automotive industry. The goal is to estimate in advance the risk of potential acoustic problems in the application of gear designs.
The available analytical methods can be divided into two main groups: the standardized methods that combine the basic calculation procedures [1], and the analytical contact analysis. The second uses a much more complex mathematical solution in describing the single stiffness, which is typically based on the deformation calculation described by Weber and Banaschek (WB) [2]. The use of this solution allows point-to-point direct determination of the stiffness during the connection. The WB formulas are defined in plain strain model. The modification of the method according to the plane stress model is detailed by Lutz's [3] for the calculation of worm wheel deformation.
It should be noted that in the analytical modelling of the mesh stiffness of the connection, it is not sufficient to determine the single stiffness. Several procedures are available to describe the determination of the cross-sectional load distribution. One of the most widely used methods is the calculation described by Kagawa [4] and Hayashi [5], it was combined by Schmidt [6] with the WB [2] method. New results on analytical calculations can be found in [7,8]. In [7] the calculation method is compared with the Fernández et al.'s [9] and Chen and Shao's [10] methods. The determination of dynamic forces during the connection is presented by Gerber [11]. These forces and the mesh stiffness are affected by the applied geometry, torques and speeds, the lubrication characteristics, the condition of the contacting surfaces, etc.
The influence of the peeling on the mesh stiffness is discussed in [12,13]. The effect of the spalling on the friction characteristics can be seen in Saxena et al.'s research [14]. The results reveal the influence of spalling location, size and shape. He et al. [15] dealt with the theoretical approximation of contact friction and its correlation with measurement results. In [15], the authors gave the correlation analysis applying as a benchmark the Rebbechi et al.'s work [16]. The determination of damping of the lubricated gear connection can be found in [17]. The impact of additional geometric errors on the system is also deserved attention. The dynamic effect of the eccentricity of helical gears and the imbalance of the rotor were examined by Zhang et al. [18]. Wang and Zhang [19] also gave results on the effect of eccentricity and profile deviations.
Other possibility for the determination of the single stiffness is the FE simulations and hybrid methods of FEM and analytical solutions which can be divided into three main groups. These are dynamic [20][21][22][23], quasi-dynamic [24] and quasi-static [25][26][27] simulations. Depending on the purpose of the tests, both two-dimensional and three-dimensional models can be used; see e.g., the works [21,22].
Our research belongs to the quasi-static calculations because most of the analytical methods are used to approximate this case. Such research can be found e.g. in the works by Hwang et al. [25] and Zhan et al. [27], they have investigated the correlation of emerging contact pressure and tooth bending stress using standard AGMA 2001-D04 through quasi-static models. The interaction effect between the tooth pairs in the contact on the single and mesh stiffness is presented by finite element method in Wanderer's work [28] by comparing the results with the analytical solution described by Ziegler's method [29].
Continuous monitoring of the stresses at stiffness calculations is also useful for reviewing the correct settings of the models since these calculations should always provide well converging absolute values for proper settings. In contrast, the results of maximum single stiffness calculation methods cannot be compared in their absolute values. The main reason for the deviation is the deformation of the hubs [30][31][32]. Therefore, the results obtained are evaluated as the effect of the change in different geometric parameters on the single stiffness. The effect of hubs on tooth deflections can be found in [33]. The impact assessment of each geometric parameter is very useful for gear designing. These can greatly facilitate the selection of a suitable design. A good example of this is the series of analysis carried out by Li [26] dealing with the effect of the addendum factor on contact pressure, tooth bending stress and mesh stiffness.
In this paper the sensitivity of the maximum single stiffness according to several parameters is presented. The results obtained by using finite element models are also compared with analytical solutions. In these tests, analytical results according to ISO 6336-1:2006 [1] and Weber and Banaschek [2] serve as benchmarks. The boundary conditions of finite element models are set according to analytical methods.
This work points out the importance of hub geometry considered. Our aim is not to approximate analytical results as much as possible, but to present the correlation of the sensitivity of the parameter at different rim thicknesses. In the paper, the full hub geometry is considered and not limited as in the ISO and Weber and Banaschek calculations. The research thus presents the topic from a particular point of view.

Profile geometry and load settings
There are several options for selecting the test version to analyze the effect of the applied geometry. The goal is to identify as comprehensive and clear-cut a series of experiments as possible. In consequence, the geometry of the meshing gears is always identical in the calculations and the single stiffness of the meshing is always evaluated in the pitch point on contact path. This definition means the comparison of the maximum single stiffness because the paired gears will show their maximum or from this slightly different stiffness in the applied load position.
The analytical determination of single stiffness of helical gear pairs can be traced back to spur gear pairs. This approximation applies to both the formulas of ISO 6336-1:2006 [1] and the equations of WB method [2]. The current calculations are also based on this approach. As a result, the value of the helix angle is always 0°.
Perhaps one of the most obvious parameters of the contact profiles is the pressure angle. Therefore first, the effect of the pressure angle used will be evaluated. The range between 15° and 25° which is the most important in practice has been selected for the analysis.
In addition to the effect of the choice of pressure angle on the single stiffness, it is also worth examining its dependence on the normal module. The tests will be carried out on gears with 1 to 5 normal modules. This investigation is especially interesting because the method used by the ISO standard ignores the effect of this parameter. WB method already takes into account the effect of this factor.
The negligibility of the load dependence of single stiffness, as in the case of the normal module, is also worth examining separately. The ISO standard considers single stiffness to be a linear function of the load if the tangential force is less than 100 N/mm, while over 100 N/mm it considers it independent. In the case of WB method, there is no simplified consideration in the effect of the load [1].
In our investigations we shall analyze the 100 N/mm minimum load range. It should be noted that although the ISO standard defines the criterion for the load dependence for tangential force, the single stiffness is calculated according to the definition using the normal force. However, since the normal load is the function of the pressure angle, the tendency of the single stiffness change is plotted as a function of the tangential force. The ISO standard defines the 300 N/mm tangential load as the base case for the formulation of the formulas. Therefore, this load case will be chosen also in this work as reference [1].
The next analyzed issue is the effect of the tooth number on the change of single stiffness. A range of tooth numbers between 35 and 105 was selected to the calculations. The gear ratio is always 1:1 in all analyzed cases. The variant used are summarized in Table 1.
It is important to note that in this case it is no longer worth relying on the comparison with analytical methods because the chosen boundary conditions radically influence the results. The main reason for this is a strong modifying effect of the selected hubs. Even if you change the normal module, it is obvious that the gears are simply enlarged on a scale, but such a clear procedure cannot be used when the tooth number varies. The theoretical considerations for the effect of the used hub geometry are discussed in more detail in Section 3. Here only specific regulations for the analysis of the impact of the tooth number are presented.
The first option for the calculations is to record the rim thickness independently of the tooth number which results in an increasingly thinner rim in relation of the gear hub by increasing the tooth number. As a result of this setup, we can predict the significant dependence of the single stiffness on the tooth number, as proved by the calculations performed. The effect of the tooth number used on single stiffness should be examined by maintaining the proportions of the hubs. This series of calculations no longer comprises the extreme change of the hubs, thus a more comprehensive comparison is possible. The SRX model case of the rim thickness is defined as: This connection indicates the notation of dimension [x m] in the abscissa on the diagrams. This approach results the direct dependence of the rim thickness on the normal module for gears with 35 teeth and it allows to be retained the proportion of the deformable rim thickness as a function of the tooth number.

Hub settings
The problem of comparability of the absolute values of each process has already been mentioned which is caused by the calculation of the deformation of the hubs. Neither the ISO nor the WB calculations consider the exact geometry of the hub. Furthermore, the deformation of the hub in the calculation of the single stiffness is interpreted only as a function of the profile geometry. The only exception is the "blank factor" of the ISO standard. Therefore, the absolute value of FE results can only be compared with knowledge of the geometry of a given hub. As a result, it is worthwhile carrying out the comparison on the tendency change of the absolute value. The aim is to show the difference in the forecast of the expected single stiffness by using the analytical calculations or the FE method at different rim thicknesses. Even though the deformation of the hubs is taken into account only in a simplified form, the analytical methods represent the relative deformation of the individual teeth well, as they accurately reflect the relative distortion of the teeth and the segment of the gear bodies that is in contact with the teeth. Therefore, if the deflection of the hubs is properly considered when determining the spatial position of the teeth, analytical methods provide a suitable tool for determining the change of connection characteristics in the contact zone.
The extra possibilities of the calculations of the deformation provided by the finite element method are of importance in the vibration analysis of the system since it is necessary to represent the specific stiffness of the mesh. During the tests performed in Section 3, the change of single stiffness of a cylindrical gear is exhibited considering the deformation of the hubs for different rim thicknesses.
Since the hub itself can be divided into a set of elements whose deformation is determined by the stress at their contact interfaces, it is worth modelling the gears as elements with different rim thicknesses. Thus, if the approximation of deformable elements as a series of springs coupled is accepted, the stiffness change of the gear pairs can be expressed as the change in the stiffness of the modified range of the gears. In the examinations it is useful to consider the relatively small area of the rim for analyzing the single stiffness of the variations of profile geometries since the modifications made will have an effect in this range. This case corresponds to the model SR4 which means that rim thickness has been chosen as the quadruple of the normal module. A further reduction of the rim thickness may influence on the tooth root stress and on the deformation of the teeth. In the model SR4, there is still a safety factor to the borderline case. The maximum thickness of rims has been taken into account during the analyses and it was determined as the model SR13. This position covers shaft diameter capable of conveying the 300 N tangential force with an adequate safety factor for the 35 teeth gears.

Results and discussion
The performed impact assessments are summarized in Section 4. The geometries used by FEM are produced with a self-made program. The MATLAB program defines the profile geometry with the modelling of machine kinematics on the base of Litvin's work [34]. Our program allows a direct control of the precision of geometry to avoid singularities. The generated geometry is then converted to Abaqus. The precision of FE mesh quality and of the applied geometry is particularly important for the precise approximation of the deformation at the contact point. The mesh is generated with quadratic quadrilateral elements and the required size of elements to the good convergent results is always controlled based on the size of the analyzed gears. The average mesh size of the contacted teeth is 0.4 % of the whole depth. The remaining areas of gears are meshed with a varied element size depending on the effect on the results. For the simulations is used the reduced integration method. The interaction between the gears is defined as a Hertzian contact. It means that a 0 friction coefficient is applied. The contact formulation is finite sliding. The boundary conditions of finite element models are set based on the analytical methods used for the comparison [1,2]. The material of gears is always steel.
The deformation results are evaluated by continuously monitoring the obtained stress values. Single stiffness is determined by the load-induced rotation angle of driving gear. For a better overview of the calculation results, Tables 1 to 6 summarizing the variants in the given evaluations are also provided. Fig. 1     • Δ FE SRX -WB -different between FE with SRX rim thickness and WB method • Δ zX-Y -different between the models with X and Y tooth number.        We remark that the slopes of the two extreme cases are very close to the relevant analytical solutions. The models SR13 show a difference of 1.5 % from the ISO results for the selected 20° pressure angle model and, in the extreme cases, show a difference of 2.3 %. The models SR4 are much closer to the results of the WB method. There is a maximum difference of 1.6 % in the 20° profile and 3.0 % in the extreme.
In addition to the variants with 35 teeth, it is also worthwhile to examine the pressure angle dependence of the variants with larger tooth numbers. The results in Fig. 3 show the single stiffness changes if the rim thickness of the models corresponds to the normal module multiplied by 9. This means that the gears have a narrowing rim thickness with increasing tooth number in relation to the full gear body. Fig. 3 shows that the finite element results do not exceed the predicted sensitivity by WB even in the model of the 105 tooth number. If the rim thickness of the gears is increased, the curves would be shifted toward the ISO results according to Fig. 2. As a result, we can conclude that the WB calculations typically overestimate the sensitivity of the single stiffness for the pressure angle change. The method described by ISO produces a constant slope in the examined range, which reflects the tendencies of the FE models with large rim thickness. Fig. 4 shows a typical change of the single stiffness obtained for normal module at different rim thicknesses.

Single stiffness changes as a function of module
To determine this, besides the models SR9, the analysis of the models SR4 and SR13 is required.
The differences in the profile angles are summarized in Fig. 5. It is observed that the finite element method shows a decrease in single stiffness by increasing the normal module. In contrast, ISO solutions are insensitive to this parameter. As a result, the finite element method shows −6.8-0 % difference in the examined cases compared to the ISO. The WB method typically overestimates the gradient of the stiffness change. This means a maximum deviation of +5.6 % for the current models.
It can be observed that the solutions provided by the finite element method are always located between the two analytical methods. As expected, the significance of the change in normal module to total deformation decreases with increasing rim thickness. Thus, ISO provides a relatively accurate trend for gears with large rim thickness.

Single stiffness changes as a function of load
The change in the value of single stiffness for the load dependence is illustrated in Fig. 6 for models SR9. Fig. 7 shows the dependence of the single stiffness change on the pressure angle at different loads. The results show less than 1.0 % difference in the stiffness change based on a 20° profile angle. This gives less than 1.8 % deviation between the extreme situations investigated. Fig. 8 represents the dependence on the normal module. The difference in the single stiffness change is less than 0.4 % in these cases.
On the base of the calculations performed, it can be stated that the effect of the change in load on the variants investigated is very small for both the chosen pressure angle and the chosen normal module. As a result, the single stiffness change as a function of the applied pressure angle and the normal module can be independent of the magnitude of the load at the variants analyzed. Of course, this does not match the condition of the load independence of the absolute values considered. With respect to  the analyzed models SR9, the single stiffness change as a function of load is observed between −5.0 and 2.5 % compared to the obtained values at 300 N tangential load.

Single stiffness changes as a function of rim thickness
The previous results indicate that the choice of the magnitude of rim thickness has a significant effect on the single stiffness. Fig. 9 depicts a typical single stiffness change as a function of rim thickness. The model SR9 served as a benchmark for the evaluation.
The stiffness change shown in Fig. 9 is not independent of the profile geometry used. Fig. 10 illustrates the pressure angle and module dependence on the effect of rim thickness. The results show a single stiffness change between +27 and −40 % relative to models SR9. It is noticeable that the variants with smaller normal modules and larger pressure angles have the most sensitive response to change of the rim thickness. The difference in stiffness change in the extreme cases is in the 6.5 % range for models SR4 and in the 5.0 % range for models SR13.

Single stiffness changes as a function of rim thickness
The importance of hub geometry was already discussed in Section 2. Fig. 11 shows the trend of the single stiffness change of the gear pairs at the same rim thickness. The selected thickness corresponds to the normal module multiplied by 9. The results show an increase in the single stiffness between 22.4 and 25.0 % in the range of variants with 35 and 70 teeth and an increase between 7.4 and 8.3 % in the range of variants with 70 and 105 teeth. This means a change from 29.8 to 33.3 % across the range.
According to the conclusions in Section 2, the effect of the tooth number used on single stiffness should be examined by using variants with the same proportion of the hubs. Fig. 12 summarizes the results of the single stiffness change of gear pairs at different tooth numbers versus the rim thickness. Fig. 12 shows a more intensive convergence of the obtained tendencies with increasing tooth number. Fig. 13 illustrates the single stiffness change associated with increasing the tooth number. It can be seen that the use of greater rim thickness cases results in a growing difference and that the difference is significantly greater for smaller tooth numbers. Accordingly, the change between the variants with 35 and 70 tooth numbers is 4.6 % in the range SR4-9 and 17.1 % in the range SR9-13. However, the models with 70 and 105 tooth numbers show only 0.5 % difference in the range SR4-9 and 3.7 % in the range SR9-13. In the latter case, a slight decrease can be seen in single stiffness at the models SR4. However, the magnitude of the difference is comparable to the expected accuracy of the calculations.

Conclusions
The investigations of the pressure angle, module and load dependency on the single stiffness indicate that the finite element simulations show a parameter sensitivity between the ISO standard and the WB calculation for the analyzed geometric variants. The correlation rate between the individual solutions depends primarily on the value of the rim thickness used.
The tests have shown that under the analyzed conditions the evolution of the single stiffness as a function of the pressure angle and the normal module can be independent of the magnitude of the load.
The key effect of the rim thickness has also been confirmed, drawing attention to the significance of the size of the gear body. The investigation of the influence of the tooth number on the single stiffness has also highlighted the effect of the size of hub.
The tests performed accurately reflect the dependence of single stiffness on different parameters. The results present the correlation between the sensitivity of each method at the change of the analyzed parameters and the limitations of comparability.