An Analytical Method for the Determination of Temperature Distribution in Short Journal Bearing Oil Film

: The aim of this paper is to derive an equation for the temperature distribution in journal bearing oil ﬁlm, in order to predict the thermal load of a bearing. This is very important for the prevention of critical regimes in a bearing operation. To achieve the goal, a partial di ﬀ erential equation of the temperature ﬁeld was ﬁrst derived, starting from the energy equation coupled with the Reynolds equation of hydrodynamic lubrication for a short bearing of symmetric geometry. Then, by solving the equation analytically, the function of temperature distribution in the bearing oil ﬁlm has been obtained. The solution is applied to the journal bearing, for which the experimental data are available in the references. Finally, the obtained results have been compared to the corresponding experimental values for two operating regimes, and a good level of agreement was achieved.


Introduction
Journal bearings are some of the oldest and most relevant machine elements. When it comes to the design, the concept of these bearings did not significantly change over the course of history, retaining its original simplicity and functionality. However, as the needs of mankind grew, machines were becoming more and more complex, and they were working in increasingly difficult conditions. To keep the bearing operation reliable even in such conditions, the materials of the bearings and their lubricants were improved over time. The bearing materials were changed in terms of improving their slip properties and embeddability, and the lubricants were changed in terms of improving their rheological properties. These improvements are the result of the studies of many researchers.
The most significant discoveries in the field of journal bearings were made by three researchers: N. P. Petrov, B. Tower and O. Reynolds [1], who recognized and formulated the phenomenon of hydrodynamic lubrication independently of each other. In the hydrodynamic lubrication regime of a journal bearing, the journal and its housing are separated by a thin layer of oil, and the friction and wear are minimized. Viscous friction of the oil causes heat generation, which in certain conditions can significantly worsen the performance of the bearing. A higher temperature of the oil film reduces the load carrying capacity of the bearing due to decrease of the oil film stiffness. In extreme cases, this can lead to a seizure of the bearing due to the rupture of the oil film, when "welding" takes place at the contact points of micro-irregularities. For this reason, the need for knowing the distribution adiabatically. A comparison is then made with the circular case to analyse the effect of this irregularity. Moreno et al. [16] transformed Reynolds equation, the energy equation, and the diffusion equation for the bearing into the corresponding finite difference equations. In order to solve these equations, they developed their own model based on an analogy with the electrical circuit. Pierre et al. [17] carried out a thermo-hydrodynamic analysis of a misaligned journal bearing to predict the size of the misalignment, pressure and temperature under steady-state conditions. The analysis is based on the numerical solution of Reynolds equation, the energy equation, and Laplace heat conduction equation for the bush material. Chauhan et al. [18] conducted a thermo-hydrodynamic analysis of an elliptical journal bearing for three different grade oils in order to determine which of them is most resistant to thermal degradation. The parameters for evaluating this degradation were oil pressure, oil temperature and load carrying capacity of the bearing. These parameters were determined by solving the Reynolds equation, the energy equation, and Laplace heat conduction equation for the bush material numerically. The reference [19] presents a model for pressure and temperature calculations in a journal bearing both for the steady and the time-dependent case, taking into account shaft and bushing thermal exchange with the external environment. While solving the problem, different numerical methods were used: the finite element method, the finite difference method, and the boundary elements method. Kornaev et al. developed a computational model of a plain fluid-film bearing with means to create artificial thermal and viscosity wedge effect. The model is based on the generalized Reynolds equation and the energy equation that are solved simultaneously using finite difference method combined with the iteration procedure [20]. Babin et al. investigated the possibility of controlled lubrication principles application and presented a complex mathematical model of an active thrust fluid-film bearing with a central feeding orifice. The model incorporates the Reynolds equation and the adiabatic energy equation, which are solved numerically using finite difference method [21]. Li et al. investigated the characteristics of water-lubricated conical hydrodynamic and hybrid bearings numerically solving the generalized Reynolds' equation and energy equation with enthalpy, and considered the effects of the gas-liquid two phases turbulence temperature, and axial flow [22].
In all the above-mentioned studies, the thermo-hydrodynamic lubrication equations are solved using different numerical methods, which can be time consuming and very demanding in terms of computer resources engagement. This can be a serious disadvantage if the real-time monitoring of a machine is necessary for diagnostic reasons. In such situations, it would be convenient to have an analytical solution for the temperature field in the oil film of the journal bearing, which could help faster diagnosis of the bearing condition. This motivated the authors of the paper to give their contribution in solving the problem. To that end, they derived an analytic expression for a two-dimensional temperature field in a journal bearing, solving simultaneously the adiabatic energy equation with a simplified Reynolds equation for a short bearing. As a criterion for assessing whether a bearing is long or short, the slenderness ratio (bearing length / bearing diameter) is used. If the ratio is greater than 2, the bearing is long, and if it is less than 1, the bearing is short [23,24]. Reynolds equation for a short bearing is a good approximation of the two-dimensional Reynolds equation if the slenderness ratio is less than 1. Since the slenderness ratio of modern bearings generally lies in the range between 0.5 and 1 [23,24], the temperature distribution equation for a short bearing would be applicable to most of the bearings currently used.
In order for the oil film temperature distribution equation to be derived analytically, it was necessary to simplify the problem by introducing a number of assumptions. For the sake of the reader, all the assumptions are listed below:

