Next Article in Journal
A Digital Twin for Friction Prediction in Dynamic Rubber Applications with Surface Textures
Next Article in Special Issue
Effect of Temperature and Surface Roughness on the Tribological Behavior of Electric Motor Greases for Hybrid Bearing Materials
Previous Article in Journal
An Identification Method for Orifice-Type Restrictors Based on the Closed-Form Solution of Reynolds Equation
Previous Article in Special Issue
Grease Performance Requirements and Future Perspectives for Electric and Hybrid Vehicle Applications
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Technical Note

Formulation to Calculate Isothermal, Non-Newtonian Elastohydrodynamic Lubrication Problems Using a Pressure Gradient Coordinate System and Its Verification by an Experimental Grease

TriboLogics Corporation, Tokyo 110-0016, Japan
Submission received: 16 February 2021 / Revised: 28 April 2021 / Accepted: 11 May 2021 / Published: 14 May 2021
(This article belongs to the Special Issue Grease)

Abstract

:
This paper presents a formulation of point contact elastohydrodynamic lubrication analysis for an isothermal, non-Newtonian flow. A coordinate system of the pressure gradient was employed herein. A Couette flow and a Poiseuille flow were considered along the directions of the zero and non-zero pressure gradients, respectively. The Poiseuille flow velocity was assumed to be represented by a 4th-order polynomial of z along the film thickness direction. The Couette flow velocity was assumed to be represented by a linear function of z. Subsequently, the modified Reynolds equation, which contains an equivalent viscosity, was obtained. Using Bauer’s rheological model, the formulation presented in this study was applied to a grease that has been previously experimented upon. The results of previous studies were compared with those of the present study and a reasonable agreement was noted. The distribution of the equivalent viscosity showed a notable difference from that of Newtonian flow. The formulation can be incorporated easily to the usual elastohydrodynamic lubrication calculation procedure for Newtonian flow. The method can be easily applied to other non-Newtonian rheological models. The equivalent viscosity can be calculated using the one-parameter Newton-Raphson’s method; as a result, the calculation can be performed within a reasonable time.

1. Introduction

Performing experiments on non-Newtonian flows, including grease flows, is considerably time-consuming and costly. Therefore, it is important to numerically analyze the phenomena corresponding to non-Newtonian flows. Numerical approaches can help obtain a variety of data that cannot be determined experimentally. Grease flows can be well defined using Bauer’s model; however, owing to the extreme complexity of this model, it is difficult to determine the exact solution as well as approximate solutions for point contact, isothermal, non-Newtonian elastohydrodynamic lubrication (EHL) analyses. Therefore, it is convenient that the non-Newtonian EHL calculation can be executed within a reasonable calculation time and without large modification to the usual Newtonian EHL calculation procedure. Kochi et al. [1] performed experiments on grease under soft EHL conditions and measured the film thickness and traction forces. The method proposed in the present study can be applied to the grease considered in Kochi et al. [1] to validate this theoretical approach.

1.1. Classification of Calculation Methods

As shown in Figure 1, the Z-direction is considered to be the film thickness direction. The flow velocities along the X- and Y-directions are denoted by u and v , respectively, which are functions of x, y, and z; however, when considering only z dependency, the velocities can be expressed as u (z) and v (z), respectively. The numerical methods for isothermal, non-Newtonian EHL analyses can be classified in terms of the accuracy of u (z) and v (z), as follows:
  • Method 1: Exact solution of u (z) and v (z) is obtained.
  • Method 2: Approximate solution of u (z) and v (z) is obtained.
Method 2 can be further classified in terms of the employed coordinate system.
  • Method 2-A: Local X-direction is the sliding direction.
  • Method 2-B: Local X-direction is the direction of the pressure gradient.

1.2. Previous Works

