Non-Newtonian mixed thermo-elastohydrodynamics of hypoid gear pairs

Abstract

Transmission efficiency is the main objective in the development of vehicular differential systems, comprising hypoid gear pairs. The overall aim is to contribute to improved vehicle fuel efficiency and thus levels of harmful emissions for modern desired eco-drive axles. Detailed predictive analysis plays an important role in this quest, particularly under realistic operating conditions, comprising high contact loads and shear rates. Under these conditions, the hypoid gear pairs are subject to mixed non-Newtonian thermo-elastohydrodynamic conditions, which is the approach undertaken in this paper. Such an approach for hypoid gear pair has not hitherto been reported in the literature.

Keywords

Differential hypoid gears non-Newtonian shear mixed thermo-elastohydrodynamics friction efficiency

Introduction

The key objectives in modern vehicle development are fuel efficiency and reduced emissions. The former has been the cornerstone of the lightweight powertrain in addition to improved combustion and transmission efficiency. These have led to vehicles of higher power-to-weight ratio, but at the expense of a plethora of noise and vibration concerns, generally referred to as NVH (noise, vibration and harshness).¹ Therefore, the drive for fuel-efficient lightweight powertrains has resulted in some major NVH concerns, for instance driveline clonk^1–4 and transmission rattle.^1,5–7 With the lightweight concept remaining the over-riding objective and increased power a desired customer attribute, the excess energy leading to impulsive action needs to be somehow palliated. Under these circumstances, friction remains an important energy sink, which in turn would adversely affect the system efficiency. Therefore, the fundamental understanding of regimes of lubrication in powertrain conjunctions under realistic operating conditions assumes an important role. This is regarded as the prerequisite in powertrain optimisation, taking into account both fuel efficiency and NVH refinement.

Some friction is essential for effective traction, transmitting the engine power to the driven wheels, whilst at the same time consuming some of the excess vibration energy, causing NVH response.⁸ Reduced friction often leads to repetitive impacts of gear teeth pairs or contact separation, thus loss of traction.^9,10 Lubricant is subjected to changes in the regime of lubrication as shown by De la Cruz et al.¹¹

In automotive transmissions and differential gearing, thin films promote mixed thermo-elastohydrodynamic regime of lubrication. Therefore, friction is generated as the result of viscous shear of lubricant as well as counterface asperity interactions. In transmission systems, lubricant shear may be assumed to be Newtonian, except for the thin adsorbed films at the summits of asperities.¹¹ However, for highly loaded contacts of hypoid gears, a non-Newtonian regime of lubrication is prevalent. Hence, a constant coefficient of friction cannot be assumed as in some reported cases,^12,13 nor the use of empirical formulae.^14–17 Another approach is proposed by He et al.,¹⁸ Velex et al.^19,20 and Kar and Mohanty,²¹ describing friction as a function of the length of the line of contact for the meshing of spur or helical gear teeth pairs.

It is more accurate to represent the meshing of gear problem as lubricated conjunctions as highlighted by De la Cruz et al.¹¹ for the case of helical gears. An analytical approach is highlighted by Karagiannis et al.,²² using the Grubin’s extrapolated oil film thickness equation as a function of contact load and sliding velocity.²³ In Karagiannis et al.,²² thin thermo-elastohydrodynamic films under non-Newtonian shear promoted asperity friction. De la Cruz et al.¹¹ showed that predicted film thickness, using the Grubin’s equation can be as much as 50% larger than that under the same operating condition when a numerical elastohydrodynamic analysis is used. Therefore, ideally a numerical analysis for hypoid gear pairs should be undertaken.

For helical gears, the contact footprint shape is elliptical with a large aspect ratio. Thus, line contact geometry is often assumed.²⁴ For the case of hypoid gear pairs, an elliptical contact footprint with angled lubricant entrainment into the contact results. This is demonstrated by Mohammadpour et al.^25,26 who provided isothermal elastohydrodynamic solutions for hypoid gear teeth pairs throughout a meshing cycle both under assumed Newtonian and non-Newtonian lubricant behaviour.

Contact friction results in generated heat, which among other issues affects the lubricant viscosity, and hence the load-carrying capacity and friction. Crook²⁷ was the first to present a detailed method for determination of generated heat and its removal in line contacts. In particular, he predicted a fall in traction as the slide–roll ratio increases. Subsequently, Cheng and Sternlicht²⁸ and Dowson and Whitaker²⁹ introduced numerical procedures for the coupled solution of Reynolds, elasticity and energy equations for line contacts.

Kim and Sadeghi^30,31 were among the first to present a full numerical solution for thermal elastohydrodynamic lubrication (TEHL) of point contacts for pure rolling under low-slip conditions. Yang et al.³² studied the general transient TEHL problem, including the first full numerical solution of the energy equation using finite differences. However, the inlet flow entrainment was assumed to be either along the minor or the major axis of the elliptical contact footprint. Their results showed the insignificant influence of rising temperature on the lubricant film thickness, but a noteworthy effect upon coefficients of traction and friction. This work was extended to the non-Newtonian investigation by Yang et al.³³ None of these contributions took into account the angled inlet flow which gives rise to side leakage of the lubricant from the contact, thus affecting the heat balance. Such conditions are prevalent in the case of hypoid gears, which is the subject of the current study. Sharif et al.³⁴ presented a thermal non-Newtonian model to obtain film thickness and traction in a variable traction drive rig. The model was presented for circular point contacts only. Habchi et al.³⁵ also presented a non-Newtonian thermal model for circular point contacts. The model did not take into account the effect of asperity interactions. The same approach was also used by Paouris et al.³⁶ for circular point contact of a sphere against a very smooth disc using a mini-traction machine (MTM). Inlet lubricant starvation and shear heating were included in the model, which agreed well with the experimental measurements of generated friction.

For hypoid gears, Kolivand et al.³⁷ provided an equation for friction of mixed elastohydrodynamic with non-Newtonian lubricant shear. However, a line contact geometry was assumed. Simon³⁸ presented a thermo-elastohydrodynamic analysis of hypoid gear pairs with elliptical point contact geometry, but with the entrainment flow vector was along the minor axis of the contact ellipse. Furthermore, a low applied input torque, not representative of vehicular differentials, was used. As already noted, the entrainment flow into the contact occurs at an angle to the minor axis of the contact ellipse. This is shown by Gohar³⁹ through optical interferometry with significant side leakage flow at an angle to the major axis of the contact ellipse.⁴⁰ Pu et al.⁴¹ further developed the work of Mohammadpour et al.²⁶ by adding an empirical flash temperature model to the analysis. Their analytical thermal model is a control volume-type model for the average temperature of the lubricant and those of the contacting surfaces. Clearly, the temperature of the lubricant alters in the contact at any point into its depth, in the direction of sliding as well as in the side leakage direction. Therefore, a more realistic model is required, using the solution of energy equation, leading to accurate estimation of lubricant viscosity and determination of film thickness. Furthermore, Pu et al.⁴¹ assume a constant coefficient of friction within the contact area, which also varies from point-to-point in the contact. Therefore, an integrated analysis using Reynolds and energy equation and heat transfer from the solid surfaces is essential, which is the approach of the current paper.

