# Fast Three-Dimensional Quasi-Static Analysis of Helical Gears Using the Finite Prism Method

M Guingand, J P de Vaujany, Y Icard

► To cite this version:

M Guingand, J P de Vaujany, Y Icard. Fast Three-Dimensional Quasi-Static Analysis of Helical Gears Using the Finite Prism Method. Journal of Mechanical Design, 2004, 126 (6), pp.10821088. 10.1115/1.1798212. hal-01499171

HAL Id: hal-01499171

https://hal.science/hal-01499171

Submitted on 31 Mar 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
# Fast Three-Dimensional Quasi-Static Analysis of Helical Gears Using the Finite Prism Method

This paper describes a model for analyzing external and internal helical gears. The finite prism method is used to study the elastic behavior of the structure. Contact deformations are also included in the model. Load sharing, pressure distribution, meshing stiffness and three-dimensional tooth fillet stresses are calculated at each instant. An experimental validation of the numerical model is also presented.

M. Guingand

J. P. de Vaujany

Y. Icard

Laboratoire de Dynamique des Machines et des Structures—CNRS UMR 5006, Institut National des Sciences Appliquées de Lyon, 20, avenue Albert Einstein, 69621 Villeurbanne Cedex, France

## 1 Introduction

Static-elastic analysis of the meshing of cylindrical gears is a challenging task as a number of modeling and numerical problems must be overcome. Two groups of problems are of primary interest. The first group is related to the computation of the contact pressure at the gear mesh. The pressures and the corresponding displacements within the contact area can be calculated by employing the theories of Hertz and Boussinesq. According to the actual geometry near the contact zone, either the results are given in the form of tabulated formulas, or a discretized contact area is used at the mesh where each grid square is subject to a certain contact pressure [1,2]. Then, geometrical conditions of nonbody penetration allow formulation of a linear problem with consistency conditions [3]. Once the pressures are known, the stress fields are calculated within the contact zone [4]. When more than one pair of teeth transmit the load, the principle used for the calculation of load distribution is based on the analysis of generalized multicontacts of two bodies [3,5]. Multitooth pair analyses of different studies subject show differences perhaps originating from the simplifying assumptions made [6–9].

The second group of problems of gear contact analysis are associated with bulk resistance. Stress calculations at the tooth root are made with the help of classical methods for continuous media. While many earlier studies used the material resistance method, general equations [10] are now used directly. Consequently, the finite element method (FEM) has been widely used in two- or three-dimensional (3D) modeling with a refined meshing [11–13]. Therefore, it is necessary to develop computational tools that are refine enough to capture the major aspects of contact phenomena, but yet global enough to take into account the entire architecture of the gears.

The aim of this paper is to present a method for analyzing cylindrical gears in an efficient and precise manner. Parametrical analysis results using this method have been presented earlier previously [14,15], however, the fundamentals of this method and its experimental validation have not been published. The objective here is to formulate an elastic model of a cylindrical gear pair in such a way that reasonably accurate computation can be made within moderate calculation times, allowing the use of this formulation as a potential design tool. The finite prism method (FPM) is chosen to reduce the size of the numerical problem and accelerate the calculations. This method is able to consider the exact geometry of the teeth as well as the global axis 3D positions of the gears altogether. With this method, the bulk deformations can be

modeled accurately as well as the stress distributions at the root of the teeth. In the following sections, FPM formulations are presented and are applied to gear pair problems.

## 2 The Finite Prism Method

Full analysis with the FEM is costly and even unnecessary for structures with a constant shape in one direction and simple boundary conditions. FPM can give results by using a much smaller number of input data [16]. FPM can be considered as a special form of FEM, based on geometrical approximations concerning the shape of the displacements along a given direction. A standard FEM uses polynomial displacement functions in every direction while FPM uses simple polynomial displacements only along two directions where the meshing is carried out and serial functions for extrusion along the third direction. Thus, the 3D problem is approximated in a two-dimensional (2D) space, providing significant savings in computing time and memory requirements.

The general form of the displacement function is given as a product of polynomial functions and series. Therefore, for a linear prism where the z axis corresponds to the generating helical line, a displacement function is written in the following form:

$$U = \sum_{m=1}^{r} N_{m}(x,y) F_{m}(z) \quad (1)$$

In this expression, the series are truncated to r terms. $N_m(x,y)$ is a polynomial expression including undetermined constants for the mth term. $F_m(z)$ is a series that satisfies the boundary conditions along the z-axis direction and also specifies the deflected shape of the structure in this direction.

The FPM used here is similar to the method proposed by Kantorovich [17], which is commonly used to reduce a partial differential equation to an ordinary differential equation. Cheung [18] developed the theory related to finite prisms. The displacement function U is often composed of two parts: the polynomial function $N_m(x,y)$, ruled by the shape of the transversal function and the series $F_m(z)$ satisfying the boundary conditions in the z direction. The transverse section of the prism is commonly modeled as an isoparametric eight-node element.

### 2.1 Choice of Displacement Functions.
The global elastic line of a gear tooth along its face width is assumed to be the curve of a bending beam. Therefore, function $F(z)$ is the solution of the spatial differential equation of a vibrating beam (see Fig. 1). The choice of this function is discussed by Olakorade [19]. The differential equation is given in the form
<figure><img src="image_1.png" /><figcaption>Fig. 1 Gear prism element</figcaption></figure>

<figure><img src="image_2.png" /></figure>

$$ \frac{d^4 F(z)}{dz^4} = \frac{\mu^4}{b^4} z \qquad (2) $$

where b is the length of the prism in consideration, and μ is defined in Refs. [20] and [21]. Therefore, the basic solution F(z) to Eq. (2) is given in the general form

$$ F(z) = C_1 \sin\left(\frac{\mu}{b}z\right) + C_2 \cos\left(\frac{\mu}{b}z\right) + C_3 \sinh\left(\frac{\mu}{b}z\right) + C_4 \cosh\left(\frac{\mu}{b}z\right) \qquad (3) $$

where coefficients Cᵢ are determined by applying the boundary conditions. They define the shape and amplitude of the displacements.

The length of the beam corresponds to the face width b of the gear teeth. In the particular case of “free-free” boundary condition, the prism is free at each end of the tooth

$$ F''(0) = F'''(0) = F''(b) = F'''(b) = 0 \qquad (4) $$

where F'' and F''' represent bending moment and transverse force terms, respectively and superscript (') denotes differentiation with respect to z. These conditions express the fact that the bending moment and the transverse loads are zero at the free ends of the teeth. By substituting Eq. (4) into Eq. (3), four equations are obtained to allow the determination of Cᵢ of Eq. (3)

$$ F_1(z) = 1, \quad F_2(z) = 1 - \frac{2}{b}z, \dots, \qquad (5) $$

$$ F_m(z) = \sin\left(\frac{\mu m}{b} z\right) + \sinh\left(\frac{\mu m}{b} z\right) - \alpha_m \left[ \cos\left(\frac{\mu m}{b} z\right) + \cosh\left(\frac{\mu m}{b} z\right) \right] \qquad (6) $$

where

$$ \alpha_m = \frac{\sin \mu_m - \sinh \mu_m}{\cos \mu_m - \cosh \mu_m} \qquad (7) $$

<figure><img src="image_3.png" /><figcaption>Fig. 2 Helical prism element</figcaption></figure>

$\mu_m \approx \frac{1}{2}(2m-3)\pi (m=3,4,5,\dots)$, and $\cos \mu_m \cosh \mu_m = 1$. In the other case of a tooth fixed at one end and free at the other end (case of a shaft pinion), the boundary conditions allow to write the conditions of a “fixed-free” prism

$$ F(0) = F'(0) = F''(b) = F'''(b) = 0 \qquad (8) $$

These equations show that the deflection and the rotation are zero at the fixed end and the bending moment and the transverse forces are zero at the free end. The corresponding serial functions for this case become

$$ F_m(z) = \sin\left(\frac{\mu m}{b} z\right) - \sinh\left(\frac{\mu m}{b} z\right) - \alpha_m \left[ \cos\left(\frac{\mu m}{b} z\right) - \cosh\left(\frac{\mu m}{b} z\right) \right] \qquad (9) $$