Several methods have been proposed to solve isothermal, non-Newtonian EHL problems. Kauzlarich and Greenwood [2] obtained the exact solution of line contact EHL problems considering the Herschel–Bulkley model. Conry et al. [3] obtained an exact solution of line contact EHL problems by considering Eyring’s model. Specifically, the velocity u (z) was represented by a function containing cosh, indicating that, in general, u (z) and v (z) cannot be precisely represented using polynomials of z. Dong and Qian [4] obtained an approximate solution of line contact EHL problems considering Bauer’s model and using the weighted residual method. Peiran and Shizhu [5] proposed a method to obtain an exact solution of point contact EHL problems for general rheological models by using an equally divided Z-direction mesh. Subsequently, the authors applied the method to a line contact EHL problem. Kim et al. [6] attempted to obtain an exact solution of point contact EHL problems by considering Eyring’s model. Ehret et al. [7] considered the X-direction to align with the sliding direction and obtained an approximate solution of u (z) and v (z) represented by the 2nd order polynomial of z by using the perturbation method. Thus, researchers have obtained two effective viscosities: one along the sliding direction and another along the perpendicular direction.
Greenwood [8] focused on the considerable amount of computation time required to obtain an exact solution and compared two approximation methods. Sharif et al. [9] considered the X-direction to align with the sliding direction and developed a method to obtain an exact solution of point contact EHL problems for an arbitrary rheological model. Using this approach, the authors obtained two effective viscosities along the X- and Y-directions. Liu et al. [10] formulated a method to obtain an exact solution of point contact thermal EHL problems considering Eyring’s model. Yang et al. [11] formulated a general Reynolds equation for point contact EHL problems by dividing the flow into Couette and Poiseuille flows. Subsequently, the authors obtained an exact solution of line contact EHL problems considering the power law model and demonstrated the effectiveness of their proposed method. Furthermore, Bordenet et al. [12] obtained an exact solution of pure rolling point contact EHL problems considering Bauer’s model for n = 1 / 2 and applied the approach to grease.

2. Overview of the Proposed Method

In this work, an isothermal, non-Newtonian EHL formulation considering Method 2-B was developed. Although this approach does not yield the exact solution of u (z) and v (z), the calculation is simple and fast. As shown in Figure 2, the local Xc-direction is considered as the direction of the pressure gradient d , and the Yc-direction is considered to be perpendicular to Xc. The flow velocities in the Xc- and Yc-directions are denoted as u c and v c , respectively. As the pressure gradient toward Yc is zero, the flow v c is assumed to be a Couette flow, which can be represented using a linear function of z. As the pressure gradient toward Xc is generally non-zero, the flow u c is assumed to be a Poiseuille flow, which can be represented using a 4th-order polynomial of z. As u c and v c cannot be precisely represented by polynomials of z, they are expanded using polynomials of z. In such cases, the 6th-order or even higher order polynomials can be considered; however, in this work, a lower 4th-order polynomial was employed. A viscosity corresponding to the Newtonian flow was obtained and termed as the equivalent viscosity. To replace the viscosity of the Newtonian flow with the equivalent viscosity, which is a typical process when evaluating EHL problems, the method proposed by Venner and Lubrecht [13] can be used without any change for the isothermal, non-Newtonian EHL calculation. Given a rheological model, any non-Newtonian isothermal point contact EHL problem can be solved using the proposed method.
The calculation procedure consists of two steps. First, the Newtonian EHL calculation is executed for the base oil that yields the Newtonian pressure distribution. Second, using the pressure distribution as the initial value, the non-Newtonian EHL calculation is executed for the non-Newtonian flow. The local coordinate system, Oc, Xc, and Yc, at a particular point, depends on the pressure gradient and it is determined simultaneously in the process to obtain the pressure distribution. The equivalent viscosity is calculated based on the local coordinate system in each iteration loop to obtain the pressure distribution. Here, the method was applied to an experimental grease characterized by Bauer’s model.

3. Calculation of Velocity Distribution as a Function of z

Figure 1 presents the XYZ coordinate system. The force balance of a fluid can be expressed as follows [10]:
τ x   z = P x
  τ y   z = P y