This paper extends the non-Newtonian mixed elastohydrodynamic analysis of hypoid gears in a vehicular differential system at realistic high loads through a combined solution of Reynolds, elasticity and energy equations. It also takes into account the inlet starvation effect as well as angled lubricant flow into the contact. Furthermore, tooth contact analysis (TCA)⁴² is used to determine the contact footprint geometry during the meshing cycle of pinion and ring gear as well as surface velocities of the contiguous solids. This approach has not hitherto been reported in literature.

The hypoid gear pair model

The hypoid gears’ model comprises simultaneous contacts of a number of teeth pairs. For the automotive differential investigated, 2–3 pairs of teeth are in simultaneous mesh. Therefore, it is necessary to calculate contact geometry, kinematics and load share for each meshing pair through tooth contact analysis (TCA).⁴² A full description of the approach is provided in Mohammadpour et al.²⁵

The meshing teeth pairs are subject to mixed-elastohydrodynamic regime of lubrication and subjected to quite high loads, and the lubricant undergoes non-Newtonian shear. Therefore, a suitable non-Newtonian model for lubricant shear behaviour should be used, as well as a suitable form of Reynolds equation.

The non-Newtonian elastohydrodynamic conjunction

Conry et al.⁴³ developed the Reynolds–Eyring equation, based on Johnson and Tevaarwerk’s shear model,⁴⁴ neglecting lubricant elastic shear for large strain rates (i.e. $D < 1$ ). Conry et al.⁴³ studied the elastohydrodynamics line contact problem, ignoring the side leakage of the lubricant from the conjunction. However, for hypoid gear teeth pairs, a two-dimensional solution is essential due to the angled flow entrainment into the contact (Figure 1).

Figure 1.

Representation of an elliptical point contact conjunction with angled entrainment flow.

Extending the approach of Conry et al.⁴³ to the case of elliptical point contacts of a hypoid gear teeth pair yields

\frac{\partial}{\partial x} [\frac{ρ h 3}{η} s (x) \frac{\partial p}{\partial x}] + \frac{\partial}{\partial y} [\frac{ρ h 3}{η} s (y) \frac{\partial p}{\partial y}] = 12 U \times {\cos θ \frac{\partial}{\partial x} [ρ h] + \sin θ \frac{\partial}{\partial y} [ρ h]}

(1)

where θ is the angle of flow entrainment, measured from the minor axis of the elliptical contact footprint (Figure 1) and U is the speed of entraining motion, where

U = \frac{1}{2} (U_{1} + U_{2})

. The form of equation (1) omits the squeeze film effect, which occurs as the result of approach and separation of gear teeth pairs. This transient effect often increases the load-carrying capacity of the contact, as noted by Gohar and Rahnejat.⁴⁵ s(x) and s(y) are dimensionless modifying parameters for lubricant viscosity, complying with the non-Newtonian hyperbolic sine model in Conry et al.⁴³

s (x) = {1 + (\frac{η Δ U \cos θ}{τ_{0} h} \frac{ξ_{x}}{\sinh ξ_{x}}) 2} 1 / 2 \times {\frac{3 (ξ_{x} \cosh ξ_{x} - \sinh ξ_{x})}{ξ_{x}^{3}}}

(2)

where

ξ_{x} = \frac{h}{2 τ_{0}} \frac{\partial p}{\partial x} and ξ_{y} = \frac{h}{2 τ_{0}} \frac{\partial p}{\partial y}

(4)

Theeffectiveviscosityis η_{eff} = \frac{η}{s (x)} or η_{eff} = \frac{η}{s (y)}

(5)

And η is the average viscosity in the direction of the film thickness updated by pressure and temperature iteratively⁴⁶

η = η_{0} \exp {(\ln η_{0} + 9.67) \times [- 1 + (1 + 5.1 \times 10 - 9 p) z_{0} (\frac{T - 138}{T_{0} - 138}) - s_{0}]}

(6)

where

z_{0} = \frac{\propto}{5.1 \times 10 - 9 (\ln η_{0} + 9.67)} and s_{0} = \frac{\bar{β} (T_{0} - 138)}{(\ln η_{0} + 9.67)}

The above viscosity relationships do not adequately deal with high pressure high shear lubricant characteristics. Recent rheological models are more suited for this purpose as highlighted by Habchi et al.³⁵ and Paouris et al.³⁶ However, the focus of this paper is on mixed non-Newtonian conditions.

Lubricant density variation with pressure and temperature is given as⁴⁷

ρ = ρ_{0} [1 + \frac{0.6 \times 10 - 9 p}{1 + 1.7 \times 10 - 9 p}] \times [1 - 0.65 \times 10 - 3 (T - T_{0})]

(7)

The speed of entraining motion of the lubricant U and the relative sliding velocity of the surfaces are obtained through TCA⁴²

The limiting shear stress is

τ_{L} = τ_{L 0} + λ' p

(8)

The combined non-Newtonian shear model is used, where shear stress is continually monitored during the meshing cycle of any pair of teeth and at any computation point in the conjunction, using equation (9). If its value exceeds that obtained from equation (8), then it is replaced by the limiting shear stress at the pressure at any given computation point

τ = η_{eff} {\overset{\cdot}{γ}}_{v} = \frac{η_{eff} Δ U}{h}

(9)

It should be noted that the pressure at any grid is assumed to represent the behaviour on the grid area. Equation (9) is the rheological model for the lubricant when the shear stress reaches its limiting value locally according to Johnson and Tevaarwerk’s shear model⁴⁴

The elastic film shape is given as

h (x, y) = h_{c 0} + s (x, y) + δ (x, y)

(10)

where the undeformed profile is:

s (x, y) = \frac{x 2}{2 R_{zx}} + \frac{y 2}{2 R_{zy}}

$R_{zx} and R_{zy}$ are the principal radii of contact of an equivalent ellipsoidal solid against a semi-infinite elastic half-space, representing the instantaneous contact of any pinion to ring gear teeth pair in the planes zx (along the minor axis of the contact ellipse) and zy (along its major axis, Figure 2)

\frac{1}{R_{zx}} = \frac{1}{r_{p}} + \frac{1}{r_{w}} and \frac{1}{R_{zy}} = \frac{1}{R_{p}} - \frac{1}{R_{w}}

(11)

Figure 2.

Equivalent contact geometry of the equivalent ellipsoidal solid.

The instantaneous radii of pinion and gear teeth are determined through TCA⁴²

The localised contact deflection $δ (x, y)$ is obtained through solution of the elasticity potential integral

δ (x, y) = \frac{1}{E_{r}} \int \int_{A} \frac{p (x_{1}, y_{1}) d x_{1} d y_{1}}{\sqrt{(x - x_{1}) 2 + (y - y_{1}) 2}}

(12)

where (x,y) represents a point where deflection of the semi-infinite elastic half-space of reduced elastic modulus E_r is determined due to any arbitrary pressure distribution

p (x_{1}, y_{1})

Boundary and initial conditions for Reynolds equation

