Effect of fluid inertia force on thermal elastohydrodynamic lubrication of elliptic contact

A non-Newtonian thermal elastohydrodynamic lubrication (TEHL) model for the elliptic contact is established, into which the inertia forces of the lubricant is incorporated. In doing so, the film pressure and film temperature are solved using the associated equations. Meanwhile, the elastic deformation is calculated with the discrete convolution and fast Fourier transform (DC-FFT) method. A film thickness experiment is conducted to validate the TEHL model considering the inertia forces. Further, effects of the inertia forces on the TEHL performances are studied at different operation conditions. The results show that when the inertia forces are considered, the central and minimum film thicknesses increase and film temperature near the inlet increases obviously. Moreover, the inertial solution of the central film thickness is closer to the experimental result compared with its inertialess value.


Introduction
As is well known, the usual Reynolds equation is derived from the Navier-Stokes equations with several simplified assumptions including negligible inertia force of the lubricant. This dealing is theoretically acceptable in the hydrodynamic lubrication analysis for the mechanical component operating at small values of Reynolds number of the lubricant [1]. Based on this assumption, the Newtonian or non-Newtonian thermal elastohydrodynamic lubrication (EHL) performances of point contacts such as ball bearings have been extensively studied theoretically and experimentally in the past years. In the case of the Newtonian thermal EHL (TEHL) study, Lee et al. [2] found that the maximum film temperature of circular contact increases with increasing entrainment velocity of the lubricant and applied load. Guo et al. [3] studied the Newtonian TEHL properties of two elliptic contact surfaces moving in opposite directions and found that there exists a dimple of the film profile in the central contact zone. Later, Liu and Yang [4] investigated the effect of the contacting body temperature on the circularcontact Newtonian TEHL performances, whose results indicated that this temperature can change the film thickness and film temperature at the low entrainment velocity of the lubricant.
At heavy-load and high-speed conditions, however, the conventional Newtonian EHL model perhaps overestimates the film pressure and film temperature due to the neglect of shear-thinning effect of the lubricant [5][6][7]. In the early stage, the point contact lubrication experiments were extensively carried out to investigate the shear-thinning phenomenon by some researchers. For example, Hirst and Moore [8] conducted the two-disc lubrication experiment to determine the relation between the shear stress and shear rate of the EHL film, and introduced the Eyring theory into the non-Newtonian lubrication model. Later, Johnson and Tevaarwerk [9] proposed a simple constitutive equation to reveal the lubricant shearing behavior based on the Eyring sinh law, which was verified with the disc machine experiment. Based on the above Ree-Eyring non-Newtonian model, the TEHL theoretical analyses of the point contact were also carried out in recent years. For example, Kaneta et al. [10] discussed the effect of a transverse ridge and groove on the elliptic contact non-Newtonian TEHL performances, and pointed out that the movement velocity and appearance time of a perturbed film are influenced mainly by the slide-roll ratio. Cui et al. [11] built a non-Newtonian TEHL model for point contacts with transverse ridges existing on contacting surface, and observed that the shear-thinning effect of the lubricant is significant and its thermal effect is weakened under the small contact radius. Recently, Xiao and Shi [12] studied the stiffness and damping performances of the point contact TEHL for the crowned herringbone gear based on the Ree-Eyring model, and found that the crown modification can significantly affect the stiffness and damping of the lubrication film.
The above studies focused on the point-contact TEHL problem, in which the fluid inertia effect is neglected. For this dealing, therefore, its numerical accuracy will drop off under operation conditions with the large Reynolds number of the lubricant. In practices some mechanical parts, such as bearings, gears, and cam-follower systems often operate at high speed conditions. In this case, the fluid inertia effect is obvious, so it is necessary to consider the fluid inertia effect in the lubrication analysis of the mechanical parts. The early researches about the inertia effect of lubricant on the EHL performances are mainly aimed at the journal bearing. Shi and Wang [13] performed a TEHL analysis of the water-based ferro-fluid lubricated sliding bearing considering the fluid inertia force, and their results show that the film thickness increases and the pressure peak decreases with the inertia force considered. Dong et al. [14] applied the inertia terms-incorporated Reynolds equation to study the EHL properties of waterlubricated journal bearings. It was found that when the fluid inertia force is considered, the film pressure at the inlet region increases, the pressure peak slightly decreases, and both the central and minimum film thicknesses increase. The TEHL analysis of the journal bearing from Liu [15] indicated that when the inertia force is considered, the film thickness increases and the maximum film temperature is reduced by about 7%. Fan [16] explored the effect of inertia force on the TEHL performances of the water-lubricated ceramic journal bearing and found that the film necking phenomenon is more obvious with the inertia force considered. Beside these, Lin et al. [17] studied the linecontact TEHL problem, and pointed out that compared with the inertialess situation, the film thickness increases slightly and the maximum mid-layer film temperature decreases by less than 2% under the inertia case.
In spite of the aforementioned studies, only the Newtonian and non-Newtonian EHL behaviors ignoring the inertia effect, along with the inertial EHL behaviors of the journal bearing and line contact, were investigated in the past years. For some mechanical parts with the elliptic contact, whose Reynolds number of the lubricant is large, the inertia force of the lubricant may has the effect on their lubrication performances. Therefore, the inertia effect should be considered in the lubrication of these mechanical parts. The past EHL studies considering the fluid inertia effect, however, aimed for either the journal bearing or the line contact, which is not applicable to the elliptic contact. In the present study, an elliptic contact non-Newtonian TEHL model considering the inertia force of the lubricant is established. In doing so, the non-Newtonian Reynolds equation considering the fluid inertia force is derived and solved using the chase-after method, which is also applied to uniformly solve the energy and heat conduction equations [4]. Besides, the DC-FFT method [18,19] is used to accelerate the elastic deformation computation of the contacting solids. Then, a film thickness test is carried out to verify the rationality of the established model based on the optical interference principle. Based on the above non-Newtonian TEHL model considering the fluid inertia force, the lubrication performances without and with the inertia effect are compared at the varied applied loads, entrainment velocities, slip-roll ratios and environment temperatures. Finally, associated conclusions are obtained.