where P is the pressure of the fluid, and τ x and τ y are the shear stresses in the X- and Y-directions, respectively. The parameters τ and γ ˙ are defined as follows:
τ = τ x 2 + τ y 2    
γ ˙ = u z 2 + v z 2
In Bauer’s model, τ is assumed to be represented as a function of γ ˙ , as follows [1,14,15]:
τ = ( τ 0 + k 1 · γ ˙ + k 2 · γ ˙ n   ) · η n η 0  
Here, τ 0 , k 1 , k 2 , and n are Bauer’s rheological parameters, and η n and η 0 represent the P dependent viscosity and ambient viscosity of the base oil, respectively. In this work, according to Dong and Qian [4], the parameters τ 0 , k 1 , k 2 , and n were assumed to be P -independent known values determined from the τ - γ ˙ curve measured at the ambient pressure. In Eyring’s model, according to Conry et al. [3] and Johnson and Tevaarwerk [16], the relationship between τ and γ ˙ can be expressed as follows:
γ ˙ = τ 0 η n sin h ( τ τ 0 )
Here, τ 0 is Eyring’s rheological parameter. Therefore, the relationship between τ and γ ˙ can be rewritten as follows:
τ = τ 0 · sin h 1 ( η n γ ˙ τ 0 )
The effective viscosity η * can be defined as follows:
η * = τ γ ˙  
The effective viscosity η * is a function of γ ˙ , which in turn is a function of z. The shear stresses τ x and τ y are assumed to be represented as
τ x = η *   u z  
τ y = η *   v z  
Substituting Equations (7) and (8) into Equations (1) and (2), respectively, yields
z ( η *   u z ) = P x
z ( η *   v z ) = P y
The pressure gradient vector d is defined as follows:
d = ( P x ,   P y )
Similar to the method employed by Yang et al. [11], this method involves the flow being divided into Couette and Poiseuille flows. As shown in Figure 2, the Xc-direction is considered to be along d , and its direction vector is n 1 . The Yc-direction is perpendicular to d , and its direction vector is n 2 . The local coordinate system, Oc, Xc, and Yc, at a given point depends on the pressure gradient and is determined simultaneously in the process to obtain the pressure distribution. In the Xc and Yc coordinate system, Equations (9) and (10) can be rewritten as follows:
z ( η *   u c z ) = d
z ( η *   v c z ) = 0
d = d · d
Furthermore, the velocities u c (z) and v c (z) satisfy the following boundary conditions:
u c ( 0 ) = U ,   u c ( h ) = U + ,   v c ( 0 ) = V ,   v c ( h ) = V +
Here, U + and U denote the velocities of the upper and lower surfaces in the Xc- direction, respectively; V + and V denote the velocities of the upper and lower surfaces in the Yc-direction, respectively. Although the velocities u c (z) and v c (z) cannot be represented by polynomials exactly [3], here they are approximated and expanded using polynomials of z so that the velocities satisfy Equation (15), as follows. In this case, the variables a 1 and a 2 are unknown.
u c ( z ) = Δ U h · z + U a 1 · z ( h z ) a 2 · z 2 ( h z ) 2
v c ( z ) = Δ V h · z + V
Here, h is the fluid film thickness, and Δ U and Δ V denote the velocity differences; specifically, Δ U = U + U and Δ V = V + V . A higher order term of z , for example, z 3 ( h z ) 3 , can also be considered; however, in this work, the lower-order approximation was chosen. As the pressure gradient toward the Yc direction is zero, v c was assumed to be a Couette flow and approximated considering a linear equation of z. Furthermore, as the pressure gradient toward the Xc-direction is generally non-zero, u c was assumed to be a Poiseuille flow and approximated using a 4th-order polynomial of z. If d = 0 , then u c is also a Couette flow, and Equation (16) can be replaced with the following equation:
u c ( z ) = Δ U h · z + U
The following equations are derived from Equations (16) and (17):
u c z = u ( z ) = Δ U h a 1 · ( h 2 z ) a 2 · 2 z ( h z ) ( h 2 z )
v c z = v ( z ) = Δ V h
The following equations are derived from Equation (19):
u ( 0 ) = Δ U h a 1 h
u ( h ) = Δ U h + a 1 h
u ( h 4 ) = Δ U h a 1 h 2 a 2 3 h 3 16
u ( 3 h 4 ) = Δ U h + a 1 h 2 + a 2 3 h 3 16
If rheological constitutive equations are given, the effective viscosity η * (z) can be calculated using Equations (4)–(6), (16) and (17). The integration of Equation (12) from z = 0 to z = h yields Equation (25), and the integration of Equation (12) from z = h / 4 to z = 3 h / 4 yields Equation (26). These equations are used to determine the values of a 1 and a 2 .
η * ( h ) · ( Δ U h + a 1 h ) η * ( 0 ) · ( Δ U h a 1 h ) = d h
η * ( 3 h 4 ) · ( Δ U h + a 1 h 2 + a 2 3 h 3 16 ) η * ( h 4 ) · ( Δ U h a 1 h 2 a 2 3 h 3 16 ) = d h 2
Equations (21) and (22) show that both η * ( h ) and η * ( 0 ) do not contain a 2 . Therefore, Equation (25) does not contain a 2 and contains only the unknown variable a 1 . The equation can be solved using the one-variable Newton–Raphson method. Although Equation (26) contains both a 1 and a 2 , a 1 has been determined using Equation (25). Consequently, Equation (26) can be considered as an equation involving only the unknown variable a 2 . Thus, it can also be solved using the one-variable Newton–Raphson method. To determine a 1 , a non-dimensional variable b 1 defined using Equation (27) and a function f 1 ( b 1 ) defined using Equation (28) are introduced. The value of b 1 can be calculated considering f 1 ( b 1 ) = 0 .
b 1 = log ( 2 a 1 η n d )
f 1 ( b 1 ) = η * ( h ) · ( Δ U h + e b 1 d h 2 η n ) η * ( 0 ) · ( Δ U h e b 1 d h 2 η n )   d · h
When b 1 is near the solution, Δ b 1 can be calculated using the following equation:
0 = f 1 ( b 1 + Δ b 1 )   f 1 ( b 1 ) + d f 1 d b 1 · Δ b 1     f 1 ( b 1 ) +   e b 1   d h 2 η n   [ η * ( h ) + η * ( 0 ) ] · Δ b 1
In other words, the new candidate b 1 new of b 1 is calculated using the iterative process of Newton–Raphson’s method, as follows:
b 1 new = b 1 + Δ b 1 = b 1 f 1 ( b 1 ) a 1 h   [ η * ( h ) + η * ( 0 ) ]
As η * ( h ) and η * ( 0 ) are originally functions of b 1 , d f 1 / d b 1 includes d η * / d b 1 ; however, in this work, the dependency was ignored, and Δ b 1 was approximated as in Equation (30). To determine a 2 , a non-dimensional variable b 2 defined using Equation (31) and a function f 2 ( b 2 ) defined using Equation (32) are introduced. The value of b 2 can be calculated considering f 2 ( b 2 ) = 0 .
b 2 = 2 a 2 h 2 η n 5 d
f 2 ( b 2 ) = η * ( 3 h 4 ) · ( Δ U h + a 1 h 2 + a 2 3 h 3 16 ) η * ( h 4 ) · ( Δ U h a 1 h 2 a 2 3 h 3 16 ) d h 2
When b 2 is near the solution, Δ b 2 can be calculated using the following equation:
0 = f 2 ( b 2 + Δ b 2 )     f 2 ( b 2 ) + d f 2 d b 2 · Δ b 2     f 2 ( b 2 ) + 15 d h 32 η n [ η * ( 3 h 4 ) + η * ( h 4 ) ] · Δ b 2
In other words, the new candidate b 2 new of b 2 is calculated using the iterative process of Newton–Raphson’s method, as follows:
b 2 new = b 2 + Δ b 2 = b 2 f 2 ( b 2 )   [ η * ( 3 h / 4 ) + η * ( h / 4 ) ] · 32 η n 15 d h
As η * ( 3 h / 4 ) and η * ( h / 4 ) are originally functions of b 2 , d f 2 / d b 2 includes d η * / d b 2 ; however, in this work, the dependency was ignored, and Δ b 2 was approximated as in Equation (34). Subsequently, in the iteration process of Newton–Raphson’s method, only η * depends on the rheological characteristics. Therefore, if the rheological equation corresponding to Equation (5) is incorporated, any isothermal, non-Newtonian EHL calculation can be performed. As per the Newton–Raphson method, the initial value for both b 1 and b 2 can be zero. As both variables a 1 and a 2 are solved using the one-variable Newton–Raphson method, the calculation can be performed within a reasonable time.