Commonly, a fully flooded inlet is assumed in most numerical analysis of EHL problems. In order to ensure this, the following conditions are set: $x \to - \infty, p \to 0$ and $y \to - \infty, p \to 0$ . This means that the inlet boundary should be set at a distance of 4–5 times the contact footprint semi-half-widths.^40,48 However, a fully flooded inlet is idealised and in practice only a partially flooded inlet occurs.^49,50 This is on account of swirl and reverse flows prior to the flow stagnation boundary.^50,51 Therefore, the inlet distance changes with the relative motion of the contiguous solid surfaces and the applied load: $m = x / b or y / a$ . Hamrock and Dowson⁵² carried out numerical analysis for starved inlets and determined the demarcation starvation boundary for just flooded condition to be

m * = 1 + 3.06 [(R / γ) 2 h_{c 0}] 0.58

(13)

where for the flow component along the minor axis,

R = R_{zx}, γ = b

, and for that along the major axis,

R = R_{zy}, γ = a

The inlet condition is established as $m = m *$ , which use the Hamrock and Dowson’s starvation boundary.

The initial estimation of the central contact oil film thickness is due to Grubin²³ with the inclusion of a side-leakage correction factor⁵³

h_{c 0} / R_{zx} = 1.212 (α η_{0} U / R_{zx}) 3 / 4 (W / E_{r} R_{zx}^{2}) - 1 / 12 \times {1 + \frac{2}{3} (R_{zy} / R_{zx})} - 2 / 3

(14)

The contact load is

W = \int_{x = x_{i}}^{x_{c}} \int_{y = y_{i}}^{y_{c}} p d x d y

(15)

The inlet boundaries given by x_i and y_i are set, using Hamrock and Dowson’s inlet boundary location.⁵² The outlet boundary conditions used are those of Swift⁵⁴ and Stieber⁵⁵: $p = \partial p / \partial x = 0 at x = x_{c} and p = \partial p / \partial y = 0 at y = y_{c}$ , where the film rupture positions $x_{c} and y_{c}$ are obtained by discarding any negatively generated pressures in the iterative numerical procedure. The outlet boundary conditions are at the exit of the computational domain in the directions of entraining motion and side leakage at fixed distances of 1.5 times the contact footprint semi-half-widths.⁵²

Thermo-elastohydrodynamics

The thermal elastohydrodynamic analysis comprises simultaneous solution of Reynolds (see The non-Newtonian elastohydrodynamic conjunction section) and the energy equations.

The energy equation for the lubricant film, neglecting heat conduction in the x- and y-directions due to the thinness of the lubricant film and convection in the z-direction (no convection by the lubricant film into the adjacent solid boundaries) becomes

ρ c (u_{f} \frac{\partial T}{\partial x} + v_{f} \frac{\partial T}{\partial y}) = k_{f} \frac{\partial 2 T}{\partial z 2} - \frac{T}{ρ} \frac{\partial ρ}{\partial T} (u_{f} \frac{\partial p}{\partial x} + v_{f} \frac{\partial p}{\partial y}) + η ((\frac{\partial u_{f}}{\partial z}) 2 + (\frac{\partial v_{f}}{\partial z}) 2)

(16)

The velocities and velocity gradients in the fluid film are given as

u_{f} = \frac{1}{2 η} \frac{\partial p}{\partial x} (z 2 - zh) + \frac{z}{h} (U_{1} - U_{2}) + U_{2}

(17)

\frac{\partial u_{f}}{\partial z} = \frac{1}{2 η} \frac{\partial p}{\partial x} (2 z - h) + \frac{(U_{1} - U_{2})}{h} + U_{2}

(18)

v_{f} = \frac{1}{2 η} \frac{\partial p}{\partial y} (z 2 - zh) + \frac{z}{h} (V_{1} - V_{2}) + V_{2}

(19)

\frac{\partial v_{f}}{\partial z} = \frac{1}{2 η} \frac{\partial p}{\partial y} (2 z - h) + \frac{(V_{1} - V_{2})}{h} + V_{2}

(20)

These velocity distributions are based on the Newtonian flow assumption. This assumption somewhat simplifies the problem and can be justified by the presented results in Figure 12. The figure shows that severe thermal thinning forces acting on the lubricant closely follow Newtonian shear. Velocity distributions, based on Newtonian and non-Newtonian assumptions, are shown in the Results and discussion section. The viscosity used in equations (16) to (20) is the average viscosity in the direction of film thickness, updated by the pressure and temperature iteratively. It is obtained from equation (6). It is noteworthy that effective viscosity is given by equations (16) to (20) if its local value reaches the limiting shear stress.

The boundary condition for the lubricant at the inlet conjunction is assumed to be that of the bulk sump temperature as

T (x_{in}, y, z) = T_{0}

(21)

The contacting solids are simplified in the form of an ellipsoidal solid and a semi-infinite plane. The energy equations for these conjunctional solid boundaries are

c_{p} ρ_{p} U_{2} \frac{\partial T}{\partial x} + c_{p} ρ_{p} V_{2} \frac{\partial T}{\partial y} = k_{p} (\frac{\partial 2 T}{\partial x 2} + \frac{\partial 2 T}{\partial y 2} + \frac{\partial 2 T}{\partial z 2}) c_{s} ρ_{s} U_{1} \frac{\partial T}{\partial x} + c_{s} ρ_{s} V_{1} \frac{\partial T}{\partial y} = k_{s} (\frac{\partial 2 T}{\partial x 2} + \frac{\partial 2 T}{\partial y 2} + \frac{\partial 2 T}{\partial z 2})

(22)

The index s denotes the ellipsoidal solid with the surface speed of the pinion tooth. Similarly, index p denotes the semi-infinite plane which is considered to be moved by the gear surface velocity. The boundary conditions for equation (27) are