•
The oil is an incompressible Newtonian fluid.

•
There is no slip of the oil at the boundaries.

•
The oil film thickness is small compared with the other dimensions.

•
Pressure and temperature are constant through the thickness of the oil film. • Circumferential pressure gradient in the oil film is neglected.

•
The oil viscosity and specific heat are constant.
• There is no thermal interaction between the oil film and the surrounding bearing structure.

•
The bearing is operating under steady-state conditions (journal speed and bearing load are constant).

•
The structural components of the bearing are rigid and smooth.

•
There is no misalignment in the bearing structure.

•
The oil supply is not taken into account.

•
No asperity contacts between the journal and the bearing.

•
The bearing geometry is symmetric.

•
The oil flow is laminar.
The rest of the paper is structured as follows. In Section 2, the initial equations are given for solving the temperature field problem in a journal bearing. Section 3 describes the procedure of how to derive and solve the differential equation of the temperature field for a short bearing. Section 4 contains the presentation and discussion of the results obtained after the derived equation is applied on a particular journal bearing. The conclusions of the research are given in Section 5.

Governing Equations
In a bearing bush with the length L and the radius R, a journal of the radius r performs rotational movement with the constant angular velocity ω under hydrodynamic lubrication conditions, as shown in Figure 1.  There is no thermal interaction between the oil film and the surrounding bearing structure.  The bearing is operating under steady-state conditions (journal speed and bearing load are constant).  The structural components of the bearing are rigid and smooth.  There is no misalignment in the bearing structure.  The oil supply is not taken into account.  No asperity contacts between the journal and the bearing.  The bearing geometry is symmetric.  The oil flow is laminar.
The rest of the paper is structured as follows. In Section 2, the initial equations are given for solving the temperature field problem in a journal bearing. Section 3 describes the procedure of how to derive and solve the differential equation of the temperature field for a short bearing. Section 4 contains the presentation and discussion of the results obtained after the derived equation is applied on a particular journal bearing. The conclusions of the research are given in Section 5.

Governing Equations
In a bearing bush with the length L and the radius R , a journal of the radius r performs rotational movement with the constant angular velocity ω under hydrodynamic lubrication conditions, as shown in Figure 1. Hydrodynamic lubrication in steady-state conditions for an incompressible fluid, is described by Reynolds equation of form [5,8]. where x -streamwise coordinate direction, θ -angular coordinate, z -axial coordinate direction, p -pressure in an arbitrary point ( x , z ), c -radial clearance, Hydrodynamic lubrication in steady-state conditions for an incompressible fluid, is described by Reynolds equation of form [5,8]. where x-streamwise coordinate direction, θ-angular coordinate, z-axial coordinate direction, p-pressure in an arbitrary point (x,z), η-oil viscosity, Energy balance of oil film in the bearing, assuming there is no heat exchange with the mating surfaces, is described by the adiabatic energy equation of the form [5] where ρ-density of oil c p -specific heat at constant pressure T-oil film temperature u-circumferential velocity w-axial velocity y-cross-film coordinate direction.
In order to determine the temperature field in the oil film, it is convenient to transform the Equation (2) to a form containing pressure gradients, i.e., the derivation procedure of which is provided in Appendix A.
In order to solve Equations (1) and (3), they need to be converted into a dimensionless form. For that purpose, the appropriate dimensionless variables θ, z * , λ, h * , p * and T * are introduced according to the expressions: Bearing in mind the dimensionless variables from the Equation (4), the Equations (1) and (3) become respectively. The dimensionless parameter λ, defined in the Equation (4), is called the slenderness ratio and is the criterion for estimating the length of the journal bearing. If it is greater than 2, then the bearing is long, and if it is less than 1, the bearing is short [23,24]. As most journal bearings in engineering applications belong to short bearings, the focus will be directed to such bearings only. In the case of a short bearing, the circumferential pressure gradient can be neglected in comparison with the axial pressure gradient. This assumption enables the coupled Equations (5) and (6) to be solved for a short bearing, in order to determine the temperature field in the bearing oil film. For a better understanding of the entire work, a block diagram of the solution algorithm is given in Figure 2.