4. Calculation of Equivalent Viscosity, Flow, and Surface Force

The flow q 1 along the n 1 direction can be defined using Equation (16), as follows. The density ρ is assumed to be independent of z.
q 1 = 0 h ρ u c   d z = 0 h ρ [ Δ U h · z + U a 1 · z ( h z ) a 2 · z 2 ( h z ) 2 ]   d z = ρ U h ρ a 1 h 3 6 ρ a 2 h 5 30
Here, U is the average velocity in the Xc-direction and can be expressed as follows:
U = U + + U 2
The equivalent viscosity η e q is defined as follows:
η e q = 5 d   2 ( 5 a 1 + a 2 h 2 )
Consequently, q 1 can be represented as
q 1 = ρ U h   ρ h 3 12 η e q · d
The flow q 2 along the n 2 direction can be derived from Equation (17), as follows:
q 2 = 0 h ρ v c   d z = ρ V h
Here, V is the average velocity in the Yc direction and can be expressed as follows:
V = V + + V 2
Hence, in the XYZ coordinate system, the flow vector q can be expressed as
q = ( ρ U h   ρ h 3 12 η e q · d ) · n 1 + ρ V h · n 2 = ρ U h · n 1 + ρ V h · n 2   ρ h 3 12 η e q · d · n 1 = ρ U h ρ h 3 12 η e q · ( P x ,   P y )
Here, U is the average velocity vector of the upper and lower surfaces, defined as follows:
U = U · n 1 + V · n 2 = ( U x ,   U y )
U x and U y denote the average velocities in the upper and lower surfaces in the XY-direction, respectively. When the mass conservation law is applied to Equation (41), the following modified Reynolds equation is obtained.
ρ U x h x + ρ U y h y   x ( ρ h 3 12 η e q · P x )   y ( ρ h 3 12 η e q · P y ) = 0
The difference in the representation of Equations (41) and (43) and that of Newtonian flow only pertains to the viscosities η e q and η n , respectively. In fact, the equivalent viscosity η e q defined by Equation (37) was determined so that Equations (41) and (43) maintain the same form as that of Newtonian flow. Therefore, the EHL calculation procedure for Newtonian flows, such as the method proposed by Venner and Lubrecht [13], can be applied to the current calculation by simply replacing η n with η e q . The shear stress τ 1 along the Xc-direction can be expressed as
τ 1 ( z ) = η *   u c z = η * [ Δ U h a 1 · ( h 2 z ) a 2 · 2 z ( h z ) ( h 2 z ) ]
Therefore, the surface forces P 10 and P 1 h acting on the lower and upper surfaces along the Xc-direction, respectively, can be defined as
P 10 = τ 1 ( 0 ) = η * ( 0 ) · ( Δ U h a 1 h )
P 1 h = τ 1 ( h ) = η * ( h ) · ( Δ U h + a 1 h )
The shear stress τ 2 along the Yc-direction is expressed as
τ 2 = η * v c z = η * Δ V h  
Therefore, the surface forces P 20 and P 2 h acting on the lower and upper surfaces along the Yc-direction, respectively, can be defined as
P 20 = τ 2 ( 0 ) = η * ( 0 ) · Δ V h
P 2 h = τ 2 ( h ) = η * ( h ) · Δ V h
In the XYZ coordinate system, the surface force vectors P 0 and P h that act on the lower and upper surfaces, respectively, can be expressed as
P 0 = P 10 · n 1 + P 20 · n 2
P h = P 1 h · n 1 + P 2 h · n 2

