Analytical Solution of Reynolds Equation: Grubin’s Approximation

Manoj Rajankunte Mahadeshwara
13.02.2017
Hertz Vs. EHL Pressure

Derivation of famous analytical solution of EHL equations is presented here. In 1949, Grubin obtained a solution for so called elasto-hydrodynamic lubrication (EHL) line contact problem with simplifications known as Grubin approximation. He was the first to combine both elastic deformation and lubricant hydrodynamic flow, including pressure-viscosity dependence. Although his solution did not satisfy strictly hydrodynamic equations of EHL, his analysis was recognized as very useful under highly loaded conditions, not only in stationary case, but also in transient conditions. Corresponding MATLAB code can be found here.

Here, the solution approach is described for the case of stationary line contact (1D Reynolds equation). The base for the validity of Grubin’s approximation is the fact that at highly loaded contacts, the hydrodynamic pressure is almost the same as Hertz contact pressure obtained for dry contacts (cylinder on plane). On top of it, the film thickness is almost constant within Hertz contact area, as shown in Fig. 1.

Hertz Vs. EHL Pressure
Fig. 1. Hertz Vs. EHL

Following system of equations is considered:

(1)   \begin{eqnarray*} \frac{d}{x} (\frac{h^3}{12\mu} \frac{dp}{dx}) = U_m\frac{dh}{dx} \\ \label{complete_sys2} \mu= \mu_0 e^{\alpha p} \\ \label{complete_sys3} p(-\infty) = 0 \end{eqnarray*}

where h, p, \mu, \alpha are hydrodynamic film thickness, pressure, viscosity, and pressure-viscosity coefficient correspondingly. Integration of Reynolds equation gives:

banner

\frac{h^3}{12\mu} \frac{dp}{dx} = h + C

where C is an integration constant. Since the pressure drops to zero in the infinity, there exist a point x' , so that \frac{dp}{dx}=0 . Hence, C=-h' , where h'=h(x') . Following expression can be obtained from the integrated Reynolds equation:

(2)   \begin{eqnarray*} \frac{dp}{dx} = \frac{h-h'}{h^3} 12\mu U_m\\ \end{eqnarray*}

Further, a following variable is introduced: p_{\mu} = \frac{1-e^{-\alpha p}}{\alpha} . Then, the equation is true:

banner

(3)   \begin{eqnarray*} \frac{dp_{\mu}}{dx} = e^{-\alpha p} \frac {dp}{dx} \\ \label{complete_sys2} => \frac{dp_{\mu}}{dx} = 12\mu_0 U_m \frac{h-h'}{h^3} \\ \end{eqnarray*}

When pressure is very high, \ p_{\mu} -> \frac{1}{\alpha} \equiv const and => \frac{dp_{\mu}}{dx} = 0 .

From this it follows that h'=h_0 and film is constant within the Hertz contact area. It is clear from Fig. 1 that outside the Hertz area, the hydrodynamic pressures are small compared to those within the contact and therefore it can be assumed that all the elastic defection is due to Herztian pressure. Therefore, the film thickness can be written in the form:

(4)   \begin{eqnarray*} h = h_0 + h_{elastic} = h_0 + \frac{W_L \delta}{E_L} \\ \end{eqnarray*}

where h_{elastic} is the combined elastic deflection of the surfaces (given by Hertz solution), W_L is the applied load per unit length of the cylinder, \frac{1}{E_L} = \frac{1}{\pi}(\frac{1-\nu_{1}^{2}}{E_1} + \frac{1-\nu_2^{2}}{E_2}) is the combined elastic modulus and \delta is given by the following expression:

(5)   \begin{eqnarray*} \delta = 2 \{\frac{x}{a} \sqrt{\frac{x^2}{a^2}-1} - ln{\frac{x}{a} + \sqrt{\frac{x^2}{a^2} - 1}} \} \\ \end{eqnarray*}

From the expression of the film thickness h-h' = h-h_0 = \frac{W_L \delta}{E_L} and hence (equation 2):

(6)   \begin{eqnarray*} => \frac{dp_{\mu}}{dx} = 12\mu_0 U_m \frac{\frac{W_L \delta}{E_L}}{(h_0 + \frac{W_L \delta}{E_L})^3} \\ => \frac{dp_{\mu}}{dx} = \frac{12\mu_0 U_m}{(\frac{W_L }{E_L})^2}\frac{\delta}{H^3} \\ => \frac{(\frac{W_L }{E_L})^2}{12\mu_0 U_m}\frac{dp_{\mu}}{dx} = \frac{\delta}{H^3}\\ \end{eqnarray*}

By introducing dimensionless variables p_{\mu}^{'} =p_{\mu} \frac{(\frac{W_L }{E_L})^2}{12\mu_0 U_m a} and x'=\frac{x}{a} , the equation reduces to the following form:

(7)   \begin{eqnarray*} \frac{dp_{\mu}^{'}}{dx'} = \frac{\delta}{H^3} \end{eqnarray*}

Integration of this equation from -\infty to 1 gives p'(-1) = \int_{-\infty}^{-1}{\frac{\delta}{H^3}} dx' . This equation was integrated by Grubin numerically for various values of H_0 and further fitted with the following simple formula:

(8)   \begin{eqnarray*} p'(-1) = 0.0986H_0^{-\frac{11}{8}} \\ \end{eqnarray*}

Now, for the line contact, a=2\sqrt{\frac{W_L R}{E_L}} with R being equivalent radius of the bodies. At x'=-1 in dimensional coordinates x=-a . Hence,

(9)   \begin{eqnarray*} p_{\mu}(x=-a) = \frac{1}{\alpha} = \frac{12\mu_0 U_m a}{\frac{W_L }{E_L})^2} p'(-1) \\ \end{eqnarray*}

Substituting the expression for dimensionless pressure one can obtain following:

(10)   \begin{eqnarray*} \frac{1}{\alpha} = \frac{12\mu_0 U_m a}{\frac{W_L }{E_L})^2} 0.0986(h_0 W_L/E_L)^{-\frac{11}{8}}\\ \end{eqnarray*}

And finally, the film thickness can be expressed as

(11)   \begin{eqnarray*} h_0^{\frac{11}{8}} = 2.3664 \frac{\mu_0 U_m \alpha \sqrt{R}}{ (W_L/E_L)^{\frac{1}{8}}}\\ \end{eqnarray*}

I am currently working as a Postgraduate Researcher at the University of Leeds, where I am actively involved in research activities. Prior to this, I successfully completed my master's degree through the renowned Erasmus Mundus joint program, specializing in Tribology and Bachelor's degree in Mechanical Engineering from VTU in Belgaum, India. Further I handle the social media pages for Tribonet and I have my youtube channel Tribo Geek.

Be the first to comment

Leave a Reply

This site uses Akismet to reduce spam. Learn how your comment data is processed.