Solution Procedure
Since the circumferential pressure gradient is neglected in short bearings, the Reynolds and the energy equations have the following form By solving the equation (7), the expression for the axial pressure gradient is obtained Introducing the equation (9) into the equation (8) it becomes or

Solution Procedure
Since the circumferential pressure gradient is neglected in short bearings, the Reynolds and the energy equations have the following form By solving the Equation (7), the expression for the axial pressure gradient is obtained Introducing the Equation (9) into the Equation (8) it becomes The solution of the Equation (11) is assumed in the form where f (θ) and g(θ) are the functions to be determined. Finding partial derivatives of the Equation (12) for θ and z * and substituting them into the Equation (11) it gives By grouping the terms containing (z * ) 2 and the terms that are functions only of θ into two separate groups, the following equations are obtained After dividing the Equation (16) with (z * ) 2 , we get The Equation (18) is a linear, non-homogeneous equation, solved by the constant variation method. For this purpose, only the homogeneous part of the differential equation is observed and it has the form The newly obtained homogeneous differential equation is separable, i.e., After introducing the substitution the Equation (20) becomes Integrating the Equation (22) gives i.e., where C is an integration constant for the homogeneous differential equation. As the non-homogeneous differential equation is solved here, then C is actually a function of θ. In order to determine C(θ), the Equation (24) should be introduced into the Equation (18), which gives a new differential equation i.e., dC dθ or The solution of the Equation (27) may be written in the general form C(θ) = 6(R/c) 2 λ 2 ε 2 · I 2 0 6 , where The integral from the Equation (29) is of the type and its solving is shown in Appendix B. The solution is After substituting the Equation (31) in the Equation (28), it can be written where D 2 is a new integration constant. The Equation (24) now becomes f (θ) = 6(R/c) 2 λ 2 ε 2 · I 2 0 6 + D 2 · (1 + ε cos θ) 2 .
After determining the function f (θ), the function g(θ) is determined by solving the Equation (17). With respect to the Equation (30), the function g(θ) may be written in the form The solving of integral I 0 0 2 is described in Appendix B and the solution is Now it is where D 3 is a new integration constant. Taking into account the Equations (33) and (36), the Equation (12) becomes The constants D 2 and D 3 are determined from the boundary condition T * (0, z * ) = T 0 * , where T 0 * is the dimensionless temperature of the oil at θ = 0. The boundary condition is based on the assumption that the oil temperature does not change in the area between the oil inlet (vertical plane) and the attitude plane defined by the angle Φ ( Figure 1). Hence, Equation (37) becomes Also, it is assumed that the temperature T 0 * is constant along the axis of the bearing. In other words, it is independent of the coordinate z * , so it must be D 2 = 0. Therefore, D 3 = T 0 * and the Equation (38) becomes By using the expression for the dimensionless temperature from the Equation (4), the Equation (39) can be transformed into the form i.e., ( 1−ε 2 ) 1/2 arccos ε+cos θ 1+ε cos θ − ε sin θ ( 1+ε cos θ ) + T 0 .
(41) The Equation (41) is the final analytic form of the temperature distribution in the oil film of short journal bearing.
While deriving this equation, it is assumed that the viscosity of oil is not changed with the temperature, which is also one of the assumptions that the Reynolds equation derivation is based on. This is justified to some extent, since journal bearings are dominantly used in internal combustion engines, where they are lubricated by so-called multi grade oils, characterized by their viscosity resistance to temperature change. On the other hand, by introducing this assumption, a considerably simpler theoretical model of temperature distribution in the oil film is derived, on the basis of which satisfactory results can be obtained quite quickly.