5. Application to a Grease

As mentioned previously, Kochi et al. [1] conducted experiments on grease under soft EHL conditions and measured the film thickness and traction forces. In the present study, the proposed method was applied to one of the greases considered in the study by Kochi et al. [1] so as to validate the theoretical approach. Grease A in the literature [1] was chosen to test the proposed method. The modified Reynolds equation, as expressed in Equation (43), was solved using a multi-level method, as reported by Venner and Lubrecht [13]. The commercial program Tribology Engineering Dynamics Contact Problem Analyzer (TED/CPA) V852 was employed. Figure 3 illustrates the calculation conditions. The upper body was a steel ball, and the lower body was a disk composed of glass or polycarbonate (PC). The rheological properties of the grease were assumed to be represented by Bauer’s model, according to existing literature [1]. Detailed properties of the steel, glass, PC, and grease are described in the previous study [1]. The pressure dependency of the density was defined using Dowson–Higginson’s formula, as follows:
ρ ( P ) = ρ 0 · p 0 + β · P p 0 + P ,   ρ 0 = 1 ,   p 0 = 590   MPa ,   β = 1.34
The pressure dependency of the base oil viscosity was assumed to be defined using Barus’ formula:
η n ( P ) = η 0 · exp ( α P ) ,   η 0 :   49.5   mPa · s ,   α :   14   GPa 1
In particular, Equation (6) diverges when γ ˙ approaches zero. It was assumed that if γ ˙ is lower than a certain value cmin=100.0 s−1, then η * varies linearly with the gradient d η * / d γ ˙ at cmin. Bauer’s parameter k 1 was assumed to be the base oil ambient viscosity η 0 . The values of Bauer’s parameters τ 0 , k 2 , and n were determined from the apparent viscosity curve of grease A, as shown in Figure 9 of Kochi et al. [1]. When P = 0 , the curve was assumed to pass through the following three points: 100 mm/s, 6.46475 Pa·s; 10,000 mm/s, 0.16712 Pa·s; and 1,000,000 mm/s, 0.06573 Pa·s.
The values of τ 0 , k 2 , and n were determined to ensure that Bauer’s curve passes through the abovementioned three points, as follows:
τ 0 = 0.000621839   MPa ,   k 2 = 6.99118 · 10 7 ,   n = 0.7248
Figure 4 presents the central film thickness of grease A and base oil as a function of the rolling velocity in the case of pure rolling and F z = 10   N . The solid lines show the calculation results and the dotted lines show the experimental results. The experimental data were read using the caliper from Figure 4 of Kochi et al. [1]. Figure 4a,b shows the cases of a PC disk and glass disk, respectively. The calculation range was set as follows: −1.0 ≤ X ≤ 0.4 and −0.6 ≤ Y ≤ 0.6.
However, in the case of grease A, a glass disk, and a velocity of less than or equal to 300 mm/s, the range was set as follows: −0.35 ≤ X ≤ 0.14 and −0.2 ≤ Y ≤ 0.2.
In these cases, the oil film was thin and the wide range calculation became hard to perform. In the case of the PC disk, the calculation results exhibited good agreement with the experimental results. In the case of the glass disk, the results of the base oil showed some difference but the other data showed reasonable agreement. Figure 5, which illustrates a sample calculation, shows the distribution of P ,   h , and η e q for the case involving a PC disk, a pure rolling velocity of 1200 mm/s, and grease A. Typically, in the case of pure rolling velocity, γ ˙ is small and η e q is large. In addition, only the appearance of the distribution of η e q differs from that of the base oil. Figure 5d shows the distribution of η e q at section Y = 0. It can be seen that η e q becomes extremely large at the center, where the pressure gradient is small and the flow volume is low.
Figure 6 presents the traction coefficient as a function of the slide roll ratio when F z = 20   N . The solid lines show the calculation results and the dotted lines show the experimental results. The experimental data were read using the caliper from Figure 8 of Kochi et al. [1]. The slide roll ratio is the difference between the upper and lower velocities divided by their average value. The traction coefficient was calculated according to the approaches proposed in the existing literature [1,17]:
T R C = T X 0 T X h 2 F z
Here, T X 0 and T X h denote the X-direction traction forces acting on the lower and upper surfaces, respectively; F z is the load. The calculation range was −0.8 ≤ X ≤ 0.5, −0.5 ≤ Y ≤ 0.5. The calculation results were in fairly good agreement with the experimental results.
Figure 7, which illustrates a sample calculation, shows the distribution of P ,   h , and η e q for the case involving grease A and a slide roll ratio of 10%. Only the appearance of the distribution of η e q differs from that of the base oil, exhibiting a figure eight in the vicinity of the contact point. Figure 7c,d shows the same distribution of η e q in different display ranges. It can be observed that η e q reduces in the rapid flow region. The figure eight phenomenon is characteristic of non-Newtonian flow, in which the apparent viscosity becomes large at a low velocity gradient. This phenomenon can be explained as follows. Let the XY coordinates of points A, B, C, and D be ( 0.1 ,   0 ) , ( + 0.1 ,   0 ) , ( 0 , 0.1 ) , and ( 0 , + 0.1 ) , respectively. The velocity gradient vector n 1 at these points are approximately ( + 1 ,   0 ) , ( 1 ,   0 ) , ( 0 , + 1 ) , and ( 0 , 1 ) , respectively, as shown in Figure 8a. When a 2 in Equation (37) is neglected, η e q can be expressed as follows:
η e q = d   2 a 1  
Calculating a 1 from Equation (25) and substituting it into Equation (56) yields:
η e q = [ η * ( h ) + η * ( 0 ) ]   2 · 1 1 + Δ U d h 2 [ η * ( 0 ) η * ( h ) ]
At points C and D, the direction of flow and that of n 1 is orthogonal, so Δ U = 0 . Consequently, the equivalent viscosity parameters η e q , C and η e q , D at points C and D are given as follows:
η e q , C = η e q , D = [ η * ( h ) + η * ( 0 ) ]   2
At point A, where n 1 directs towards + X and the velocity gradient at Z = h is greater than that at Z = 0 (as shown in Figure 8b), the following equation is satisfied:
η * ( 0 ) > η * ( h ) ,   Δ U = U x + U x > 0
At point B, where n 1 directs towards X , and the velocity gradient at Z = h is smaller than that at Z = 0 (as shown in Figure 8c), the following equation is satisfied:
η * ( 0 ) < η * ( h ) ,   Δ U = ( U x + U x ) < 0
In any case, the equivalent viscosity parameters η e q , A and η e q , B at points A and B become smaller than η e q , C and η e q , D by the effect of the second term of Equation (57). In the pure rolling case, where Δ U is 0, Equation (57) results in Equation (58). It can be understood that in such a case, no figure-eight-shaped distribution appears as is shown in Figure 5c.