where

$$ \alpha_m = \frac{\sin \mu_m + \sinh \mu_m}{\cos \mu_m + \cosh \mu_m} \qquad (10) $$

$\mu_m \approx \frac{1}{2}(2m-1)\pi (m=1,2,3,\dots)$, and $\cos \mu_m \cosh \mu_m = -1$. Concerning the serial part of the displacement, $F_m$ and $F'_m$ are both used in the analysis of the bidimensional elasticity problem. The equation given later is based on the assumption of small deviations of the beam, in which the transverse deviation u is related to the longitudinal displacement F(z)

$$ F(z) = A \frac{du}{dz} \qquad (11) $$

<figure><img src="image_4.png" /></figure>

Direction of the prism

P = 254.51 N/mm²

v = 0.3

R1 = 10 mm

R2 = 50 mm

L = 80 mm

E = 210000 MPa

<figure><img src="image_5.png" /><figcaption>Fig. 3 Pipe model</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 4 A prismatic solid</figcaption></figure>

A general form of the serial function is given by

$$F_m^u = F_m, \quad F_m^v = F_m, \quad F_m^w = F_m' \quad (12)$$

In the case of a spur gear pair, FPM can be used directly. For helical gears, a transformation of coordinates must be applied [5]. First, a curvilinear system (O,R,θ,Z) is used as shown in Fig. 2 such that

$$X = R \cos \theta, \quad Y = R \sin \theta, \quad R = [X^2 + Y^2]^{1/2}, \\ \theta = \text{Arc tan} \left[ \frac{Y}{X} \right], \quad Z = z \quad (13)$$

with $R>0$ and $\theta_0 \le \theta \le \theta_f$. Next, an oblique system (O,r,s,t) is applied in Fig. 2 as

$$R=r, \quad \theta=s+t \sin \beta, \quad Z=t \cos \beta \quad (14)$$

## 2.2 Validation of FPM.
A pipe pressurized internally constitutes the first validation example. The pipe is fixed rigidly externally and an uniform pressure is applied inside as shown in Fig. 3. The serial function is defined along the axial direction. All the results obtained with a prism match well with the analytical predictions within an error nearly 4% on average. The larger differences are observed on the inner surface of the pipe.

In order to simulate the loading of the gear teeth, a clamped prismatic solid was considered as shown in Fig. 4 and the FPM results were compared with those obtained by using a classical 3D FEM. The main objective was to validate the results with different boundary conditions (fixed-free and free-free) and also to investigate the appropriate number of terms that must be considered in the series. The block is fixed rigidly in y direction for all z, the nodes are either free on the two small sides (case free-free), or fixed at boundary at z=0 (case of fixed-free). The stresses are analyzed on the dotted line aa' (see Fig. 4) close to the rigidly fixed side. The displacements are observed on the bb' line and the load is applied on top of the tooth. Figures 5 and 6 show the results for the free-free boundary condition along the length of the beam. The displacements (Fig. 5) and the values of the principal stress (Fig. 6) are given for different numbers of terms in the series (T1 is one-term, T2 is two-term and soon). The results are compared to those obtained with the 3D finite element method.

<figure><img src="image_7.png" /><figcaption>Fig. 5 Displacements study (x direction)</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 6 Maximum main stress study</figcaption></figure>