{T (x_{in}, y, z) = T_{0} I T (x, y, + d) = T_{0} II T (x, y, - d) = T_{0} III

(23)

The first boundary condition (I) assumes a constant bulk housing temperature of the lubricant. This assumption is based on a thermal balance having been reached for the differential unit. The remaining boundary conditions II and III assume temperature penetration depths of d and –d into the solid boundaries. These are ﬁnite depths, sufficiently large to constitute a zero temperature gradient into the solid boundaries. The value of d in the current analysis is taken as three times the Hertzian contact semi-half-width in the direction of entraining motion.^32,33

The thermal conductivity and heat capacity of the lubricant has been assumed to be constant. In reality, both of these quantities change with the pressure and temperature. The effect of these properties on the film thickness and friction of a circular point contact is presented by Habchi et al.⁵⁶ They showed that the variable thermal conductivity and heat capacity could change the friction as the main tribological quantity up to 8%.

To satisfy the continuity of heat flux between the solid surfaces and the lubricant film, the following interfacial boundary conditions should be satisfied

k \frac{\partial T}{\partial z} {z = 0 + = k_{p} \frac{\partial T}{\partial z} {z = 0 - k \frac{\partial T}{\partial z} {z = h - = k_{s} \frac{\partial T}{\partial z} {z = h +

(24)

Figure 3 shows the thermal boundary conditions and the employed coordinate system.

Figure 3.

Boundary conditions and the coordinate system used in the thermal model.

Prediction of friction

Hypoid gears usually operate with quite thin lubricant films under partially lubricated non-Newtonian conditions. Thus, the generated friction is the combined result of non-Newtonian viscous shear of the lubricant film as well as some degree of boundary friction due to interaction of asperities of counter face surfaces. Therefore, the prevailing conditions can be termed as mixed-thermo-elastohydrodynamic. The extent of direct solid boundary interaction is determined through calculation of oil film parameter, λ at any instant through the meshing cycle for any pair of contacting teeth, where

λ = \frac{h}{σ}

(25)

where the composite root mean square roughness of the counter face surfaces is

σ = \sqrt{σ_{1}^{2} + σ_{2}^{2}}

(26)

The surface of the hypoid gear pair, subject of the current study, is lapped and super-finished, and thus fairly smooth,^37,57 where boundary interactions are expected to follow the method proposed by Greenwood and Tripp.⁵⁸ This method assumes a Gaussian distribution of asperity heights, which is not always the case in practice. However, this is the case considered here (see Appendix 1).

The real area of interacting asperity peaks,A_a, is a function of bulk separation and is described by the statistical function $F_{2} (λ)$ in terms of Stribeck’s oil film parameter⁵⁸

A_{a} = π 2 (ξ β σ) 2 {AF}_{2} (λ)

(27)

The proportion of load carried by the asperities is obtained as

F_{a} = \frac{8 \sqrt{2}}{15} π (ξ β σ) 2 \sqrt{\frac{σ}{β}} E * {AF}_{5 / 2} (λ)

(28)

In these, formulae A is the Hertzian contact area (apparent area of contact for perfectly smooth counter face surfaces),ξ is the density of asperity peaks per unit area and β is the average radius of curvature of asperity peaks. In addition, $E * = \frac{2 E_{r}}{π}$ and σ is the root mean square roughness of the counter faces (equation (26)). $F_{5 / 2} (λ)$ is another representative statistical function of surface topography⁵⁸

According to Teodorescu et al.⁵⁹

F_{5 / 2} = - 0.1922 λ 3 + 0.721 λ 2 - 1.0649 λ + 0.6163

(29)

F_{2} = - 0.116 λ 3 + 0.4862 λ 2 - 0.7949 λ + 0.4999

(30)

According to Greenwood and Tripp,⁵⁸ the roughness parameter $ξ β σ$ is usually in the range 0.03–0.05; the ratio $σ / β$ is a measure of asperity slope and is usually in the range 10⁻⁴–10⁻². These parameters were measured using a white light interferometer (Alicona with the vertical/height resolution of 10 nm and lateral resolution of 0.175 µm) to be $ξ β σ$ = 0.055; therefore, in this case

A_{a} = 0.0298 {AF}_{2} (λ)

(31)

Also, using the measured topography: $σ / β$ = 0.001, and the load carried by the asperities becomes

F_{a} = 0.000227 E * {AF}_{5 / 2} (λ)

(32)

The conjunctional friction is the result of combined viscous shear of a thin lubricant film, F_v and direct interaction of asperities on the counterface surfaces (boundary friction, F_b), thus

F_{f} = F_{b} + F_{v}

(33)

The viscous friction becomes

F_{v} = \int \int | \pm \frac{h}{2} \overset{⇀}{\nabla} p - Δ \overset{⇀}{V} \frac{η}{h} | d x d y

(34)

where

\overset{⇀}{\nabla}

is the Laplacian operator and

Δ \overset{⇀}{V}

is the vector of relative sliding velocity of contiguous surfaces in the x and y directions.

The viscosity is the average viscosity in the direction of film thickness, updated by the pressure and temperature using equation (6). It should be noted that viscosity assumes its effective value as the local viscosity promotes limiting shear stress. The double integral is calculated over the apparent contact area.

Boundary friction is as the result of direct interaction of a very thin film adsorbed to the summit of counterface asperities, which is subject to non-Newtonian shear. It is assumed that limiting shear stress at the asperity pressure governs these interactions, where

F_{b} = τ_{b} A_{a}

(35)

where boundary shear is obtained as

τ_{b} = τ_{L 0} + λ' P_{a}

(36)

where P_a is the asperity pressure

P_{a} = \frac{F_{a}}{A_{a}}

(37)

Method of solution

Coupling between Reynolds and energy equations is not explicit as temperature affects Reynolds equation indirectly through changes in the lubricant rheological state (density and viscosity). In the conjunctional inlet, flow and thermal effects are strongly coupled due to the viscous nature of the flow. On the other hand, in the high pressure region of the contact, viscosity increases significantly with pressure, whilst the lubricant film thickness decreases. Therefore, the Poiseuille component of flow becomes negligible and the convective heat transfer is no longer significant. This is also shown by Kim et al.^60,61

Reynolds equation is discretised using finite differences in the same manner as that described by Jalali-Vahid et al.⁴⁰ Solution for pressure at any nodal position (i,j) within the computational grid is based on low relaxation effective influence Newton–Raphson method, fully described by Jalali-Vahid et al.^40,62

A full numerical solution using finite difference form of the energy equations is required to calculate the temperature distribution in the fluid film. Lubricant film temperature alters three-dimensionally, unlike the pressure distribution (where: $\frac{\partial p}{\partial z} = 0$ ). Thus, the temperature distribution is calculated as $f (x, y, z)$ .

The following iterative procedure is used:

At any pinion angle, ϕ (commencing at the beginning of a meshing cycle), an initial guess is made for the central lubricant film thickness using equation (14). The required data for contact geometry, load and speed of entraining motion are obtained through TCA⁴²

Using the film thickness obtained in step 1, the inlet boundary condition is set as $m = m *$ for both the inlet boundaries ahead of the major and minor axes of the Hertzian contact ellipse. Now a grid of $n_{x} \times n_{y}$ is used to cover the domain bounded by the inlet menisci along the x (minor axis) and y (major axis) and extend beyond the elliptical contact to include the conjunctional exit boundaries at the lubricant film rupture points x_c and y_c.

The pressure distribution and the corresponding lubricant film contour are obtained by simultaneous solution of equations (1) to (7) and (10) to (12) in an iterative manner, where two convergence criteria should be satisfied.

The first criterion seeks to compute nodal pressures within a specified limit

\sum_{i} \sum_{j} | \frac{p_{i, j}^{k} - p_{i, j}^{k - 1}}{p_{i, j}^{k}} | \leq ɛ_{p}

(38)

where

10 - 5 \leq ɛ_{p} \leq 10 - 4

. If the criterion is not satisfied, the generated pressures are under-relaxed as

p_{i . j}^{k} = p_{i, j}^{k - 1} + Ω Δ p_{i, j} \in i, j

. The under-relaxation factor is usually

0.01 \leq Ω \leq 0.8

and the steps 3–4 are repeated.

The second criterion seeks to obtain the temperature distribution. Using the obtained pressure distribution and film thickness, the temperature distribution is obtained as described in the Thermo-elastohydrodynamics section. The average temperature at any cross section through the lubricant film is used to update the lubricant viscosity and density (equations (6) and (7)). It should be noted that the temperature variation can be directly integrated into Reynolds equation. Such an approach has been reported in Peiran and Shizhu.⁶³

The convergence criteria for temperature are similar to that for pressure, but the summation is carried out volumetrically in the (x, y and z) domain as

\sum_{i} \sum_{j} \sum_{k} \frac{| T_{new} - T_{old} |}{T_{old}} \leq ɛ_{T}

(39)

where usually

10 - 5 \leq ɛ_{T} < 10 - 4

New pressure distribution is obtained using the updated values of viscosity and density, repeating steps (1)–(4). After new converged thermo-elastohydrodynamic pressures are obtained, an instantaneous load balance is sought between the lubricant reaction (equation (15)) and the applied contact load (determined through TCA). Therefore, the third convergence criterion is

| \frac{F - W}{F} | \leq ɛ_{w}

(40)

where

0.001 \leq ɛ_{w} \leq 0.05

. If the criterion is not met, the central lubricant film thickness

h_{co}

is adjusted and the entire iterative process is repeated

h_{c 0}^{l} = h_{c 0}^{l - 1} (\frac{F}{W}) ς

(41)

where the damping factor is usually in the range:

- 0.1 \leq ς \leq 0.1

In the above process, the indices i and j refer to a computational grid position, k denotes the pressure convergence iteration counter and l the load convergence iteration counter.

When all the convergence criteria are met, the pinion angle is advanced within the meshing cycle and the entire process is repeated. To observe the contact conditions for any pair of meshing teeth, the meshing cycle is sub-divided into 20 discrete steps of the pinion rotation angle.

Results and discussion

The current analysis investigates the contact conditions for a moderate to highly loaded hypoid gear pair of a C-segment medium front wheel drive transaxle vehicle. The results presented are for a typical gear teeth pair meshing cycle. Tables 1 and 2 provide the specifications of the pinion and ring gear, respectively. Figure 18 in Appendix 2 provides a graphical representation of the parameters in Tables 1 and 2. Table 3 lists the lubricant data and mechanical and thermal properties of the surfaces. The engine conditions considered correspond to an output torque of 144 Nm at 2524 r/min. The transmission is engaged in the fifth gear with a ratio of 0.702:1.

Table 1.

Pinion specification.

Pinion parameters
Number of teeth	13
Face-width	33.851 (mm)
Face angle	29.056 (°)
Pitch angle	29.056 (°)
Root angle	29.056 (°)
Spiral angle	45.989 (°)
Pitch apex	−9.085 (mm)
Face apex	1.368 (mm)
Outer cone distance	83.084 (mm)
Pinion offset	24.0000028 (mm)
Pinion hand	Right

Table 2.

Gear specification.

Gear parameters
Number of teeth	36
Face-width	29.999 (mm)
Face angle	59.653 (°)
Pitch angle	59.653 (°)
Root angle	59.653 (°)
Spiral angle	27.601 (°)
Pitch apex	8.987 (mm)
Face apex	10.948 (mm)
Outer cone distance	95.598 (mm)
Gear offset	24 (mm)

Table 3.

Lubricant and solids data.

Pressure viscosity coefficient (α)	1.39e-8 (Pa⁻¹)
Temperature viscosity coefficient (β)	0.026 (K⁻¹)
Atmospheric dynamic viscosity ( $η_{0})$	0.0499 (Pa.s)
Eyring shear stress $τ_{0}$ ⁵³	5 (MPa)
Limiting shear stress at atmospheric pressure $τ_{L 0}$ ⁵³	2.3 (MPa)
Atmospheric density $ρ_{0}$ of lubricant	829.3 (kg/m³)
Thermal conductivity of lubricant	0.14 (W/mK)
Heat capacity of lubricant	2000 (J/kgK)
Modulus of elasticity of contacting solids	210 (GPa)
Poisson’s ratio of contacting solids	0.3 (−)
Density of contacting solids	7850(kg/m³)
Thermal conductivity of contacting solids	46 (W/mK)
Heat capacity of contacting solids	470 (J/kgK)
Surface roughness Ra	0.5 (µm)
Bulk temperature	330 (K)

Simulation studies were carried out on a 2.93 GHz Pentium Intel i7 machine with the total CPU time of 5000 s.

Using TCA, contact kinematics, load and geometry are calculated. Table 4 lists the speed of entraining motion and the instantaneous radii of curvature of mating teeth during a meshing cycle.

Table 4.

Speed of entraining motion and radii of curvature of a teeth pair through mesh.

Pinion angle ϕ (Rad)	Magnitude of entraining velocity (m/s)	Principal radius $R_{zx}$ (m)	Principal radius $R_{zy}$ (m)
0.5027	12.78	0.0157	1.0067
0.5341	12.66	0.0158	1.0297
0.5812	12.48	0.0160	1.0626
0.6283	12.30	0.0162	1.0937
0.6754	12.13	0.0164	1.1228
0.7226	11.97	0.0166	1.1501
0.7697	11.80	0.0168	1.1754
0.8168	11.64	0.0171	1.1988
0.8639	11.49	0.0174	1.2204
0.9111	11.34	0.0177	1.2400
0.9582	11.20	0.0180	1.2578
1.0053	11.05	0.0183	1.2736
1.0524	10.92	0.0186	1.2876
1.0996	10.78	0.0190	1.2996
1.1467	10.66	0.0194	1.3098
1.1938	10.53	0.0198	1.3180
1.2409	10.42	0.0202	1.3243
1.2881	10.30	0.0206	1.3288
1.3352	10.19	0.0211	1.3313

The hypoid gear pair is subject to high contact loads during a typical meshing cycle under the simulated conditions. This is shown in Figure 4. There is also significant side leakage of the lubricant from the contact along the major axis of the elliptical contact footprint. This is also shown in the same figure, indicating that the side leakage flow accounts for 50% of the inlet entrained lubricant flow. Therefore, solutions for hypoid gear pair meshing should always take into account the angled flow entrainment into the contact and embody two-dimensional solution of Reynolds equation irrespective of the large aspect ratio of the elliptical contact footprint during parts of the meshing cycle. High contact loads and significant side leakage flow out of the contact domain yield thin thermo-elastohydrodynamic films with high generated pressures and viscous shear, thus promoting non-Newtonian conditions, which significantly deviate from the usual isothermal elastohydrodynamic analyses reported in earlier literature.

Figure 4.

Contact load and percentage side leakage flow in a meshing cycle.

Solutions are obtained for 10 instantaneous contacting locations of the hypoid gear teeth pair during the simulated meshing cycle, using a 60 × 60 grid across the elliptical contact footprint. The number of grid points is chosen after a grid sensitivity study. The result of this test is shown in Figure 5 at the centre of the meshing cycle. Results are presented for 60 × 60 and 120 × 120 grid points. The latter predicts 0.57 µm and 0.60 µm for the central and minimum exit film thickness positions respectively. These results are 0.56 µm and 0.58 µm, when using a 60 by 60 computational grid. This difference translates to 1.75% and 3.3% absolute differences for the central and the minimum exit oil film thickness, respectively. The temperature distribution at any grid location (cross section into the depth of the lubricant film and the adjacent solid contacting surfaces) is obtained at 25 equally spaced grid locations in each case (for the lubricant film thickness and solid surface layers).

Figure 5.

Pressure distribution and film thickness along the minor axis for different grid numbers.

Figure 6 shows the pressure distribution for an isothermal, as well as a thermal non-Newtonian analysis for the pinion angle of 1.3352 rad (this marks the end of the meshing cycle). The main difference in the results is the diminution of the secondary pressure peak at the contact outlet in the latter case, caused by reduced lubricant viscosity in a thermal contact.

Figure 6.

3D pressure distribution at the end of meshing cycle.

Figure 7 shows the cross-sectional pressure profile along the minor axis of the contact footprint at the mid-point of the major axis of the 3D distributions in Figure 6. It also shows the corresponding film shape. In general, there is an insignificant difference in the lubricant film thickness between the two analyses, except for a small change in the minimum film thickness due to reduced lubricant viscosity when thermal analysis is included.

Figure 7.

Pressure profile and corresponding film shape at the contact mid-section along the minor axis of the elliptical contact footprint.

The results in Figure 7 correspond to the instance at the end of the meshing cycle. Figure 8 shows the variation in the minimum film thickness and the maximum primary pressure peak (corresponding to the equivalent maximum Hertzian pressure) during a typical steady state meshing cycle. The maximum pressure occurs approximately at mid-meshing cycle, where the absolute minimum film thickness occurs as well. The thermo-elastohydrodynamic minimum film thickness is significantly lower than that predicted by the idealised isothermal analysis due to reduced lubricant viscosity with temperature, thus its load-carrying capacity.

Figure 8.

Variation of maximum primary peak pressure and minimum film thickness during a meshing cycle.

An important point to note in Figures 7 and 8 is that the film thickness is barely in excess of the root mean square roughness of the contiguous surfaces (measured to be 0.5 µm). Therefore, a mixed regime of lubrication would be expected, which is in fact the case throughout the meshing cycle. This is the reason for the current mixed thermo-elastohydrodynamic analysis.

Figure 9 shows the 3D average (through the depth of the film) lubricant temperature distributions at the mid and at the end of the meshing cycle, where maximum and minimum loading occurs. In this transition, the contact load decreases from 4 kN to 750 N, whilst the side leakage flow is reduced by 20% (see Figure 4). At the same time, the maximum lubricant temperature is reduced by 80 K. This reduction in temperature not only shows the significance of compressive heating, but also demonstrates the important role of shear heating. The latter is influenced by viscosity variation under high pressures.^35,36 The contours of lubricant temperature at mid-section of the major axis of the elliptical contact footprint, along its minor axis are also provided in the figure. These show the temperature variation in the lubricant, as well as the adjacent solid contacting surfaces; the equivalent ellipsoidal solid (Figure 2) and the semi-infinite elastic half space. Figure 9(c) also shows the velocity distribution, based on the equations (17) to (20), as well as that based on a non-Newtonian assumption provided by Conry et al.⁴³ These velocity distributions are for different positions along the minor axis of the contact footprint for the highly loaded condition at the centre of the meshing cycle. This condition corresponds to the most severe non-Newtonian behaviour and is the most appropriate condition for ascertaining the validity of using a Newtonian velocity distribution in the energy equation.

Figure 9.

(a) 3D temperature distributions and contours for maximum and minimum loaded conditions along the meshing cycle; (b) velocity distributions.

The results show that for a realistic prediction of film thickness, it is necessary to undertake thermo-elastohydrodynamic analysis. Figure 10 shows the results for the average lubricant temperature profile (into the depth of the film) along the minor axis of the contact ellipse at two instances along the meshing cycle. The corresponding temperature of the surfaces (equivalent ellipsoidal solid and the semi-infinite plane, see Figure 2) is also shown. As expected, the lubricant temperature is higher than the adjacent solid surfaces due to its shear heating. The lubricant conducts the heat into the solid boundaries. The temperature of the semi-infinite elastic half-space (represented by the lower surface velocity of the gear surface) remains higher than that of the equivalent ellipsoidal solid which is represented by the higher pinion surface velocity. At the start of the meshing cycle, the lubricant temperature is that of the bulk and ascends to the temperature of the solid surfaces, assumed to be the same. In physical terms, the hotter surfaces convect heat into the incoming lubricant at the bulk temperature of the differential bath. At the end of the meshing cycle, the temperature of lubricant and the surfaces converge to that of the lubricant bath. The maximum temperature occurs at the maximum contact pressure, not the maximum sliding velocity. This is due to two effects; first compressive heating, and second the rise in lubricant viscosity due to increased pressures and subsequent higher shear heating.

Figure 10.

Cross-sectional temperature profile for maximum and minimum loaded conditions along the minor axis of the contact ellipse during a meshing cycle.

Figure 11 shows the maximum lubricant temperature in the contact as well as its maximum exit temperature. At any instant in the meshing cycle, the difference between these two temperatures is indicative of the lubricant heat loss mainly into the adjacent solid surfaces for the predicted thin lubricant films. The temperature variation of the surfaces closely follows that of the lubricant with lower values for the moving surface (semi-infinite plane) and the assumed stationary surface (the equivalent ellipsoidal solid). The surfaces conduct the heat to the assumed penetration depth of d (see boundary conditions (23)).This means that the lubricant exit temperature would always remain below that of the solid surfaces, as would be expected.

Figure 11.

Transient variation in the lubricant and solid surface temperatures for a meshing cycle.

As already noted in the Introduction section, reduction of friction (thus, improved efficiency) is the primary objective in transmission engineering. Therefore, prediction of friction, comprising viscous shear and asperity interactions, is important. Figure 12 shows the viscous shear stress distribution along the minor axis of the contact ellipse through the mid major axis location. The Eyring shear stress limit is shown in the figure, together with the limiting shear stress of the lubricant (equation (8)) for given local pressure, p. The figure also includes the predicted shear stress for isothermal Newtonian and non-Newtonian analyses. It can be seen that in both cases, the shear stress distributions exceed both the Eyring and limiting shear stress, particularly excessively in the case of the former. This indicates that assumption of isothermal Newtonian behaviour cannot be upheld. When thermal effects are taken into account, more realistic shear distributions are predicted, because under isothermal non-Newtonian condition the shear stress is allowed to exceed its limiting value. The shear stress falls in the region between its Eyring value and the limiting shear. Therefore, the appropriate viscous shear stress characteristics is thermal non-Newtonian. At the end of the meshing cycle for lower contact load and sliding velocity, both thermal Newtonian and non-Newtonian shear stress distributions coincide and fall below Eyring shear stress. This indicates that for low to medium contact loads and sliding velocity, a thermal Newtonian traction regime may be assumed.

Figure 12.

Shear stress distribution along the minor axis of contact ellipse.

With predicted thin films, some of the load is carried by the interaction of asperities on the counter face surfaces. Figure 13 shows the asperity contact area as a proportion of the apparent area as well as the percentage of load carried by direct boundary interactions. These account for several percentage of the total contact area and contact load. The same is also true in the case of generated friction (Figure 14). The boundary friction contribution accounts for a larger proportion of generated friction at the beginning and at the end of the meshing cycle. This behaviour is as a result of complex interactions between several effects. Lower sliding velocities at the beginning and the end of meshing cycle leads to reduced shear rate and shear stress. On the other hand, this also leads to reduced shear thinning and a higher effective lubricant viscosity. Hence, a rise in shear stress as well as thicker film thickness would be expected, leading to reduced asperity interactions. A balance between these effects is found by the iterative numerical solution which is presented here.

Figure 13.

Asperity contribution to contact conditions.

Figure 14.

Friction variation during a meshing cycle.

Finally, the variation of coefficient of friction is shown in Figure 15. For purely elastohydrodynamic conditions, one would expect a coefficient of friction < 0.05.³⁴ It also shows that a constant value cannot be assumed with such transient meshing conditions.

Figure 15.

Variation of coefficient of friction during a meshing cycle.

It is noteworthy that the isothermal Newtonian model yields high values of coefficient of friction, given higher predicted values of shear stress in Figure 14. Of course, this is not a realistic representation of the prevailing conditions.

Concluding remarks

The results show that under realistic driving conditions, the differential hypoid gear of a C-segment vehicle is subjected to high contact loads and thin elastohydrodynamic films of the order of the composite roughness of the contacting surfaces. Angled lubricant entrainment flow occurs into the elliptical contact footprint with significant side leakage flow out of the contact, resulting in thinning of the lubricant film. Furthermore, high rate of shear for the investigated conditions promotes non-Newtonian lubricant behaviour, as well as the generated temperature due to both compressive and shear heating. Thus, a realistic analysis should be mixed non-Newtonian thermo-elastohydrodynamics, such as the one reported in this paper. In fact, the analysis shows that asperity friction can account for up to 7–8% of the total friction. Nevertheless, the generated heat has a small effect on the film thickness relative to an isothermal analysis, except for the minimum film thickness as the result of reduced lubricant viscosity with the generated temperature and shear thinning. The effect of temperature is more pronounced on the viscous shear stress distribution, indicating lubricant shear stress residing within the Eyring regime of traction, unlike unreasonable predictions with an isothermal analysis exceeding the lubricant limiting shear stress. It should be noted that the current analysis does not take into account the heat generated as the result of local interaction of counterface asperities. The inclusion of this constitutes a future development of the model. One should also note that the presented model can be enhanced by incorporating more recent rheological models as described in literature.^35,36

Footnotes

Acknowledgement

The authors would like to express their gratitude to Advanced Numerical Solutions Inc. for access to the TCA tool CALYX.

Declaration of Conflicting Interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

Appendix

References

Rahnejat H. (ed.) Tribology and dynamics of engine and powertrain. Cambridge, UK: Woodhead Publishing Ltd., 2010.

Krenz R. Vehicle response to throttle tip in-tip out. SAE Trans., Paper no. 850967, 1985.

Menday

Rahnejat

Ebrahimi

. Clonk: an onomatopoeic response in torsional impact of automotive drivelines. Proc IMechE, Part D: J Automobile Engineering 1999; 213: 349–357.

Vafaei

Menday

Rahnejat

. Transient high-frequency elasto-acoustic response of a vehicular drivetrain to sudden throttle demand. Proc IMechE, Part K: J Multi-body Dynamics 2001; 215: 35–52.

Tangasawi

Theodossiades

Rahnejat

. Lightly loaded lubricated impacts: Idle gear rattle. J Sound Vibrat 2007; 308: 418–430.

Ottewill

Neild

Wilson

. Intermittent gear rattle due to interactions between forcing and manufacturing errors. J Sound Vibrat 2009; 321: 913–935.

Rocca

Russo

. Theoretical and experimental investigation into the influence of the periodic backlash fluctuations on the gear rattle. J Sound Vibrat 2011; 330: 4738–4752.

Rahnejat

. Multi-body dynamics: vehicles, machines and mechanisms, Bury St Edmunds, UK: Professional Engineering Publishing, 1998.

Brancati

Rocca

Russo

. A gear rattle model accounting for oil squeeze between the meshing gear teeth. Proc IMechE, Part D: J Automobile Engineering 2005; 219: 1075–1083.

10.

De la Cruz

Theodossiades

Rahnejat

. An investigation of manual transmission drive rattle. Proc IMechE, Part K: J Multi-body Dynamics 2010; 224: 167–181.

11.

De la Cruz

Chong

WWF

Teodorescu

et al.

Transient mixed thermo-elastohydrodynamic lubrication in multi-speed transmissions. Tribol Int 2012; 49: 17–29.

12.

Denny CM. Mesh friction in gearing. AGMA, Technical Paper No. 98FTM2, 1998.

13.

Michlin

Myunster

. Determination of power losses in gear transmissions with rolling and sliding friction incorporated. Mech Mach Theory 2002; 37: 167–174.

14.

Misharin YA. Inﬂuence of the friction condition on the magnitude of the friction coefﬁcient in the case of rollers with sliding. In: Proceedings international conference on gearing. London, UK: IMechE, 1958, pp.159–164.

15.

Benedict

Kelly

. Instantaneous coefﬁcients of gear tooth friction. Trans ASLE 1960; 4: 59–70.

16.

O’Donoghue

Cameron

. Friction and temperature in rolling sliding contacts. Trans ASLE 1966; 9: 186–194.

17.

Drozdov

Gavrikov

. Friction and scoring under the conditions of simultaneous rolling and sliding of bodies. Wear 1968; 11: 291–302.

18.

Gunda

Singh

. Effect of sliding friction on the dynamics of spur gear pair with realistic time-varying stiffness. J Sound Vibrat 2007; 301: 927–949.

19.

Velex

Cahouet

. Experimental and numerical investigations on the influence of tooth friction in spur and helical gear dynamics. Trans ASME J Mech Des 2000; 122: 515–522.

20.

Velex

Sainsot

. An analytical study of tooth friction excitations in errorless spur and helical gears. Mech Mach Theor 2002; 37: 641–658.

21.

Kar

Mohanty

. An algorithm for determination of time-varying frictional force and torque in a helical gear system. Mech Mach Theor 2007; 42: 482–496.

22.

Karagiannis

Theodossiades

Rahnejat

. On the dynamics of lubricated hypoid gears. Mech Mach Theor 2012; 48: 94–120.

23.

Grubin

. Contact stresses in toothed gears and worm gears. Central Scientific Research Institute for Technology and Mechanical Engineering 1949; 30: 337–337.

24.

Kahraman

. A transient mixed elastohydrodynamic lubrication model for spur gear pairs. Trans ASME J Tribol 2010; 132: 011501 (9p.)>–011501 (9p.).

25.

Mohammadpour

Theodossiades

Rahnejat

. Elastohydrodynamic lubrication of hypoid gears at high loads. Proc IMechE, Part J: J Engineering Tribology 2012; 226: 183–198.

26.

Mohammadpour

Theodossiades

Rahnejat

et al.

Non-Newtonian mixed elastohydrodynamics of differential hypoid gears at high loads. Meccanica 2014; 49: 1115–1138.

27.

Crook

. The Lubrication of Rollers III. A theoretical discussion of friction and the temperatures in the oil film. Phil Trans R Soc Lond 1961; 254: 237–258.

28.

Cheng

Sternlicht

. A numerical solution for the pressure, temperature, and film thickness between two infinitely long lubricated rolling and sliding cylinders under heavy loads. Trans ASME J Basic Eng 1965; 87: 695–707.

29.

Dowson

Whitaker

. A numerical procedure for the solution of the elastohydrodynamic problem of rolling and sliding contacts lubricated by a Newtonian ﬂuid. Proc IMechE 1966; 180: 57–71.

30.

Kim

Sadeghi

. Three-dimensional temperature distribution in EHD lubrication. Part I: circular contact. Trans ASME J Tribology 1992; 114: 32–41.

31.

Kim

Sadeghi

. Three-dimensional temperature distribution in EHD lubrication. Part II: point contact and numerical formulation. Trans ASME J Tribology 1993; 115: 36–45.

32.

Yang

Cui

Jin

et al.

Transient elastohydrodynamic analysis of elliptical contacts. Part 2: thermal and Newtonian lubricant solution. Proc IMechE, Part J: J Engineering Tribology 2005; 219: 187–200.

33.

Yang

Cui

Jin

et al.

Transient elastohydrodynamic analysis of elliptical contacts. Part 3: Newtonian lubricant solution under iso-thermal and thermal conditions. Proc IMechE, Part J: J Engineering Tribology 2007; 221: 63–73.

34.

Sharif

Evans

Snidle

et al.

Modeling of film thickness and traction in a variable ratio traction drive rig. J Tribol 2004; 126: 92–104.

35.

Habchi

Vergne

Bair

et al.

Influence of pressure and temperature dependence of thermal properties of a lubricant on the behaviour of circular TEHD contacts. Tribol Int 2010; 43: 1842–1850.

36.

Paouris

Rahmani

Theodossiades

et al.

An analytical approach for prediction of elastohydrodynamic friction with inlet shear heating and starvation. Tribol Lett 2016; 64: 10–10.

37.

Kolivand

Kahraman

. Prediction of mechanical gear mesh efﬁciency of hypoid gear pairs. Mech Mach Theor 2010; 45: 1568–1582.

38.

Simon

. Inﬂuence of machine tool setting parameters on EHD lubrication in hypoid gears. Mech Mach Theor 2009; 44: 923–937.

39.

Gohar

. Oil film thickness and rolling friction in elastohydrodynamic point contact. Trans ASME, J Lubn Tech 1971; 93: 371–382.

40.

Jalali-Vahid

Rahnejat

Gohar

et al.

Prediction of oil-film thickness and shape in elliptical point contacts under combined rolling and sliding motion. Proc IMechE, Part J: J Engineering Tribology 2000; 214: 427–437.

41.

Wang

Yang

et al.

Mixed elastohydrodynamic lubrication with three-dimensional machined roughness in spiral bevel and hypoid gears. J Tribol 2015; 137: 041503–041503.

42.

Litvin

Fuentes

Fan

et al.

Computerized design, simulation of meshing, and contact and stress analysis of face-milled format generated spiral bevel gears. Mech Mach Theor 2002; 37: 441–459.

43.

Conry

Wang

Cusano

. A Reynolds-Eyring equation for elastohydrodynamic lubrication in line contacts. Trans ASME J Tribol 1987; 109: 648–658.

44.

Johnson

Tevaarwerk

. Shear behaviour of elastohydrodynamic oil films. Proc R Soc Lond A 1977; 356: 215–236.

45.

Gohar

Rahnejat

. Fundamentals of tribology, London, UK: Imperial College Press, 2008.

46.

Roelands CJA. Correlation aspects of viscosity–temperature–pressure relationship of lubricating oils. PhD Thesis, Delft University of Technology, Delft, The Netherlands, 1966.

47.

Dowson

Higginson

. Elastohydrodynamic lubrication, Oxford, UK: Pergamon Press, 1966.

48.

Hamrock

Dowson

. Isothermal elastohydrodynamic lubrication of point contacts, Part I – Theoretical formulation. Trans ASME J Lubn Technol 1976; 98: 223–229.

49.

Birkhoff

Hays

. Free boundaries in partial lubrication. J Math Phys 1963; 32: 2–2.

50.

Tipei

. Boundary conditions of a viscous flow between surfaces with rolling and sliding motion. Trans ASME J Lubn Technol 1968; 90: 8–16.

51.

Mohammadpour

Johns-Rahnejat

Rahnejat

et al.

Boundary conditions for elastohydrodynamics of circular point contacts. Tribol Lett 2014; 53: 107–118.

52.

Hamrock

Dowson

. Isothermal elastohydrodynamic lubrication of point contacts, Part IV – Starvation results. Trans ASME J Lubn Tech 1977; 99: 15–23.

53.

Gohar

. Elastohydrodynamics, London, UK: Imperial College Press, 2001.

54.

Swift

. The wetability of lubricating films in journal bearings. Minutes Proc J Inst Civil Engrs 1932; 223: 267–288.

55.

Stieber

. Das schwimmlager: Hydrodynamische theorie des gleitlagers, VDI-Verlag, 1933.

56.

Habchi

Vergne

Bair

et al.

Influence of pressure and temperature dependence of thermal properties of a lubricant on the behaviour of circular TEHD contacts. Tribol Int 2010; 43: 1842–1850.

57.

Masseth J and Kolivand M. Lapping and superfinishing effects on hypoid gears surface finish and transmission errors. In: ASME conference, international design engineering technical conferences and computers and information in engineering conference (IDETC), Las-Vegas, Nevada, USA, 2007, pp.759–765.

58.

Greenwood

Tripp

. The contact of two nominally flat rough surfaces. Proc IMechE 1970–1971; 185: 625–633.

59.

Teodorescu

Balakrishnan

Rahnejat

. Integrated tribological analysis within a multi-physics approach to system dynamics. Tribol Interf Eng Ser 2005; 48: 725–737.

60.

Kim

Ehret

Dowson

et al.

Thermal elastohydrodynamic analysis of circular contacts. Part 1: Newtonian model. Proc IMechE, Part J: J Engineering Tribology 2001; 215: 339–352.

61.

Kim

Ehret

Dowson

et al.

Thermal elastohydrodynamic analysis of circular contacts. Part 2: non- Newtonian model. Proc IMechE, Part J: J Engineering Tribology 2001; 215: 353–362.

62.

Jalali-Vahid

Rahnejat

Jin

et al.

Transient analysis of isothermal elastohydrodynamic circular point contacts. Proc IMechE, Part C: J Mechanical Engineering Science 2001; 215: 1159–1172.

63.

Peiran

Shizhu

. A generalized Reynolds equation for non-Newtonian thermal elastohydrodynamic lubrication. J Tribol 1990; 112: 631–636.