6. Conclusions

In this study, an isothermal, non-Newtonian EHL formulation of Bauer’s model was performed using the local coordinate system of the pressure gradient. The flow toward the pressure gradient was assumed to be a Poiseuille flow and was approximated using a 4th-order polynomial of z. The flow along the direction of the zero pressure gradient was assumed to be a Couette flow and was approximated using a linear function of z. The following results were obtained.
(1)
A modified Reynolds equation, which contains an equivalent viscosity, was obtained.
(2)
The EHL calculation procedure for Newtonian flows can be applied to non-Newtonian flows by simply replacing the Newtonian viscosity with the equivalent viscosity.
(3)
If rheological equations are incorporated, any isothermal, non-Newtonian EHL calculation can be performed easily.
(4)
As the equivalent viscosity is calculated using the one-variable Newton–Raphson method, the EHL calculation can be performed within a reasonable calculation time.
(5)
Using Bauer’s model, the formulation was applied to a grease that was evaluated experimentally by Kochi et al. [1]. The results obtained using the proposed method and the experimental results were compared, and reasonable agreement was noted.
(6)
In the case of sliding velocity, the equivalent viscosity shows a figure-eight-shaped distribution in the vicinity of the contact point.
However, the proposed method yields an approximate solution. If the Poiseuille flow and Couette flow cannot be approximated using a 4th-order polynomial of z and a linear function of z, respectively, the obtained results may be inaccurate. The application limits of the current formulation are not clear. Therefore, future work should be focused on determining these limitations.