Governing equations
The TEHL solution of the elliptic contact considering the inertia force of the lubricant is based on the elliptic contact model. As shown in Figure 1a, the elliptic contact is simplified so as to achieve the match between an elastic ellipsoid B with an infinitely flat plate A. The contact center of the plate and ellipsoid is located at the coordinate origin o. Symbols x and y are separately the coordinates along and perpendicular to the lubricant's flow direction, and z is the coordinate along the film thickness direction. u a and u b are the velocities of the lower and upper surfaces in the x-direction. Q denotes the load applied on the ellipsoid B. Moreover, the calculation domain of the ellipsoid-plate contact is illustrated in Figure 1b, in which x in and x out separately represent the inlet and outlet positions of the lubricant in the x-direction, and y in and y out separately denote the edges of the solution domain along the y-direction. a and b separately denote the semi-length and semi-width of the contact ellipse in the x-and y-directions, in which case the ellipticity ratio k is defined as k ¼ b=a.
The present study aims at exploring the non-Newtonian steady-state TEHL performances of the elliptic contact under the consideration of the inertia forces of the lubricant. In doing so, the lubricant is assumed to be continuous, isotropic, and its body force and variation in the film pressure across film thickness are ignored. In this case, the steady Navier-Stokes equations used in the analysis can be simplified as Here, I x and I y separately denote the inertia forces of any fluid cell in the x-and y-directions, which are labelled in Figure 1, and defined as I x =r(u ∂ u/∂ x +v ∂ u/∂ y +w ∂ u/∂ z) and . For other variables, p and h separately denote the film pressure and film thickness, r and h* separately denote the density and effective viscosity for the non-Newtonian lubricant, u and v separately represent the velocity components of the lubricant in the x-and y-directions. In addition, the lubricant velocity in the z-direction w can be obtained according to the continues equation, that is, Integrating equation (1) twice with respect to z, the velocity components of the lubricant are The mass flux terms in the x-and y-directions are written as Substituting equation (2) into equation (3) yields where the entrainment velocity is defined as u r = (u a + u b )/2.
Besides, the detailed expressions of r=h ð Þ e , h e , h 0 e , r e , r 0 e and r 00 e are given by Yang and Wen [20,21].
According to the mass conservation law, the continuity equation of the lubricant is given as Integrating equation (5) with respect to z over the intervals 0 to h, using the Newton-Leibnitz integral rule, and substituting equation (4) into the equation (5), the generalized Reynolds equation including the inertia terms in steady-state elliptic contact TEHL can be obtained as To solve equation (6) for the film pressure, the following boundary conditions must be satisfied The film thickness is written as Here, h 0 is the normal approach of the ellipsoid surface and flat plate, whose equivalent curvature radii in the xoz and yoz planes are denoted by R x and R y , respectively. The composite elastic deformation of the contacting solids can be expressed by the Boussinesq integration where E 0 is the composite elastic modulus with an expression of 2=E 0 ¼ ð1 À y 2 a Þ=E a þ ð1 À y 2 b Þ=E b . E a and E b are the elastic moduli of the plate A and ellipsoid B, whose Poisson's ratios are y a and y b , respectively.
Directly solving equation (9) for the elastic deformation requires a quadruple loop in programming, in which case the computation is usually time-consuming. To save the computation time, the DC-FFT method employed by Meng et al. [18,19] and Liu et al. [22,23] is applied here.
The density-pressure-temperature relation proposed by Dowson and Higginson is employed as [24] r ¼ r 0 1þ 0:6 Â 10 À9 p 1 þ 1:7 Â 10 À9 p where r 0 is the ambient lubricant density, T 0 is the reference temperature, and T f is the film temperature. Based on the above equations, the lubricant velocities u along the x-direction and v along the y-direction can be obtained through equation (2), and the lubricant density r can be obtained through equation (10). Combining the expression of the z-direction fluid velocity w, the inertia forces I x and I y can be solved according to their expressions. It is noteworthy that in solving the lubricant velocities u and v, the inertia-related terms I x , I y , I x 0 , and I y 0 need to be obtained in advance, which means the solution of lubricant velocities and inertia forces are nested with each other. In the present study, the initial lubricant velocities are deemed as ones obtained by the TEHL analysis without the inertia effect of the fluid under the same conditions.
To obtain the viscosity of the non-Newtonian lubricant, one of Newtonian lubricant is obtained in advance, whose viscosity-pressure-temperature relation is evaluated according to the Roelands equation [24], that is, where h 0 is the ambient lubricant viscosity, Z 0 = a/[5.1Â10 À9 (lnh 0 +9.67)], and S 0 = b(T 0 À 138)/(lnh 0 + 9.67). Symbols a and b are the viscosity-pressure coefficient and viscosity-temperature coefficient of the lubricant, respectively. For the Ree-Eyring non-Newtonian lubricant applied in the present study, its effective viscosity h * is defined as where t 0 is the characteristic shear stress of the lubricant. t e is the module of the shear stress vector with an , wherein t zx and t zy are separately the shear stresses of the film in the x-and y-directions. The applied load Q is balanced through the pressure integration in the solution domain V Regardless of effects of body force and heat radiation, and ignoring the heat conduction in the x-and y-directions, the energy equation of the lubricant is used to obtain film temperature T f , which is given below where c f and k f are the specific heat and thermal conductivity of the lubricant, respectively. The first and second terms on the left hand of the equation (14) separately represent the convective heat dissipation and compression work. The two terms on the right hand of the equation (14) successively stand for heat conduction and shearing heat. Since the energy equation ignores the heat conduction in the x-direction, the convection is the only heat transfer mode in the x-direction. Thus, at the lubricant inlet, no temperature boundary condition is required for the countercurrent region exists, but T f = T 0 for the downstream region. Similarly, temperature boundary conditions are not required for the film outlet and two longitudinal sides.
The heat conduction equations of the plate A and ellipsoid B are given by where c a (c b ), r a (r b ), and k a (k b ) separately stand for the specific heat, density, and thermal conductivity of solid A (B). T a and T b separately denote the temperatures of solids A and B, which, along with the film temperature T f , are solved simultaneously through equations (14) and (15) based on the following boundary conditions that are equations (16)- (18). For the plate A and ellipsoid B, the following temperature boundary conditions in z-direction must be satisfied Here, T 0 is the ambient temperature and h t is the heat permeating thickness of the contacting solids.
At the starting position along the x-direction of the contacting solids, the temperature boundary conditions are given as The interfacial heat flow continuity conditions for the plate A, the film, and the ellipsoid B are given by k Treating the temperature of the plate A, the film, and the ellipsoid B as a whole, the solution domain of the temperature to be solved is By solving the above equations simultaneously, the non-Newtonian TEHL performances of elliptic contacts without or with inertia force considered can be obtained.
In what follows, the three kinds of temperatures are all denoted by symbol T to be convenient for the analysis.