Results and Discussion
In order to validate the mathematical model of temperature distribution in a journal bearing oil film (Equation (41)), the experimental results given in the reference [25] were used. In that sense, the bearing parameters were also taken from the reference [25] and are presented in Table 1. The operating point of the bearing is defined by the force F and journal speed n and the two parameters significantly affect the eccentricity ratio ε. In order to apply the temperature distribution Equation (41) to the data in Table 1, it is necessary to relate the parameters F, n and ε to each other. The parameters are related by the equation of the load capacity that is given in the form [8] By introducing the data from Table 1 into the Equation (42) and applying the "trial and error" method to the Equation (42) it is found that the value ε = 0.44 corresponds to the operating point (F = 4000 N, n = 2000 rpm), and the value ε = 0.38 corresponds to the operating point (F = 6000 N, n = 4000 rpm). The corresponding values of attitude angle Φ are 58.04 • and 62.39 • , respectively. Furthermore, by substituting the data from Table 1 and the corresponding ε values into the Equation (41), the theoretical oil film temperature distribution in the midplane cross-section of the bearing has been obtained at the two operating points.
The first operating point is defined by the journal speed of n = 2000 rpm and the load of F = 4000 N. The theoretical temperature distribution corresponding to that operating point is represented by the blue line in Figure 3. In Figure 3, the oil temperature values measured under the same operation conditions are also plotted and marked by the red circles. Similarly, Figure 4 shows the theoretical temperature distribution and the corresponding experimental values [25], at n = 4000 rpm and F = 6000 N.
The first operating point is defined by the journal speed of n = 2000 rpm and the load of F = 4000 N. The theoretical temperature distribution corresponding to that operating point is represented by the blue line in Figure 3. In Figure 3, the oil temperature values measured under the same operation conditions are also plotted and marked by the red circles. Similarly, Figure 4 shows the theoretical temperature distribution and the corresponding experimental values [25], at n = 4000 rpm and F = 6000 N.  Looking at the both diagrams, one can notice that in the area where the bearing angle is less than Φ , the measured oil temperature is slightly higher than at the oil inlet. That is the result of mixing the oil that enters the bearing with the oil that is already inside the bearing. This phenomenon is not taken into account while deriving the temperature distribution equation, hence the deviation of about 1-2 °C in relation to the previous assumption about the temperature constancy in the bearing area between the oil supply and the attitude plane. However, this The first operating point is defined by the journal speed of n = 2000 rpm and the load of F = 4000 N. The theoretical temperature distribution corresponding to that operating point is represented by the blue line in Figure 3. In Figure 3, the oil temperature values measured under the same operation conditions are also plotted and marked by the red circles. Similarly, Figure 4 shows the theoretical temperature distribution and the corresponding experimental values [25], at n = 4000 rpm and F = 6000 N.  Looking at the both diagrams, one can notice that in the area where the bearing angle is less than Φ , the measured oil temperature is slightly higher than at the oil inlet. That is the result of mixing the oil that enters the bearing with the oil that is already inside the bearing. This phenomenon is not taken into account while deriving the temperature distribution equation, hence the deviation of about 1-2 °C in relation to the previous assumption about the temperature constancy in the bearing area between the oil supply and the attitude plane. However, this Looking at the both diagrams, one can notice that in the area where the bearing angle is less than Φ, the measured oil temperature is slightly higher than at the oil inlet. That is the result of mixing the oil that enters the bearing with the oil that is already inside the bearing. This phenomenon is not taken into account while deriving the temperature distribution equation, hence the deviation of about 1-2 • C in relation to the previous assumption about the temperature constancy in the bearing area between the oil supply and the attitude plane. However, this assumption is justified by the fact that the experimental data show a slight change in temperature in the mentioned area.
In the area of temperature rise (the area between the angle Φ and Φ + 180), there is a very good agreement between the experimental and theoretical values with a deviation of up to 2 • C, at both operating points considered. Furthermore, the maximum oil temperature calculated is slightly higher than the maximum temperature measured, which can be useful if one wants to predict critical operating points of the bearing.
In the area between the bearing angle corresponding to the maximum temperature and the bearing angle corresponding to the oil inlet, temperature decreases. A fairly good agreement between theoretical and experimental results is clearly seen when the bearing operates under light-duty conditions (Figure 3). However, for heavier operating conditions (Figure 4), the deviation is slightly higher (up to about 3 • C).
A more complete view of temperature distribution in the oil film of the bearing is enabled by a 3D-diagram showing the temperature change along both the bearing circumference and length. Such a diagram for the bearing considered at n = 4000 rpm and F = 6000 N, is shown in Figure 5. It is noted that the temperature variation along the axis of the bearing has a slightly parabolic shape that is symmetrical in relation to the bearing midplane (z * = 0). Furthermore, the temperature decreases from the side planes to the midplane of the bearing. A similar temperature distribution was obtained experimentally by Mitsui et al. [26].
Symmetry 2020, 12, x FOR PEER REVIEW 12 of 19 assumption is justified by the fact that the experimental data show a slight change in temperature in the mentioned area.
In the area of temperature rise (the area between the angle Φ and Φ + 180), there is a very good agreement between the experimental and theoretical values with a deviation of up to 2 °C, at both operating points considered. Furthermore, the maximum oil temperature calculated is slightly higher than the maximum temperature measured, which can be useful if one wants to predict critical operating points of the bearing.
In the area between the bearing angle corresponding to the maximum temperature and the bearing angle corresponding to the oil inlet, temperature decreases. A fairly good agreement between theoretical and experimental results is clearly seen when the bearing operates under light-duty conditions (Figure 3). However, for heavier operating conditions (Figure 4), the deviation is slightly higher (up to about 3 °C).
A more complete view of temperature distribution in the oil film of the bearing is enabled by a 3D-diagram showing the temperature change along both the bearing circumference and length. Such a diagram for the bearing considered at n = 4000 rpm and F = 6000 N, is shown in Figure 5. It is noted that the temperature variation along the axis of the bearing has a slightly parabolic shape that is symmetrical in relation to the bearing midplane ( 0 * = z ). Furthermore, the temperature decreases from the side planes to the midplane of the bearing. A similar temperature distribution was obtained experimentally by Mitsui et al. [26]. The good agreement of the derived equation with the experimental data is encouraging when it comes to the possibility of its application. The authors hope that such an equation could serve to estimate the temperature of the oil film quickly, which may be useful in a diagnostic tool designed to monitor the bearing conditions in real time. A detailed analysis of the possibility of applying the derived equation in a diagnostic tool, as well as its potential implementation, could be a topic for future research. In the text to follow, only a hint is given of how this equation could work in real time.
Most of the input data in equation (41) are constants, while the others depend on operating parameters of the bearing. The constants are: -characteristics of the oil used (density ρ , viscosity η and specific heat at constant pressure p c ) and -geometric parameters of the bearing (bearing radius R , radial clearance c and bearing length L ).
These constants do not affect the real time monitoring concept and their values are stored in the computer's memory. The good agreement of the derived equation with the experimental data is encouraging when it comes to the possibility of its application. The authors hope that such an equation could serve to estimate the temperature of the oil film quickly, which may be useful in a diagnostic tool designed to monitor the bearing conditions in real time. A detailed analysis of the possibility of applying the derived equation in a diagnostic tool, as well as its potential implementation, could be a topic for future research. In the text to follow, only a hint is given of how this equation could work in real time.
Most of the input data in Equation (41) are constants, while the others depend on operating parameters of the bearing. The constants are: -characteristics of the oil used (density ρ, viscosity η and specific heat at constant pressure c p ) and geometric parameters of the bearing (bearing radius R, radial clearance c and bearing length L).
These constants do not affect the real time monitoring concept and their values are stored in the computer's memory.
On the other hand, the operating parameters (oil inlet temperature T 0 , journal speed n and bearing load F) are measured using appropriate sensors. On the basis of the measured values, the angular velocity ω and eccentricity ratio ε are calculated. Such a scenario could be possible for steady-state conditions of the bearing operation. After each transition from one steady-state condition to another, the temperature distribution should be recalculated.