Funding

This research was funded by the TriboLogics Corporation.

Acknowledgments

The author sincerely thanks TriboLogics Corporation (www.tribology.co.jp/indexEng.htm, accessed on 13 May 2021) for permitting the use of TED/CPA. The author is also very grateful to Hiroyuki OHTA of the Nagaoka University of Technology, Japan, for his constant guidance, encouragement, and valuable suggestions.

Conflicts of Interest

The author declares no conflict of interest.

References

  1. Kochi, T.; Ichimura, R.; Yoshihara, M.; Dong, D.; Kimura, Y. Film Thickness and Traction in Soft EHL with Grease. Tribol. Online 2017, 12, 171–176. [Google Scholar] [CrossRef] [Green Version]
  2. Kauzlarich, J.J.; Greenwood, J.A. Elastohydrodynamic Lubrication with Herschel-Bulkley Model Greases. ASLE Trans. 1972, 15, 269–277. [Google Scholar] [CrossRef]
  3. Conry, T.F.; Wang, S.; Cusano, C. A Reynolds-Eyring Equation for Elastohydrodynamic Lubrication in Line Contacts. J. Tribol. 1987, 109, 648–654. [Google Scholar] [CrossRef]
  4. Dong, D.; Qian, X. A theory of elastohydrodynamic grease-lubricated line contact based on a refined rheological model. Tribol. Int. 1988, 21, 261–267. [Google Scholar] [CrossRef]
  5. Peiran, Y.; Shizhu, W. A Generalized Reynolds Equation for Non-Newtonian Thermal Elastohydrodynamic Lubrication. J. Tribol. 1990, 112, 631–636. [Google Scholar] [CrossRef]
  6. Kim, K.H.; Sadeghi, F. Non-Newtonian Elastohydrodynamic Lubrication of Point Contact. J. Tribol. 1991, 113, 703–711. [Google Scholar] [CrossRef]
  7. Ehret, P.; Dowson, D.; Taylor, C.M. On lubricant transport conditions in elastohydrodynamic conjuctions. Proc. R. Soc. A Math. Phys. Eng. Sci. 1998, 454, 763–787. [Google Scholar] [CrossRef]
  8. Greenwood, J.A. Discussion: Two-dimensional flow of a non-Newtonian lubricant. Proc. Inst. Mech. Eng. Part J J. Eng. Tribol. 2001, 215, 121–122. [Google Scholar] [CrossRef]
  9. Sharif, K.J.; Kong, S.; Evans, H.P.; Snidle, R.W. Contact and elastohydrodynamic analysis of worm gears Part 1: Theoretical formulation. Proc. Inst. Mech. Eng. Part C J. Mech. Eng. Sci. 2001, 215, 817–830. [Google Scholar] [CrossRef] [Green Version]
  10. Liu, X.; Jiang, M.; Yang, P.; Kaneta, M. Non-Newtonian Thermal Analyses of Point EHL Contacts Using the Eyring Model. J. Tribol. 2005, 127, 70–81. [Google Scholar] [CrossRef]
  11. Yang, Q.; Huang, P.; Fang, Y. A novel Reynolds equation of non-Newtonian fluid for lubrication simulation. Tribol. Int. 2016, 94, 458–463. [Google Scholar] [CrossRef]
  12. Bordenet, L.; Dalmaz, G.; Chaomleffel, J.-P.; Vergne, F. A study of grease film thicknesses in elastorheodynamic rolling point contacts. Lubr. Sci. 1990, 2, 273–284. [Google Scholar] [CrossRef]
  13. Venner, C.; Lubrecht, A. Foreword. Veh. Tribol. 2000, 37, 7–8. [Google Scholar] [CrossRef]
  14. Bauer, W.H.; Finkelstein, A.P.; Wiberley, S.E. Flow Properties of Lithium Stearate-Oil Model Greases as Functions of Soap Concentration and Temperature. ASLE Trans. 1960, 3, 215–224. [Google Scholar] [CrossRef]
  15. Palacios, J.; Palacios, M. Rheological properties of greases in ehd contacts. Tribol. Int. 1984, 17, 167–171. [Google Scholar] [CrossRef]
  16. Johnson, K.L.; Tevaarwerk, J.L. Shear behaviour of elastohydrodynamic oil films. Proc. R. Soc. Lond. Ser. A Math. Phys. Sci. 1977, 356, 215–236. [Google Scholar] [CrossRef]
  17. De Vicente, J.; Stokes, J.R.; Spikes, H.A. Rolling and sliding friction in compliant, lubricated contact. Proc. Inst. Mech. Eng. Part J J. Eng. Tribol. 2006, 220, 55–63. [Google Scholar] [CrossRef]