Solution scheme
The above governing equations are dimensional. In order to reduce the numerical error and make the numerical results generalized, the dimensional governing equations need to be nondimensionalized with the following dimensionless parameters: the dimensionless coordinates x ¼ x=a, y ¼ y=b, and z ¼ z=h. In detail, z ¼ z=h aims for the plate A (lower solid) when Àh t z 0, for the lubricant film when 0 z h, and for the ellipsoid B (upper solid) when h z h + h t . Here, h t is equal to 3.15a [20], in which the boundary locations are regarded as sufficiently far from the oil-solid interfaces. Below this value, the solid temperature is not change along the film thickness direction. The dimensionless film thickness h ¼ h=a, the dimensionless film pressure p ¼ p=p H , and the maximum Hertzian pressure p H = 3Q/(2pab). The dimensionless film velocities in the x-and y-directions are u ¼ u=u r and v ¼ v=u r , respectively. Besides, the dimensionless temperature is T ¼ T =T 0 , which represents the temperature distribution of the lubricant film or the contacting bodies.
After nondimensionalizing the governing equations based on the above dimensionless parameters, the columnby-column is used to solve the Reynolds equation to obtain the film pressure [11], and the same method is applied for solving the energy equation of the lubrication film and the heat conduction equations of the contact solids to obtain the temperature field. The composite elastic deformation can be separately solved using the DC-FFT method. The numerical calculation is conducted in the solution domain of À4:5 x 1:5 and À1:8 y 1:8. Such a large solution domain can avoid numerically the starved lubrication as pointed out by Liu et al. [25].
Too sparse a grid system results in an inaccurate evaluation of TEHL performances, while too dense a grid system brings about a long computation time. In the present study, the grid number is NX = 128 in the x-direction and NY = 512 in the y-direction. In this case, grid intervals in x-and y-directions are almost the same, which contributes to the fast and accurate calculation, as pointed out by Kaneta et al. [10]. In the z-direction, the grid number across the film is NZ = 10, NZA = 6 within the plate A, and NZB = 6 within the ellipsoid B. From the solid-liquid interface to the interior of the upper or lower solid, the grid interval between two closely neighboring nodes increases by two times with increasing serial number of the grid. The starting grid interval of the solid is set as 0.05a. Such a grid distribution can guarantee to achieve accurate and efficient convergent solutions for the film pressure and temperature. The numerical calculation flowchart of the elliptic-contact TEHL considering the inertial effect is shown in Figure 2. The specific calculation process is described as follows.
Step 1: Give the initial film pressure p, film temperature T and velocity components (u and v), which can be derived from TEHL analysis ignoring the inertia force using the same input parameters.
Step 2: Calculate the lubricant density r and its viscosity h according to equations (10) and (11), respectively, and further obtain the film thickness h considering the elastic deformation based on equation (8).
Step 3: Calculate the coefficient matrices of the Reynolds equation (6), including terms r=h ð Þ e , r Ã x , V x , and V y .
Step 4: Calculate the lubricant velocities u and v along the x-and y-directions using equation (2).
Step 5: Determine whether the fluid velocities converge according to the following convergence criteria.
where the superscript "new" denotes the current iteration, and "old " indicates the previous iteration. e u and e v separately represent the convergence precisions for the lubricant velocities in the x-and y-directions, whose values are both taken as 1.0 Â 10 À4 . If the velocities not satisfy the above criteria, update them with the relaxation factor of 0.0005, and then return to Step 4. It is notable that steps 3-5 constitute the velocity iteration process considering the lubricant inertia forces, as shown in the dotted box of Figure 2.
Step 6: If the lubricant velocities convergence accuracy is met, calculate the film pressure p based on the Reynolds equation expressed in equation (6), and then the film temperature T based on the energy equation of the lubrication film expressed in equation (14) and heat conduction equations of the solids expressed in equation (15). In doing so, the corresponding boundary conditions such as ones in equations (16)- (19) are used.
Step 7: Check whether the film pressure and temperature converge simultaneously. If the pressure does not satisfy the convergence criterion in equation (21) or the temperature does not satisfy the convergence criterion in equation (22), the pressure and temperature are updated with the relaxation factors of 0.01 and 0.05, respectively, and then return to Step 2.
In equations (21) and (22), the pressure convergence precision e p is taken to be 1.0 Â 10 À4 , and temperature convergence precision e T is set as 1.0 Â 10 À5 in the present study. Such convergence precisions guarantee the computation accuracy.
Step 8: Further determine whether the load is balanced with the following convergence criterion.
Here, Q T is the dimensionless load-carrying capacity, which is obtained through the integration of the dimensionless film pressure over the whole computation solution.
Q R is the dimensionless referenced applied load (Q R ¼ 2p=3). e Q is the load convergence precision, set as 1.0 Â 10 À3 . If the load-carrying capacity of the lubricant is balanced by the applied load, terminate the whole calculation. Otherwise, update the normal approach h 0 in equation (8) with the relaxation factor of 0.003 and continue the calculation by returning to Step 2.

Experimental verification
To validate the correctness of the elliptic-contact TEHL model considering the inertia force, a central film thickness experiment is carried out with the TFM-150 film thickness measuring apparatus with a repeatability error less than 0.5%. The applied load and entrainment velocity ranges of the apparatus are 0-100 N and 0-10 m/s, respectively, which meet the requirements of the experiment. Figure 3, the central film thickness between a disc and a ball is measured in the experiment. The disc is made of K9 optical glass, whose elastic modulus, Poisson's ratio, density, specific heat and thermal conductivity are separately 81 GPa, 0.208, 2510 kg/m 3 , 840 J/(kg°C) and 1.11 W/(m°C). The ball with the radius of 22.225 mm is made of SUJ2 steel, whose elastic modulus, Poisson's ratio, density, specific heat and thermal conductivity are separately 207 GPa, 0.29, 7810 kg/m 3 , 533 J/(kg°C) and 40.11 W/(m°C). The lubricant used is 4050 synthetic aviation lubricant with the density of 971.2 kg/m 3 , viscosity of 0.0165 Pa · s, viscosity-pressure coefficient of 1.88 Â 10 À8 Pa À1 and viscosity-temperature coefficient of 0.035°C À1 . Other parameters for this oil are as follows: specific heat of 1910 J/(kg°C), thermal conductivity of 0.152 W/(m°C), and characteristic shear stress of 10 MPa.

As illustrated in
The feeding temperature of the lubricant is kept constant as 50°C with an oil heating tank. To reduce the repeatability error, the experiment under each working condition is repeated three times. Unlike the elliptic contacts in Section 2, the ball-on-disc contact problem is a circular-contact problem. Therefore, in the corresponding simulation in this section, the solution domain is defined as À2:5 x 1:5 and À2:0 y 2:0, the grid numbers NX and NY are separately set as 256 and 128, the temperature convergence precision e T is taken to be 1.0 Â 10 À4 , and other input parameters are the same as those in Section 2. Figure 4 compares the central film thickness obtained by the experiment and the simulation at the varied entrainment velocity. In the comparison, the applied loads are 30 and 60 N and slide-roll ratio s = 0.3. Moreover, the simulation parameters for this comparison are consistent with those of the steel ball-glass disc contact pair in the experiment. As seen from Figure 4, the greater the entrainment velocity, the larger the central film thickness, which is due to the enhanced hydrodynamic effect of the lubricant. More importantly, compared with the central film thickness under the inertialess situation, the central film thickness under the inertial situation is closer to the experimental result. The maximum error between the inertia and inertial solutions is 7.76% in Figure 4a with Q = 30 N, while it is 3.82% for Figure 4b with Q = 60 N. The above discrepancies may be caused by the impurities of the lubricant and wear debris of the ball disk. These influencing factors will cause a slight reduction in the refractive index of the lubricant, in which case the measured central film thickness is slightly higher than the simulation result. However, the inertial solution is in good agreement with the experimental data in whole, indicating that consideration of the inertia force in the non-Newtonian TEHL analysis is reasonable.

Results and discussion
Based on the above-mentioned TEHL model considering the inertia force effect, working conditions on the elliptic contact non-Newtonian TEHL performances are studied. In the numerical calculation, FAG B7014AC angular contact ball bearing with an ellipticity ratio of 5.23 is used. This ratio is far less than 8.0, and thus the obtained result is different from the line contact case [26]. The materials of solids A and B are both GCr15, whose material and thermal properties are listed in Table 1.
To form the lubrication film in the elliptic contact region, the high viscous lubricant PAO6 [27] is applied, whose rheological properties along with the working condition are listed in Table 2. In the subsequent simulation, these parameters are unchanged until specified.

At varied applied load
The elliptic contact TEHL performances considering the inertial force are first studied at different applied loads Q. The values of Q are taken to be 65, 104 and 155 N in turn, whose corresponding contact pressure ranges from 1.2, 1.4 and 1.6 GPa. Owing to the symmetry of film pressure and film thickness along the y-axis, their profiles along the xaxis only are discussed. Figure 5a compares the film pressure without and with inertia forces. Regardless of the inertia or inertialess effect, the larger the applied load, the greater the peak pressure, and the closer the secondary pressure peak is to the outlet. With the inertia forces considered, the pressure peak remains little changed. Figure 5b presents the load effect on the film thickness with and without the inertia forces. As seen from this figure, the larger the applied load, the thinner the film thickness, and the closer the location of the minimum film thickness is to the outlet. Moreover, the inertia force makes the central and minimum film thicknesses increase.
The above film thickness change is further illustrated in Figure 6. At the applied loads of 65, 104 and 155 N, the central film thicknesses h c are increased separately by 2.42%, 5.14% and 5.20%, and the minimum film thicknesses h min are increased by 1.96%, 4.08% and 4.46%. Such increases are due to the phenomenon that the inertia force can promote the centrifugal motion of the lubricant, thus the film thickness increasing. Figure 7 gives a comparison between the representative mid-layer film temperatures T with and without inertia force effects at the varied applied load Q. The temperatures are taken at the middle cross section of y = 0. As seen from Figure 7a, regardless of the inertia or inertialess effect,     increasing the load results in a rise in the maximum film temperature. This is due mainly to the phenomenon that at the relative large load, the more work is done by the loadcompressed lubrication film. It should be pointed out that the maximum film temperature position depends on the operation conditions. It sometime appears near the inlet due to the phenomenon that there may exist the relative large film pressure at this position [28,29].
When the inertia forces are considered, the maximum film temperature in the contact center increases slightly. The unobvious temperature rises are 0.52%, 0.56% and 0.42% at Q = 65, 104 and 155 N, respectively. A similar rule was also discovered in the journal bearing work by Liu [15]. When the inertia forces are considered, however, the temperature rise zone expands towards the inlet. Near the inlet, the obvious temperature rise can be found, as shown in Figure 7b. An illustrative temperature rise up to 9.9% appears at the position of x/a = À3.375 when Q = 104 N.
The above phenomena can be explained through the illustrative Figure 8, which compares the lubricant velocities in the x-direction for the inertialess and inertial situations at Q = 104 N. In Figure 8, the negative velocity means the occurrence of the lubricant countercurrent. The comparison between countercurrents in Figures 8a and 8b shows that the range of the countercurrent region in the inertia case is obviously larger and closer to the inlet, compared with that in the inertialess case. This shows that the inertia force can argument the lubricant countercurrent to a degree. The augmented countercurrent resists the free flow of the lubricant in the contact region and make the molecular interaction between different lubricant layers stronger, thus the inlet temperature increasing obviously. Figure 9 gives the friction coefficients m f with and without inertia force effects at different applied loads. With the increase of the applied loads Q, the lubricant viscosity and the shear stress of the lubricant are enhanced, which cause the increment of the friction force and friction coefficient. Meanwhile, as seen from Figure 9, the friction coefficient with the inertial effect is obviously lower than that without the inertialess effect at the same load. The reason is that the film thickness increases and the film's shear stress is weakened when the inertial forces are considered, which brings out reductions of the friction force and friction coefficient.

At varied entrainment velocity
The entrainment velocity of the lubricant is another significant factor affecting the elliptic contact TEHL performances. The entrainment velocity u r used varies from 2 to 11 m/s, which corresponds to the modified Reynolds number Re Ã ranging from 1.48 to 16.31. Figure 10 illustrates the film thickness ignoring and considering the inertia forces at the varied entrainment velocity. As shown in Figure 10a, the larger the entrainment velocity u r , the thicker the film thickness and the narrower the flat segment of the film thickness curve. As shown in Figure 10b, the central and minimum film thicknesses with the inertia force are larger than the corresponding values without the inertia force. Besides, the maximum relative percentage increases in the central film thickness h c at u r = 2, 5, 8 and 11 m/s are separately 4.95%,  1.55%, 0.28% and 0.65%, while increases in the minimum film thickness h min for the four velocities are separately 4.35%, 1.34%, 0.21% and 1.19%. With the increase in the value of u r , the Reynolds number Re of the lubricant linearly increases according to the expression Re = r 0 u r a/h 0 . Further, as Re increases, the film pressure nonlinearly varies due to the two-order partial differential terms on the right hand of the Reynolds equation, which results in the nonlinear variations of the film thickness. Figure 11a compares the mid-layer film temperature T at different entrainment velocities u r. When the velocity increases, the maximum film temperature increases due to the enhanced hydrodynamic action of the lubricant. Moreover, the film temperatures increases when the inertia force is considered. At u r = 2, 5, 8 and 11 m/s, the increments in the maximum film temperature at the contact center are separately 2.25%, 0.56%, 0.08% and 1.31%, and the maximum relative percentages of the temperature rise near the inlet are separately up to 5.01%, 8.46%, 5.12% and 2.85%, as shown in Figure 11b.
The above temperature rise can be explained through the variation in velocity fields (z = h/2) with and without inertia forces. As shown in Figure 12, when inertia force is considered, the amplitude of the lubricant velocity increases to a degree and the countercurrent velocity becomes stronger. Moreover, the reflow region of the inertia force is closer to the inlet. The stronger countercurrent arguments the interaction between molecules at different fluid layers; therefore, the inlet temperature increases. Figure 13 shows the effect of entrainment velocity u r on the friction coefficient u f without and with the inertial forces of lubricant. With the increment in the entrainment velocity, the enhanced film temperature appears due to the strong hydrodynamic action of the lubricant. In this case, the lubricant viscosity decreases so that the friction force and friction coefficient decrease due to the weakened shear stress of the lubrication film. Moreover, the friction coefficient under the inertial situation is lower than the inertialess solution at the same entrainment velocity. The reduction is up to about 4% at u r = 9 m/s. Figure 14 shows the film thicknesses h with and without the inertia forces at the varied slide-roll ratio s. The slip ratio s used varies from 0.05 to 1.8.

At varied slide-roll ratio
As shown in Figure 14, the large slide-roll ratio leads to a thin film thickness. As the slide-roll ratio increases, both  the shear heating and shear thinning of the lubricant are enhanced, which reduces the lubricant viscosity and eventually leads to a decrease of the film thickness.
Another phenomenon observed is that in the case of relatively large slide-roll ratio (e.g., s = 1.5), the film thickness curve inclines towards the outlet shown in Figure 14a. In the case of small slide-roll ratios (s = 0.02-0.5), this curve basically maintains flat on the bottom, as shown in Figure 14b. Similar inclination phenomenon can also be found in literature [20]. The larger slide-roll ratio implies larger relative sliding velocity between the matching surfaces, in which case the higher film temperature rise occurs at the outlet. Due to the thermal-thinning effect, the higher film temperature induces a lower lubricant viscosity at the outlet and weakens the film's ability to resist the external applied load. In order to balance the applied load, the film thickness at the outlet drops naturally and the film thickness curve becomes inclined towards the outlet. Thus, the larger slide-roll ratio causes a slight inclination due to the more significant thermal-thinning effect.
It can be also found from Figure 14 that when the inertia force is considered, both the central and minimum film thicknesses increase. Moreover, this increment becomes obvious at the small value of slide-roll ratio, in which case the relative roll state of matching bodies is dominate. In Figure 14a with the relative large value of s, the central film thicknesses at s = 0.5, 1.0, 1.5 and 1.8 are separately increased by 4.11%, 2.60%, 2.12% and 1.93%, and the corresponding minimum film thicknesses are separately increased by 6.62%, 3.59%, 2.51% and 2.02%. In Figure 14b with the relative small value of s, the maximum percentage increment of the central film thicknesses at s = 0.02,0.05 and 0.08 are found as 5.31%, 5.14% and 4.66%, while their corresponding values are separately 4.13%, 4.08% and 3.72% for the minimum film thicknesses. Figure 15 shows the comparisons of the maximum midlayer film temperatures at the contact center with and without inertia forces. When the inertia force is considered,  the maximum film temperature increases. Such a temperature rise is more obvious at the small slide-roll ratio s, as illustrated in Figure 15a. The small slide-roll ratio means smaller relative sliding velocity between contacting bodies. In this case, the shearing heat caused by the inertia force of the lubrication film is not easily to be taken away from the contact region, thus the temperature rise near the inlet position becomes more obvious. The temperature rise discrepancy illustrated above shows that the inertia effect is significant in the case of small slide-roll ratio. Figure 16a shows the whole comparison of the mid-layer film temperature T with and without the inertia forces. The maximum film temperature increases slightly, but the obvious temperature rise near the inlet can be found.    Figure 16b shows the maximum relative percentages of the temperature rise near the inlet, whose values are separately up to 5.43%, 6.37%, 5.77% and 4.96% at s = 0.5,1.0,1.5 and 1.8. Such a nonlinear film temperature increase is due to the two-order partial differential terms on the right hand of the Reynolds equation. Figure 17 compares the friction coefficients u f at the varied slide-roll ratio s. With the increase of the slide-roll ratio s, the film temperature is enhanced, and the lubricant viscosity and the shear stress of the lubrication film decrease. In this case, both the friction force and friction coefficient decrease. Moreover, the friction coefficient under the inertial situation is lower than that under the inertialess case at the same slide-roll ratio. The reduction is up to about 3.4% at s = 0.4 due to the weakened shear stress of lubricant.