For the meshing, hexahedron elements with 20 nodes are used in a commercially available finite elements package (two elements in X, four elements in Y, and ten elements in Z. Different meshing conditions are also considered to analyze the effect of the meshing smoothness on the results of both between FEM and FPM. The results show that, in this example case, eight terms are necessary to achieve acceptable accuracy in displacement predictions and five terms for stress computations.

In summary, the FPM uses a 2D mesh on a transversal section of an object, and an implicit approximation of the displacements is performed along the third dimension. This method gives satisfying results, compared to those obtained by other conventional methods. In addition, the calculation time is much shorter (80 times shorter according to Olakorede [5]), and the memory space required is much smaller (nearly six times less) compared to conventional finite elements. Therefore, this method can be useful for solving gear design problems.

# 3 Application of FPM to Gear Analysis

## 3.1 Description of the Process.
A realistic geometry of the gear distinguishes three levels of definition: the architecture (web and rim), the position of the gear (misalignment, eccentricity, instantaneous meshing position), and the teeth (profile, crowning, surface waviness, pitch, and distortion errors). Obviously, for each meshing position, the load is shared among the meshing teeth, contact pressure distribution is established and a stress field spreads inside the considered elements. The complete process, integrated in the design software called PRINCE (French acronym for: Program Interactif de Calcul d'Engrenage), can be summed up as follows (Fig. 7).

(1) After determining the active profiles of the teeth, tooth contact analysis (TCA) simulation without load is performed for several meshing positions.

(2) This preliminary TCA allows the definition of the initial contact points and the unloaded transmission error caused, for example, by gear eccentricity, modified profile or undesirable profiles errors.

(3) A 2D mesh of the gear is then produced to introduce the FPM method. Influence coefficients are used to calculate load

<figure><img src="image_9.png" /><figcaption>Fig. 7 Computational process flowchart</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Fig. 8 Contact zone</figcaption></figure>

sharing. A unit force is applied on the tooth contact surface and produces a displacement. Then, by considering unit forces at different points, a set of bulk influence coefficients can be deduced. A model of the contact calculates the surface effects between the teeth. At this stage of the process, the contact pressures and the static loaded transmission error are obtained. Finally, once the load sharing is known, the FPM is used again to calculate the 3D stresses at the root of the teeth at the last step of the process.

**3.2 Load Analysis.** The determination of tooth load distribution is a problem of contacts. Once potential contact zones are established, the load sharing in these areas must satisfy the conditions of displacement compatibility in every point *k* in Fig. 8. Inside the contact zone

$$W_{1k} + W_{2k} + \epsilon_k - \alpha = 0 \quad \text{and} \quad P_k \ge 0, \quad (15)$$

and outside the contact zone

$$W_{1k} + W_{2k} + \epsilon_k - \alpha \ge 0 \quad \text{and} \quad P_k = 0 \quad (16)$$

$W_{1k}$ and $W_{2k}$ are the normal displacements of bodies 1 and 2 at point *k*, respectively, $P_k$ is the contact pressure at point *k*, $\epsilon_k$ is the normal initial deviation at point *k* (computed in the unloaded TCA), and $\alpha$ is the global adjustment between the gear bodies (transmission error expressed on the line of action).

The normal displacement of body *i* at point *k* can be written as

$$W_{ik} = \sum_{j=1}^{N} C_{kj}^{i} P_{j} \quad (17)$$

For two bodies, $C_{ik} = C_{kj}^{1} + C_{kj}^{2}$ and Eq. (17) yields

$$\sum_{i=1}^{N} C_{ik} P_{k} + \epsilon_{k} = \alpha, \quad (18)$$

$$\sum_{k=1}^{N} P_k S_k = F \quad (19)$$

inside the contact zone. The fixed-point algorithm is used to solve the system of equations. Coefficients $C_{ik}$ are the influence coefficients obtained with a uniform pressure $P_k$ applied on a small surface $S_k$ around point *k* of the contact zone of body *i*.

<figure><img src="image_11.png" /><figcaption>Fig. 9 Local meshing</figcaption></figure>

<table><caption>Table 1 Parameters of the example helical gear pair</caption><thead><tr><th></th><th>Pinion</th><th>Gear</th></tr></thead><tbody><tr><td>Number of the teeth</td><td>20</td><td>31</td></tr><tr><td>Module (mm)</td><td>6</td><td>6</td></tr><tr><td>Pressure angle (deg)</td><td>20</td><td>20</td></tr><tr><td>Helix angle (deg)</td><td>15</td><td>15</td></tr><tr><td>Face width (mm)</td><td>70</td><td>70</td></tr></tbody></table>

<figure><img src="image_12.png" /><figcaption>Fig. 10 Potential contact lines</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Fig. 11 Tooth load sharing</figcaption></figure>

Two terms define $C_{ik}$, $C_{ik} = C_{ik}^{f} + C_{ik}^{s}$. Term with superscripts *f* and *s* represents the bulk deformations and the surface deformations, respectively. $C_{ik}^{f}$ is calculated using the FPM method and $C_{ik}^{s}$ is found by using the well-known Boussinesq theory [22]. In order to implement these calculations, a local square mesh is produced in the tangent plane of the contact zone as shown in Fig. 9. The influence coefficients are obtained by applying unit loads at each node that are likely to be loaded. Two types of result are obtained: (i) load sharing and (ii) contact pressure distribution. Load sharing refers to the sharing of the torque among meshing tooth pairs. The contact pressure distribution is found along the face width for the different teeth in contact.

<figure><img src="image_14.png" /><figcaption>Fig. 12 Contact pressure for the kinematics position 1</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>Fig. 13 Contact pressure during meshing on one helical tooth</figcaption></figure>

**3.3 Model Results.** In this section, results for an example helical gear pair are presented. The design parameters of the helical gear pair as listed in Table 1. The input torque is 2000 Nm and the involute gear tooth profiles are not modified. Figure 10 presents potential contact lines on the teeth of the driving pinion for three different kinematic positions. For position 1, only teeth 2 and 3 are in contact. For an intermediate position (position 5), three pairs of teeth are in contact (teeth 2, 3, and 4). For tooth 3, the contact line on tooth 3 extent is along the entire tooth width. Tooth 2 experiences only a contact at the tip of tooth. At the last kinematic position, tooth 2 is not in contact while tooth number 5 enters the mesh.

The load sharing characteristics are presented in Fig. 11 as a result of computation of the elastic behavior. Tooth 1 is not loaded for this range of kinematic positions. Tooth 4 enters in the contact at position 2 and tooth 5 at position 8 while tooth 2 leaves the contact at position 6. The contact pressure distribution along the flanks of teeth 2 and 3 is presented in Fig. 12 for the kinematic position 1. Peaks of pressure on the edge of the teeth can be noted. Similarly, Fig. 13 presents the field of contact pressure on the square grid defined around the contact zone for several kinematic positions, from the beginning of a tooth meshing, to the exit of the tooth from the contact zone. The transmission error and the gear mesh stiffness are also predicted. An example of stiffness is shown in Fig. 14 as a function of the angular position of the driving pinion. The stiffness is calculated with the ratio of the input load to the transmission error ($\alpha$). The stresses in the tooth roots are computed by using the distribution of the load. Figure 15 presents the maximal principal stresses in four consecutive roots at kinematic position 5. Roots 2 and 3, located between loaded teeth experience compressive stresses.

<figure><img src="image_16.png" /><figcaption>Fig. 14 Meshing stiffness variation</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 15 Principal stresses in tooth roots (kinematic position 5)</figcaption></figure>

## 4 Experimental Validation

Any numerical tool must be validated experimentally before it can be used in confidence as a design and analysis tool. This section aims at validating the gear analysis model presented in earlier sections experimentally.

Tests were performed on an existing back-to-back type fatigue test machine shown in Fig. 16. Table 2 specifies the parameters of the helical gears used in the experimental test. The gear manufacturing quality is ISO4 and tooth profiles are not modified. The schematic layout of the machine is shown in Fig. 17. The shafts

<figure><img src="image_18.png" /></figure>

Nominal speed: 1610 tr/mn
Allowed torque: 200 to 600 daNm
Maximal dissipative power: 100KW

<caption>Fig. 16 Picture of the test machine</caption>

<table><caption>Table 2 Geometrical characteristics of the gears</caption><thead><tr><th></th><th>Pinion<br>(driving and driven)</th><th>Gear<br>(driving and driven)</th></tr></thead><tbody><tr><td>Number of teeth</td><td>159</td><td>32</td></tr><tr><td>Working module (mm)</td><td>3.56</td><td>3.56</td></tr><tr><td>Width (mm)</td><td>86</td><td>80</td></tr><tr><td>Pitch diameter (mm)</td><td>566.039</td><td>113.92</td></tr><tr><td>Helix angle (deg)</td><td>19.967</td><td>19.967</td></tr><tr><td>Helix hand</td><td>Left</td><td>Right</td></tr><tr><td>Working pressure angle (deg)</td><td>25.347</td><td>25.347</td></tr><tr><td>Tip diameter (mm)</td><td>572.1</td><td>121.2</td></tr><tr><td>Root diameter (mm)</td><td>556.1</td><td>105.1</td></tr></tbody></table>

<figure><img src="image_19.png" /><figcaption>Fig. 17 Layout of the test machine</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>Fig. 18 Location of the gauges</figcaption></figure>

were mounted on high-capacity tapered roller bearings. A rotating hydraulic torque device applies different torque values on the inner test loop and the system is driven by an electric motor.

During the tests, the angular rotation of the pinions is kept very low (5 rpm) so that a comparison with the quasi-static model is possible with no dynamic effects. The gears were instrumented with strain gages. A total of six gauges were mounted on each gear, as shown in Fig. 18 at a 120 deg tooth fillet position and two locations along the tooth fillet width. The gauges are located near the peak bending stress location. An optical encoder was also used to mark the angular position as gears rotate. The tests were performed for one rotation of the larger gears that corresponded to

<figure><img src="image_21.png" /><figcaption>Fig. 19 Meshing and boundary conditions</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>Fig. 20 Stresses for the driven pinion, face width position 1</figcaption></figure>

nearly five rotations of the pinion, so that the strain-stress variations could be monitored for each strain gage. When the strain gauged tooth fillet entered the meshing zone, the signal varied from positive to negative values (see Fig. 20). The results given for the circular positions 1, 1' or 1" are also plotted in Fig. 20.

The tested gear was modeled numerically. Figure 19 shows the FPM of the test gears and the boundary conditions used. The rim thickness is greater than 3.5 times the module, a five-tooth segment is sufficient to perform a reasonably good numerical simulation. The input torque value was set at 600 N·m. In Fig. 20, example predictions and the corresponding stresses during meshing are compared. The measurements for the three gauges situated at an angular distance of 120 deg agree very well with the simulation results (solid line in Fig. 20). Differences between the numerical and measured results can be explained by secondary effects such as pitch error, shaft misalignments, and helix angle errors [23].

## 5 Conclusions

The present paper describes a numerical tool that aids the design of cylindrical gears. The elastic behavior of gear contact was studied with the finite prism method. This method has been shown to correlate to the results of conventional FEM well. It was shown that the numbers of terms in the series functions influence of the accuracy of the predictions. The integration of this method to calculate the gear bending coefficients and the application of Boussinesq theory allows to computation of the load sharing and the tooth root stresses. The FPM gives results in three dimensions with significantly less computational demand. An experimental test results were also shown to match the predictions of this model well.

## Acknowledgments

The initial research program was financed by the EUROCOP-TER FRANCE Company mainly during the Ph.D of A. A. Olakorede. Additional works were performed in the framework of the Ph.D of H. C. Kim and J. P. de Vaujany. The authors thank D. Remond and J. Mahdoudh for their help during the experiments.

## References

[1] Von Eiff, H., Hirschmann, K. H., and Lechner, G., 1990, "Influence of Gear Tooth Geometry on Tooth Stress of External and Internal Gears," ASME J. Mech. Des., **112**, pp. 575-583.

[2] Shuting, L., 2002, "Deformation and Bending Stress Analysis of a Three-Dimensional, Thin-Rimmed Gear," ASME J. Mech. Des., **124**, pp. 129-135.

[3] Conry, T. F., and Seireg, A., 1971, "A Mathematical Programming Method for Design of Elastic Bodies in Contact," J. Appl. Mech., **38**, pp. 387-392.

[4] Wagaj, P., and Kahraman, A., 2002, "Influence of Tooth Profile Modification on Helical Gear Durability," ASME J. Mech. Des., **124**, pp. 501-510.

[5] Olakorede, A. A., and Play, D., 1991, "Load Sharing, Load Distribution and Stress Analysis of Cylindrical Gears by Finite Prism Method," *Design Productivity International Conference, Honolulu, Hawaii*, pp. 921-927.

[6] Stegemiller, M. E., and Houser, D. R., 1993, "A Three-Dimensional Analysis of the Base Flexibility of Gear Teeth," ASME J. Mech. Des., **115**, pp. 186-192.

[7] Vijayakar, S. M., Rusby, H. R., and Houser, D. R., 1987, "Finite Element Analysis of Quasi-Prismatic Bodies Using Chebyshev Polynomials," Int. J. Numer. Methods Eng., **24**, pp. 1461-1477.

[8] Elkholy, A. H., 1989, "Load Distribution on Contact Lines of Helical Gear Teeth," *Proc. International Power Transmission and Gearing Conference, ASME, 1*, pp. 135-142.

[9] Vijayakar, S. M., and Houser, D. R., 1993, "Contact Analysis of Gears Using a Combined Finite Element and Surface Integral Method," Gear Technology, pp. 26-33.

[10] Wilcox, L., and Coleman, W., 1973, "Application of Finite Elements to the Analysis of Gear Tooth Stresses," ASME J. Eng. Ind., **95**, pp. 1139-1148.

[11] Frater, J. L., and Kasuba, R., 1981, "Extended Load-Stress Analysis of Spur Gearing: Bending Strength Considerations," *International Symposium on Gearing & Power Transmissions*. Tokyo, pp. 147-152.

[12] Sun, H., Mavriplis, D., Huston, R. L., and Oswald, F. B., 1989, "Comparison of Boundary Element and Finite Element Method in Spur Gear Root Stress Analysis," *International Power Transmissions and Gearing Conference, ASME, 1*, pp. 163-166.

[13] Baud, S., and Velex, P., 2002, "Static and Dynamic Tooth Loading in Spur and
Helical Geared Systems-Experiments and Model Validation," ASME J. Mech. Des., **124**, pp. 334–346.

[14] Kim, H. C., de Vaujany, J. P., Guingand, M., Bard, C., and et Play, D., 1995, "Effects of Rim, Web and Constraint Conditions on Stresses of External and Internal Cylindrical Gears," *International Congress-Gear Transmission'95, Sofia, Bugaria*, Vol. 2, pp. 164–171.

[15] de Vaujany, J. P., Kim, H. C., Guingand, M., and et Play, D., 1996, "Effect of Rim and Web on the Stresses of Internal Cylindrical Gears," *7th International Power Transmission and Gearing Conference, ASME International, San Diego CA*.

[16] Olakorede, A. A., and Play, D., 1988, "Development of the Finite Prism Method in Computer: Aided Design of Spur Gears," *Proc. of the fourth SAS, World Conf. FEMCAD, Paris, 17–19, Oct.* 1988, Gournay-Sur-Marne (Fr.): IITT-International., Vol. 1, pp. 384–391.

[17] Kantorovich, L. V., and Krylov, V. I., 1964, *Approximate Method of Higher Analysis*, Interscience Publishers New York, p. 681.

[18] Cheung, Y. K., 1976, *Finite Strip Method in Structure Analysis*, Pergamon, Oxford, p. 233.

[19] Olakorede, A., 1990, "Répartition de Charges et Résistance en Conception d'Engrenages Cylindriques, Application de la Méthode des Prismes Finis en C A O," Thèse de Doctorat. Institut National des Sciences Appliquées de Lyon, p. 157.

[20] Bishop, R. E. D., and Johnson, D. C., 1956, *Vibration Analysis Tables*, Cambridge University Press, London, p. 59.

[21] Meirovitch, L., 1986, *Elements of Vibration Analysis*, 2nd ed., McGraw-Hill, New York, p. 560.

[22] Boussinesq, J., 1885, *Application des Potentiels à L'étude de L'équilibre et du Mouvement des Solides Élastiques*, Gauthiers-Villard, Paris.

[23] Kim, H.-C., 1996, "Distribution des Charges et des Contraintes Dans les Engrenages Cylindriques Avec Voile," *Thèse Doct. Institut National des Sciences Appliquées de Lyon*, p. 151.