Conclusions
In this paper, the equation of temperature field in the oil film of short journal bearing is derived analytically. The results of applying this equation to a journal bearing have been compared to the corresponding experimental results of other researchers for the same bearing, showing a small deviation. The good agreement of the derived equation with experimental data, as well as its explicit analytical form, allows for temperature distribution in the oil film to be determined quickly, with satisfactory accuracy but without demanding computer resources. This indicates that the derived equation might be useful in a computer diagnostic tool as a part of a more complex procedure for the real time monitoring of the bearing condition.
Another contribution of this research is that the energy equation is transformed into a form, which allows it to be coupled with the Reynolds equation for short journal bearing. As a result, a partial differential equation of the temperature field for a short bearing is obtained, which is in the form that allows it to be solved analytically and, as far as the authors are aware, it is not encountered in any previous literature.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Conversion of the Energy Equation into the Form Containing Pressure Gradients
The initial form of the energy equation is [5] In order for this equation to be coupled with the Reynolds equation, it is necessary to express flow velocities u and w through the pressure gradients ∂p ∂x and ∂p ∂z . All these quantities are contained in the Navier-Stokes equations (from which Reynolds equation is derived), and it is these equations that are used to establish the relationship between the energy equation and the Reynolds equation. By applying the Navier-Stokes equations to the oil flow in the journal bearing, their simplified form is obtained The Equation (A3) shows that the pressure does not change along the oil thickness, i.e., p = const. The solutions of the quadratic equation are and y 2 = y us2 = 1 12ah 6ah 2 − 6U + 2 (3a 2 h 4 + 9U 2 ) .
After introducing these solutions into the Equation (A16) and arranging it, we get and du dy (y = y us2 ) = √ 3 3h (3U 2 + a 2 h 4 ). The solutions of the quadratic equation are and After introducing the solutions into the Equation (A17) and arranging it, we get dw dy (y = y ws1 ) = −