At varied environment temperature
The environment temperature affects theological properties of the lubricant and further elliptic contact TEHL performances. In the present section, the environment temperatures T 0 are separately set as 25, 35 and 45°C, whose corresponding viscosity are separately 0.048, 0.029 and 0.018 Pa · s, and corresponding densities are separately 866.0, 860.6 and 855.3 kg/m 3 . Figure 18a, the higher environment temperature brings out flatter film thickness at the contact region. This is due to the phenomenon that the higher temperature leads to lower lubricant viscosity, and thus more lubricant is squeezed out of the contact region. When the environment temperature increases from 25 to 45°C, the inertial value of the film thickness is larger than its inertialess solution. As shown in Figure 18b, the maximum relative percentage increases in the central film thickness h c are 1.55%, 2.93% and 3.97% in turn, while such increases in the minimum film thickness h min are separately 1.34%, 2.56% and 3.56%. Figure 19 shows the rise in the mid-layer film temperature DT at the varied environment temperature T 0 . Figure 19a shows that increasing the value of T 0 leads to a decline in maximum film temperature due to the weakened hydrodynamic action of the lubricant. When the inertia force is considered, the maximum film temperatures at T 0 = 25, 35 and 45°C increase separately by 0.57%, 1.10% and 1.57%. Moreover, the film temperature near the inlet position significantly increases due to the enhanced countercurrent of the lubricant. As shown in Figure 19b, the maximum relative percentages of the temperature rise near the inlet are subsequently 9.9%, 7.1% and 4.4% for the above three environment temperatures.  By combining Figures 18 and 19, the reasonable reduction of the environment temperature is helpful to forming the lubrication film. Too low an environment temperature, however, can result in higher temperature rise near the inlet due to the inertia force effect of the lubricant.

Conclusions
A non-Newtonian thermal elastohydrodynamic lubrication (TEHL) model considering the inertia forces of the lubricant is established for the elliptic-contact problem. Meanwhile, a central film thickness experiment is performed to prove the rationality of the proposed TEHL model. Further, effects of the applied load, entrainment velocity, slide-roll ratio, and environment temperature on the TEHL performances are investigated numerically. The main conclusions drawn are as follows: -When the inertia force is considered, the film pressure peak in the contact region is little changed. -With the inertia force considered, the central and minimum film thicknesses increase because the inertia force overcomes the stickiness between the lubricant molecules. -When the inertia force is considered, the countercurrent zone of the lubricant expands towards the inlet, the maximum film temperature increases slightly, and the significant temperature rise near the inlet appears. -The numerical result of the central film thickness under the inertial situation is closer to the experimental result, which shows that it is reasonable and necessary to consider the inertia forces of the lubricant in the elliptic contact non-Newtonian TEHL analysis.