Figure 1. Global coordinate system O, X, Y, Z.
Figure 1. Global coordinate system O, X, Y, Z.
Lubricants 09 00056 g001
Figure 2. Local coordinate system Oc, Xc, Yc.
Figure 2. Local coordinate system Oc, Xc, Yc.
Lubricants 09 00056 g002
Figure 3. Calculation condition.
Figure 3. Calculation condition.
Lubricants 09 00056 g003
Figure 4. Comparison of film thickness between the calculation and the experiment.
Figure 4. Comparison of film thickness between the calculation and the experiment.
Lubricants 09 00056 g004
Figure 5. Distribution of P, h, and η e q of grease A at a rolling velocity of 1200 mm/s.
Figure 5. Distribution of P, h, and η e q of grease A at a rolling velocity of 1200 mm/s.
Lubricants 09 00056 g005
Figure 6. Comparison of the traction coefficient between the calculation and the experiment.
Figure 6. Comparison of the traction coefficient between the calculation and the experiment.
Lubricants 09 00056 g006
Figure 7. Calculation results of grease A at a 10% slide ratio.
Figure 7. Calculation results of grease A at a 10% slide ratio.
Lubricants 09 00056 g007
Figure 8. Explanation of the figure eight phenomenon.
Figure 8. Explanation of the figure eight phenomenon.
Lubricants 09 00056 g008
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Kakoi, K. Formulation to Calculate Isothermal, Non-Newtonian Elastohydrodynamic Lubrication Problems Using a Pressure Gradient Coordinate System and Its Verification by an Experimental Grease. Lubricants 2021, 9, 56. https://0-doi-org.brum.beds.ac.uk/10.3390/lubricants9050056

AMA Style

Kakoi K. Formulation to Calculate Isothermal, Non-Newtonian Elastohydrodynamic Lubrication Problems Using a Pressure Gradient Coordinate System and Its Verification by an Experimental Grease. Lubricants. 2021; 9(5):56. https://0-doi-org.brum.beds.ac.uk/10.3390/lubricants9050056

Chicago/Turabian Style

Kakoi, Kunihiko. 2021. "Formulation to Calculate Isothermal, Non-Newtonian Elastohydrodynamic Lubrication Problems Using a Pressure Gradient Coordinate System and Its Verification by an Experimental Grease" Lubricants 9, no. 5: 56. https://0-doi-org.brum.beds.ac.uk/10.3390/lubricants9050056

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop