Abstract

This research focuses on the multiphase oil film tribology between the piston pin and the connecting rod in an internal combustion engine and establishes a new computational approach for thin-film lubrication with unsteady flow channel variation. First, the pin and the connecting rod are considered as rigid bodies, and 3D numerical analysis of the cavitating lubricating oil flow is performed when combustion load is applied to the pin. We find that dynamic pressure does not increase around the connecting rod edge and that pressure is potentially insufficient to support the load. In the second numerical analysis, the pin and the connecting rod are considered to be elastically deformable structures, and coupled 3D multiphase fluid–structure interaction simulation is performed. The boundary lubrication area is detected using a statistical Greenwood–Tripp model as unevenness of the contacted metal surface. The results show that pressure distribution spreads more widely than in the result for rigid bodies and that the film was thicker as well. Also, the pin deformed like a bow, but the deformation of the connecting rod was quite small, suggesting a potential mechanical contact at the edge of the connecting rod with the pin. By comparison with an actual operationally used piston pin, we find that the fluid–structure coupled analysis qualitatively predicted the seizure location.

1 Introduction

There is a growing need for further improvements to the efficiency of internal combustion engines due to their impact on environmental and sustainability problems [1]. In internal combustion engines for automobiles, mechanical energy loss accounts for 10–20% of the energy produced, so there are increasing efforts toward its reduction [2]. Furthermore, about 40% of mechanical loss is due to friction in the piston system, so it is important to clarify tribological phenomena there [3]. Reducing the weight and size of parts and reducing the viscosity of lubricating oil effectively realizes low friction, but at the same time, there is an increased risk of internal combustion engine failure due to contact between parts.

In a reciprocating engine, power generated by combustion is extracted as reciprocating motion and converted to rotary motion for driving force [4]. Reciprocating engines have many uses but are primarily used to power automobiles. It is important to efficiently transmit a large amount of power from combustion to the drive system, but high durability is required for parts, which must operate under extremely harsh conditions. The most common cause of failure in a reciprocating engine is seizure of its parts, a phenomenon that occurs when there is a break in the lubricating, thus allowing metal parts to come into contact and causing damage or adhesion (Fig. 1). Seizure generally makes engines inoperative, often causing such extensive damage that the entire engine must be replaced. Lubricating oil plays a major role in preventing seizure, so reciprocating engines must be designed to maintain an oil film between parts.

Fig. 1
A piston pin mechanically contacting the small end of the connecting rod seizes. Adapted from DEEP FIELD website.
Fig. 1
A piston pin mechanically contacting the small end of the connecting rod seizes. Adapted from DEEP FIELD website.
Close modal

In a reciprocating engine, reciprocating motion of a piston is converted into rotational motion of a crankshaft through use of a crank mechanism [5]. The rod connecting the piston and the crankshaft is called the connecting rod, and its structure is divided into a small end connected to the piston side and a large end connected to the crankshaft side. Both ends are connected by pins, and relative rotational movement occurs at a sliding part. In particular, the small end is connected by a pin called the piston pin, which is one of the most important lubrication points in a reciprocating engine. The small end supports the explosive pressure generated in the combustion chamber, transmitting that force to the connecting rod and crankshaft. In contrast to the large end, the relative motion of the small end is not rotational motion that keeps sliding in the same direction but oscillating motion that changes the direction of rotation at regular intervals. In other words, at the moment when the direction of rotation switches, the relative velocity becomes zero, and there is a very high possibility that metal parts will come into contact, causing seizure. From the above, lubrication between the piston pin and connecting rod is bearing lubrication under the harshest conditions in an internal combustion engine.

Since the piston system repeats up–down reciprocating motion at high speed, large inertial forces are applied. Weight reductions in the piston system are thus indispensable for improving fuel efficiency and durability. The piston pin is generally made of iron so that it can withstand deformation due to explosive pressure, so its mass is disproportionately large compared with that of the entire piston system. Currently, the piston pins used in reciprocating engines are mostly what is called the fully floating type, in which there is clearance between the piston pin and the connecting rod. Filling this clearance with lubricating oil prevents the seizure of these parts. The piston pin has a hollow cylindrical shape for weight reduction, but its inner diameter is optimized to minimize deformation even under heavy loads. Automakers and others are making various efforts toward reducing the weight of piston pins. Measures to reduce their diameter, for example, include distributing the applied load.

As described above, there is increasing demand for internal combustion engines with higher efficiency, and research on piston system lubrication for reducing mechanical energy loss is actively conducted. It is difficult to visualize and quantitatively analyze lubrication phenomena in experiments, so it is essential to establish numerical analysis methods that can be applied on computers. There have been some interesting studies involving both simulations and experimental methods related to piston pin lubrication problems. Allmaier and Sander [6] found that most piston pins rotate at very slow rotational speeds and even change rotational directions according to operating conditions. They also investigated influencing effects on this dynamic behavior such as clearance and pin surface roughness to see how the lubrication of this crucial part could be improved. Fang et al. [7] developed a numerical model to analyze lubricated full-floating pin bearings in planar multibody systems, which couple a mixed lubrication model with multibody motion equations. To consider the effects of surface roughness and asperity contact, they integrated the average Reynolds equation and the GT asperity contact model to establish the mixed lubrication model. Ferretti et al. [8] performed a tribological study showing that two different asperity contact models have been implemented, the former based on standard Greenwood–Tripp theory and the latter based on a complementarity formulation of the asperity contact problem, and compared the results. They showed that both approaches can provide similar results if the model input data are calibrated.

There have also been many papers on piston pin and connecting rod tribology, focusing mainly on mechanisms and multibody dynamics [7], but few studies have focused on multiphase flow behavior inside narrow thin-film channels, which is the subject of this research.

Even in studies of the state between connecting rods and the small end of piston pins, analyses commonly use the Reynolds equation [9,10], which is a simplification of the Navier–Stokes equation, to perform numerically stable analyses. Furthermore, the large loads applied to the piston pin deform the piston pin–connecting rod structure. Conventional lubrication theory underestimates the gap width because it assumes that the structure is a rigid body, so it is necessary to perform analyses based on elastohydrodynamic lubrication (EHL) theory, which treats the structure as an elastic body [1116]. In general EHL theory, however, only one structure is treated as an elastic body [17,18], so there are limits to considering deformation of both the piston pin and connecting rod, which is the subject of this research.

Ozasa et al. [19] used EHL theory to perform calculations for the large end of the connecting rod, which experiences less harsh lubrication conditions compared with those of small end, confirming the validity of calculation results regarding film thickness and axial locus by comparison with experimental results from an actual machine. Lin et al. [20] performed numerical analyses of lubrication between the small end of the connecting rod and the piston pin. They used the Reynolds equation as the governing equation for lubricating oil and attempted to elucidate seizure phenomena of the piston pin by considering structural deformation. These numerical analyses demonstrated that deformation of the piston pin could cause seizure with the small end of the connecting rod, but satisfactory results were not obtained for lubricating oil flow, so this work serves only for discussion.

Thermal effects in the lubricant are required to formulate the governing equations because temperature distribution due to frictional heat is formed in the thin-film layer where unsteady fluctuating pressure is applied. Also, cavitation occurs with phase change, which influences the behavior and seizure location. In addition, although many previous studies have considered the effects of cavitation in the tribology field [2127], few studies have evaluated interactions between EHL effects and cavitation.

This research focuses on the multiphase narrow gap flow between the small end of piston pins and connecting rods with the unsteady variation of channel width. Furthermore, it established a numerical method based on the fluid–structure interaction (FSI) model that considers the momentum and energy exchange at the piston pin-lubricant interface and the exchange at the conrod-lubricant interface. To that end, we develop a 3D multiphase fluid–structure coupled analysis that considers elastic deformation of the structure and changes in the flow channel, along with an analysis method that reproduces the complex phenomena of lubricating oil with cavitation under harsh conditions. Our objective is to make proposals for optimizing part shapes and predicting seizure points using the newly developed analysis method. This study also focuses on elastic deformation of piston connecting rods due to flow in a narrow gap with phase change and fluid–structure interactions rather than multibody and kinematic investigations such as the motion of piston pin rotation. We also performed a coupled multiphase fluid–structure analysis with 3D unsteady narrow-channel deformation, which is difficult to analyze with TEHL theory [28], so a simple comparison with TEHL theory is difficult. We perform a comparison of rigid and elastic models to demonstrate the effectiveness of this analysis.

2 Rigid-Body Analysis of Lubricating Oil Flow Between the Piston Pin and the Small End of the Connecting Rod

We first perform a coupled analysis using rigid-body approximation, then a comparative study to clarify the effect of elastic deformation on the EHL characteristics. The numerical analyses in this section are used to perform composite fluid analysis of the lubricating oil and rigid bodies, thereby reproducing flow between the piston pin and connecting rod. These numerical analyses consider the connecting rod as a rigid body and apply load to reproduce a flow channel that unsteadily changes. This section describes the fluid model used in the analyses, the mesh morphing method, the analysis conditions, and the analytical results. To consider the numerical analysis results, we also describe an analytical method for pressure distributions in the bearing, considering that lubrication between the piston pin and connecting rod can be simulated by the bearing lubrication model.

2.1 Basic Fluid Phase Equations and Numerical Method.

This analysis considers lubricating oil as a compressible fluid by using the continuity equation and the Navier–Stokes equation, where are respectively given by
ρt+(ρU)=0
(1)
and
ρUt+(ρUU)=p+τ
(2)
where ρ is the density, U is the flow velocity vector, p is the pressure, and τ is the viscous stress tensor represented as
τ=μ[U+(U)T23(U)I]
(3)
where μ is the viscosity and I is the unit tensor.

2.2 Cavitation Model.

When a large load is applied to and displaces the piston pin, a high-pressure region is generated in the direction of load application, but a large negative pressure is generated on the opposite side. If this negative pressure is below the saturated vapor pressure, the lubricating oil evaporates in parts of that region, causing cavitation. Identification of the gas–liquid phase interface is generally unimportant for cavitation in bearing lubrication, so a homogenous equilibrium model (HEM) [2934] with high numerical stability is used in such numerical analyses. The region of cavitation generation in an HEM has a single density value, with density of the lubricating oil expressed as
ρ=αρv+(1α)ρl
(4)
where v and l subscripts indicate vapor and liquid phases, respectively, and the volume fraction of vapor α is
α=ρρl,satρv,satρl,sat
(5)
where the subscript “sat” indicates a value at saturated vapor pressure.

2.3 Lubrication Compressibility.

Extremely large loads are applied to the lubricating oil between the piston pin and connecting rod, so we must consider the compressibility of both the liquid and vapor phases. The liquid phase density is calculated by the Dowson–Higginson equation [14,3539], which is widely used in EHL theory
ρl=ρl,0C1+C2pC1+p
(6)
Here, ρl,0 is the liquid phase density under atmospheric pressure, and coefficients C1 = 5.9 × 108 (Pa) and C2=1.34().
Vapor phase density is given by
ρv=ψvp
(7)
where ψv is the compressibility of vapor phase [33,40].

2.4 Rheology Model.

Viscosity of the lubricating oil, the target of this numerical analysis, is known to change under conditions where the pressure becomes extremely high. This is called the piezo-viscosity effect, and many studies have aimed to predict the correct viscosity under high shear-rate conditions. As the load increases and pressure in the lubricating oil rises, the density changes due to compressibility, decreasing the kinematic viscosity coefficient. However, as pressures increase further, the viscosity coefficient of the lubricating oil sharply increases due to the piezo-viscosity effect, significantly affecting the change in kinematic viscosity. To predict the viscosity coefficient of lubricating oil under extremely high shear-rate conditions, it is therefore necessary to use an appropriate rheological model [14,15,35,41,42].

To model the pressure dependence of the lubricating oil in this numerical analysis, we use the Barus equation
ηBarus=η0exp(αpvp)
(8)
where η0 is the viscosity coefficient under atmospheric pressure, and αpv is the pressure–viscosity coefficient, which is a characteristic value specific to the fluid. Lubricating oils are known to show non-Newtonian properties under high shear-rate conditions. In this numerical analysis, we use the Eyring model, usefulness of which has been confirmed in EHL research
γ˙=τ0ηBarussinhττ0
(9)
Here, ηBarus is the fluid viscosity as calculated from Eq. (8), γ˙ is the shear rate, τ is the shear stress, and τ0 is the Eyring stress. Shear stress τ is calculated as
τ=ηBarusγ˙
(10)
By substituting Eq. (10) into Eq. (9), we obtain the Eyring equation
ηEyring=τ0γ˙sinh1(ηBarusγ˙τ0)
(11)
In this numerical analysis, when the shear rate γ˙ exceeds 1.0 × 10−8 1/s we use Eq. (11) to find the viscosity of the lubricating oil, and Eq. (8) otherwise.

2.5 Mesh Morphing Model.

This numerical analysis predicts the flow of lubricating oil in the clearance between the piston pin and connecting rod, which changes when a load is applied to the piston pin. We therefore perform mesh morphing to deform the original mesh describing the flow channel. Mesh morphing follows the Laplace equation
(γd)=0
(12)
where γ is the displacement diffusion coefficient and d is the amount of mesh displacement [43]. Displacement diffusivity is a parameter that determines the mesh distribution in the calculation area and is a function of the distance r from the rigid-body interface, given by
γ=γ(r)
(13)
Mesh morphing has a significant effect on mesh quality. In particular, when the mesh displacement is large or when the mesh is displaced in a narrow region, it is necessary to set the displacement diffusion coefficient so that orthogonality is maintained.

2.6 Analytical Method for Bearing Lubrication.

In this computation, numerical constraints necessitate adding boundary conditions for oil inflow. An approximate analytical solution can be obtained for such bearing lubrication. Therefore, we obtain an analytical solution by applying the Sommerfeld condition [44] to the Reynolds equation [10,45], which is a simplification of the Navier–Stokes equation.

2.7 Half-Sommerfeld Condition.

In the Sommerfeld condition solution shown above, the pressure is negative below atmospheric pressure in the range where the gap width becomes large. In that negative pressure region, the lubricating oil falls below the saturated vapor pressure and evaporates, causing cavitation, so the assumption that the bearing gap is filled with only lubricating oil does not hold under these conditions.

The various available cavitation models are extensively compared in Ref. [46], and both the Reynolds (Swift-Stieber) condition and the Jakobsson–Floberg–Olsson (JFO) cavitation model [4749] have been used in conventional research to investigate the performance of laser-textured piston rings or sliders [5052].

However, it still needs to be experimentally validated if multiphase cavity flow with narrow channel deformation can be accurately modeled with a thin-film approximation (e.g., the Reynolds equation) employing a saturation parameter (e.g., the JFO cavitation model) [5357].

However, mass conservation is not satisfied, so the load capacity is overestimated in cases where mass flux limits the pressure increase. Still, behavior of the asperity parameters was qualitatively identical to that in the JFO model, resulting in the same optimal asperity parameters [58].

The half-Sommerfeld condition proposed by Gümbel [5961] is a model that considers the occurrence of such cavitation. In this model, the negative pressure part in the pressure distribution obtained under the Sommerfeld condition is set to 0 pressure. For these reasons, we applied the half-Sommerfeld condition in this study.

Figure 2 shows a schematic of the bearing lubrication area considering cavitation, and Fig. 3 shows the circumferential pressure distribution obtained under the half-Sommerfeld condition.

Fig. 2
Bearing lubrication area considering cavitation
Fig. 2
Bearing lubrication area considering cavitation
Close modal
Fig. 3
Circumferential pressure distribution obtained under the half-Sommerfeld condition: (a) Reynolds equation model and (b) current FSI model
Fig. 3
Circumferential pressure distribution obtained under the half-Sommerfeld condition: (a) Reynolds equation model and (b) current FSI model
Close modal

2.8 Analysis Model and Numerical Analysis Conditions.

This section describes the analysis model and numerical analysis conditions used for flow analysis of lubricating oil in the narrow flow channel between the piston pin and connecting rod. We performed this numerical analysis and computing using cavitatingTribologyFoam, a customized cavitation solver based on our previously developed solver cavitatingLesInterFoam [33] in the open-source computational fluid dynamics (CFD) code openfoam [62].

Lubricating oil flows in from an oil reservoir at the top of the connecting rod, filling the gap between the piston pin and the connecting rod. In this analysis, therefore, we extract the gap region between the piston pin and the connecting rod and create a fluid mesh. Figure 4 shows the calculation region and the mesh used for analysis, with the z-axis showing the axis of rotation, the x-axis showing the horizontal radial direction, and the y-axis showing the vertical radial direction. Figure 5 shows dimensions in the numerical model.

Fig. 4
Fluid-phase calculation region and mesh
Fig. 4
Fluid-phase calculation region and mesh
Close modal
Fig. 5
Dimensions in the numerical model
Fig. 5
Dimensions in the numerical model
Close modal

Both the piston pin and the connecting rod are treated as rigid bodies. The initial clearance width between the two parts is extremely narrow at 11.75 μm, but we created a mesh with seven layers to obtain a sufficiently accurate velocity distribution between the piston pin and connecting rod (Fig. 6). The number of mesh cells is approximately 200,000. We apply a load in the y-axis direction and treat the piston pin as a moving rigid body. In an actual reciprocating engine, such as a four-stroke engine, each cycle comprises four steps: intake, compression, combustion, and exhaust. The total inertial force due to vertical reciprocating motion and explosive pressure received from the combustion chamber is applied to the piston pin. Figure 7 shows the load history on a piston pin over the four cycles of a reciprocating engine. Based on actual engine operating conditions, the load history on the piston pin data and the connecting rod velocity data were obtained by flexible multibody dynamic analysis using RecurDyn [63,64]. The highest load of 4672 N occurs during the combustion process. To reduce calculation times in this analysis, instead of unsteady analysis through four cycles, we perform calculations by extraction from the middle of the compression process to the middle of the exhaust process, where the possibility of seizure is the highest. Figure 8 shows the load history applied to the piston pin by numerical analysis.

Fig. 6
Mesh in the thin-film region between the piston pin and connecting rod
Fig. 6
Mesh in the thin-film region between the piston pin and connecting rod
Close modal
Fig. 7
Load on the piston pin over four cycles
Fig. 7
Load on the piston pin over four cycles
Close modal
Fig. 8
Load on the piston pin during the combustion process
Fig. 8
Load on the piston pin during the combustion process
Close modal

Table 1 and Fig. 9 show the boundary conditions used for numerical analysis. Both the upper oil reservoir and the gap between the piston pin and connecting rod are set to free inflow and outflow conditions, and the pressure boundary is set to atmospheric pressure. In other words, lubricating oil flows in or out depending on the pressure difference between the inside and outside of the calculation region. Both the piston pin and the connecting rod wall surface have no-slip boundary conditions. Oscillatory movement of the connecting rod provides velocity according to the swing angle shown in Fig. 10 as the velocity boundary condition on the piston pin wall surface.

Fig. 9
Schematic of boundary conditions for the analysis model
Fig. 9
Schematic of boundary conditions for the analysis model
Close modal
Fig. 10
Angle history of the connecting rod due during oscillatory movement
Fig. 10
Angle history of the connecting rod due during oscillatory movement
Close modal
Table 1

Boundary conditions

BoundaryPressureVelocity
InletConstant atmospheric (101.325 kPa)Continuous in–out flow
OutletConstant atmospheric (101.325 kPa)Continuous in–out flow
Piston pinNeumann conditionNo slip
Connecting rodNeumann conditionNo slip
BoundaryPressureVelocity
InletConstant atmospheric (101.325 kPa)Continuous in–out flow
OutletConstant atmospheric (101.325 kPa)Continuous in–out flow
Piston pinNeumann conditionNo slip
Connecting rodNeumann conditionNo slip

Since the gap between the piston pin and connecting rod changes unsteadily, we apply an adaptive step and perform the unsteady calculation so that the Courant number is one or less. We also assume that the engine is operating at 100C, with kinematic viscosity ν = 1.0 × 10−5 m2/s, and density ρ = 823 kg/m3 for the lubricating oil. In addition, since the Reynolds number Re is assumed to be about 100 at maximum, the flow is assumed to be laminar. Since this study deals with multiphase flow with cavitation and a free surface based on phase change, a high-precision CFD model is necessary to accurately capture the interfacial flow in a narrow-channel thin-liquid film. Table 2 shows physical characteristics of the lubricating oil used in the analysis. We assume SAE-viscosity grade 10W–30 [65] as physical properties of the oil.

Table 2

Physical characteristics of the lubricating oil [65]

LubricantVaporUnit
Viscosity1.0 × 10−58.97 × 10−6(m2/s)
Density8500.029(kg/m3)
Compressibility4.9 × 10−75.7 × 10−6(1/kPa)
Saturated vapor pressure5000(Pa)
Pressure–viscosity index0.689(–)
Eyring stress5.0 × 106(Pa)
LubricantVaporUnit
Viscosity1.0 × 10−58.97 × 10−6(m2/s)
Density8500.029(kg/m3)
Compressibility4.9 × 10−75.7 × 10−6(1/kPa)
Saturated vapor pressure5000(Pa)
Pressure–viscosity index0.689(–)
Eyring stress5.0 × 106(Pa)

2.9 Numerical Analysis Calculations and Discussion.

The following describes the results of flow analysis of the lubricating oil between the piston pin and connecting rod. Figures 1116 respectively show pressure distribution contours (left) and graphs in which the perimeter of the z = 0 vertical cross section is plotted on the horizontal axis and the pressure is plotted on the vertical axis (right). Additionally, Figs. 1416 respectively show the lubricating oil velocity vector distribution (left) and the volume fraction distribution of vapor phase due to cavitation (right) at each time.

Fig. 11
Pressure distribution at t = 5.0 × 10−4 s (left) and z = 0 cross-sectional pressure distribution (right)
Fig. 11
Pressure distribution at t = 5.0 × 10−4 s (left) and z = 0 cross-sectional pressure distribution (right)
Close modal
Fig. 12
Pressure distribution at t = 2.0 × 10−3 s (left) and z = 0 cross-sectional pressure distribution (right)
Fig. 12
Pressure distribution at t = 2.0 × 10−3 s (left) and z = 0 cross-sectional pressure distribution (right)
Close modal
Fig. 13
Pressure distribution at t = 2.5 × 10−3 s (left) and z = 0 cross-sectional pressure distribution (right)
Fig. 13
Pressure distribution at t = 2.5 × 10−3 s (left) and z = 0 cross-sectional pressure distribution (right)
Close modal
Fig. 14
Velocity vector distribution at t = 5.0 × 10−4 s (left) and volume fraction distribution of cavitation (right)
Fig. 14
Velocity vector distribution at t = 5.0 × 10−4 s (left) and volume fraction distribution of cavitation (right)
Close modal
Fig. 15
Velocity vector distribution at t = 2.0 × 10−3 s (left) and volume fraction distribution of cavitation (right)
Fig. 15
Velocity vector distribution at t = 2.0 × 10−3 s (left) and volume fraction distribution of cavitation (right)
Close modal
Fig. 16
Velocity vector distribution at t = 2.5 × 10−3 s (left) and volume fraction distribution of cavitation (right)
Fig. 16
Velocity vector distribution at t = 2.5 × 10−3 s (left) and volume fraction distribution of cavitation (right)
Close modal

In Fig. 11, we can see pressure rising (to a maximum of approximately 650 MPa) in the direction of load application, that is, at the bottom of the analysis model. We adjusted maximum and minimum values for the color legend range to more clearly visualize the pressure distribution in the region (left figure). We can also confirm that because the piston pin moves rapidly, large negative pressure is generated in the upper part of the model, and cavitation occurs. The load applied to the piston pin increases up to t = 1.5 × 10−3 s, and the wall surface of the piston pin and connecting rod rapidly approach, causing the squeeze effect [66] in the Reynolds equation to appear prominently, increasing pressure in the lubricating oil. The increase in the pressure magnitude is also due to the decreased local film thickness caused by the load increase. We can also see that as the load increases, the pressure distribution increases locally in the central part of the model. This is likely caused by a local increase in the viscosity of the lubricating oil due to the piezo-viscosity effect as modeled by the Barus equation, in which viscosity exponentially increases as the pressure increases. This pressure increase allows the lubricating oil to support the load applied to the piston pin. However, the high-pressure region in the lubricating oil is limited to the center of the lower part of the model, and pressures around the edge are near atmospheric pressure, so no large pressure increases are observed there. As Fig. 14 shows, this occurs because lubricating oil is pushed out by the sudden movement of the piston pin. We can therefore expect the piston pin, which is more easily deformed than the connecting rod, to undergo a bow-like deformation.

Figures 12 and 13 show the analysis results for pressure distributions when load on the piston pin reaches the peak value of 4672 N, and the load gradually decreases. As the load decreases and the movement speed of the piston pin in the y direction decreases, the pressure increase effect of the lubricating oil due to the squeeze effect decreases, and the local pressure increase dulls the sharply pointed distribution. Because the eccentricity nearly constantly rotates, the wedge effect becomes stronger, and a pressure distribution similar to that under the half-Sommerfeld condition arises at t = 2.5 × 10−3 s. Figures 15 and 16 show that the amount of lubricating oil pushed out was small and that the oscillatory motion of the connecting rod formed a velocity field that rotated in the circumferential direction. We can also confirm that the cavitation region occurring in the upper part shrank with the load phenomenon and that cavitation sparsely occurred in the lower part of the model.

Figure 17 shows the initial state (at t = 0 s) and the gap width between the piston pin and connecting rod at t = 1.5 × 10−3 s, where the flow channel was narrowest. The flow channel width was 11.75 μm at the initial state but became as narrow as 1.8 μm due to negative displacement of the piston pin in the y-direction under loading. The aspect ratio was approximately 50 at the initial state and increased to approximately 290 at t = 1.5 × 10−3 s, but we confirmed the possibility of calculations that maintain numerical stability.

Fig. 17
Initial state of the lower part of the analysis model (left) and distance between the piston pin and connecting rod at t = 1.5 × 10−3 s (right)
Fig. 17
Initial state of the lower part of the analysis model (left) and distance between the piston pin and connecting rod at t = 1.5 × 10−3 s (right)
Close modal

3 Coupled Fluid–Structure Interaction Analysis of Lubricating Oil Flow Between the Piston and Small End of the Connecting Rod

In the previous section, we performed a rigid-body analysis on lubrication flow between a piston pin and the small end of a connecting rod to predict the location where seizure will occur if the structure does not deform. In Fig. 3, we also showed the validity of the current numerical analysis method through qualitative comparison with analytical solutions by the Reynolds equation. In this section, we extend the solver used in the previous section to perform numerical analysis that considers elastic deformation of the structure. Namely, we conduct coupled FSI analysis [6769] of lubricating oil flow between the piston and the small end of the connecting rod. By adding an energy equation and implementing a temperature coupling between the structure and the fluid [70], we also consider the thermal effects, which is an essential factor in seizure phenomena [36,44,7173]. Finally, using the results obtained from this solver, we propose a method for determining piston pin shapes that satisfy the required balance between seizure resistance and weight reduction.

3.1 Governing Equation and Numerical Method in the Fluid Model.

In this section, we reconsider lubricating oil as a compressible fluid and use the continuity Eq. (14) and the Navier–Stokes Eq. (15), which are respectively
ρt+(ρU)=0
(14)
and
ρUt+(ρUU)=p+τ
(15)
Furthermore, because this analysis also considers thermal effects, we introduce the energy equation
D(ρh)Dt=kT+τ:U+DpDt
(16)
In these equations, h is the enthalpy, k is the thermal conductivity, T is the temperature. The second and third terms on the right side of Eq. (16) respectively represent frictional heat (energy transport due to viscous dissipation) and compressive heat. These generated heats are transferred to the solid side via conjugate heat transfer (CHT) [74], resulting in thermal deformation of the solid.

3.1.1 Rheology Model Considering the Thermal Effect.

In the rigid-body analysis in the previous section, we created a rheological model for lubricating oil by using the pressure dependence of the viscosity coefficient by the Barus equation [14] and a non-Newtonian fluid model by Eyring. However, because temperature change affects the viscosity of lubricating oil, in this analysis, we extend the Barus equation and use the Houpert model [7578], which allows consideration of temperature. Viscosity in the Houpert model is
ηHoupert=ηBarusexp[β*(TT0)]
(17)
Here, ηBarus is the viscosity according to the Barus equation, T is the temperature, T0 is the basis temperature, and β* is defined as
β*=(ln(η0)+9.67)(1+5.1×109p)ZS0T0138
(18)
where η0 is the viscosity under atmospheric pressure, Z is defined as
Z=αpv5.1×109(ln(η0)+9.67)
(19)
and S0 is defined as
S0=β(T0138)ln(η0)+9.67
(20)
In the above equations, αpv is the pressure–viscosity coefficient, and β is the temperature–viscosity coefficient, which are characteristic values specific to the fluid.

3.2 Governing Equation in the Structural Model.

This analysis also considers deformation of the structure due to large loads. Assuming that the deformation is insufficient to cause plastic deformation, we consider only elastic deformation. Elastic deformation can be obtained by solving the linear momentum equation
ρs2Dt2+σ=fs
(21)
To consider the effect of thermal expansion [73,77] due to heat generated by the fluid, we also solve the following energy equation:
ρcv,sTt=(ksT)+σ:Dt
(22)
In these equations, subscript s denotes the structure side. ρs is the density of the structure, D is the displacement vector, σ is the stress tensor, fs is the body force, cv,s is the specific heat, and ks is the thermal conductivity. Using the linear elastic stress σLe and thermal stress σT, we can find the stress tensor as
σ=σLe+σT
(23)
The linear elastic stress σLe is calculated as
σLe=2μsϵ+λstr[ϵ]I
(24)
and the thermal stress σT is calculated as
σT=2(2μs+3λs)αs(TsT0)I
(25)
where ϵ is the strain tensor, αs is the volume expansion rate, Ts is the solid phase temperature, and T0 is the solid phase reference temperature
ϵ=12[D+(D)T]
(26)
μs and λs are Lame’s constants
μs=E2(1+ν(poi))
(27)
λs=ν(poi)E(1+ν(poi))(12ν(poi))
(28)
where E is the Young’s modulus, ν(poi) is the Poisson’s ratio.

3.3 Fluid–Structure Coupling.

Since the fluid and the structure are coupled, it is necessary to show that the boundary surface is continuous. In this analysis, we set continuous conditions for the boundary surface movement velocity, u, as
ui=us
(29)
for displacement amount as
Di=Ds
(30)
and for temperature as
Ti=Ts
(31)
where subscript i denotes the interfacial value. We also respectively define the exchange of stress and heat flux, q at the interface between the fluid and the structure as
niσi=nsσs
(32)
and
niqi=nsqs
(33)
where subscript n is the unit vector normal to the surface.

3.3.1 Solution by Strongly Coupled Partitioning.

There are two main methods for finding solutions in fluid–structure coupled analysis: the partitioned method [79] and the monolithic method (Fig. 18). The former method calculates fluids and structures using different solvers, while the latter uses a single matrix to simultaneously solve for fluids and structures. The monolithic method is generally considered to provide higher numeric stability, but there are many issues related to its extension to parallel computing with regional segmentation. In this analysis, the total number of mesh cells for the fluid and the structure is about 330,000, the aspect ratio is large, and the time-step cannot be increased, so computational loads are large. Parallel computing is thus indispensable, so this analysis is performed using the partitioned method.

Fig. 18
Fluid–structure coupled analysis
Fig. 18
Fluid–structure coupled analysis
Close modal

Partitioned solutions can be further subdivided into weakly and strongly coupled solutions. Weakly coupled analysis is a method that uses values at the boundary surface obtained from fluid analysis as the boundary condition on the structure side. However, this method is disadvantageous in that interactions between the fluid and the structure are considered in only one direction. This is inappropriate for weakly coupled analysis of phenomena such as EHL, where the structure is subjected to large forces. In this analysis, therefore, we attempt numerical analysis using strong coupling.

In strongly coupled analysis, we perform iterative calculations until the amount of deformation of the fluid and the structure converges between timesteps. In this analysis, we use the Aitken relaxation method, which was originally developed to simultaneously handle FSI and CHT [74]. The residual Resfsi for convergence testing of the fluid–structure interface is
Resfsi=DfbsolidDfb
(34)
where the amount of displacement at the boundary surface of the fluid region Dfb and the amount of displacement at the boundary surface of the solid region Dfbsolid are interpolated on the fluid side. We use the obtained boundary surface residual Resfsin and the immediately preceding residual Resfsi0 to find the relaxation factor
frelaxn=|frelaxrm0ΣResfsi0(ResfsinResfsi0)Σ|ResfsinResfsi0|2|
(35)
Using the obtained relaxation factor frelaxn, we can find the displacement amount of the fluid interface as follows:
Dfb=Dfb+frelaxnResfsi
(36)
Figure 19 shows the flow of analysis using the strongly coupled partitioning method.
Fig. 19
Analysis flow for strongly coupled partitioning
Fig. 19
Analysis flow for strongly coupled partitioning
Close modal

3.3.2 Boundary Lubrication Model.

When no load is applied to the piston pin, the wall surfaces of the piston pin and connecting rod are sufficiently separated, with the two parts separated by only lubricating oil. As the load increases, however, the gap between the piston pin and the connecting rod narrows. Then, microscale asperities of the metal surface begin to cause mechanical contact; this lubrication state is called boundary lubrication (Fig. 20). Mixed lubrication may also occur in the system, wherein hydrodynamic and asperity contact coexist. Predicting boundary lubrication is vital for predicting seizure in lubrication where large loads are generated, as in this analysis.

Fig. 20
Transition from fluid lubrication to boundary lubrication
Fig. 20
Transition from fluid lubrication to boundary lubrication
Close modal

3.3.3 The Greenwood–Tripp Model.

In this analysis, the region of the transition to boundary lubrication is determined according to the Greenwood–Tripp model [18,80,81]. In this model, asperities on the metal surface are assumed to be spherical, and by statistically handling the asperity height distribution, force due to mechanical contact can be obtained without using a mesh to reproduce the asperity shape in detail [82,83]. In this analysis, the distribution of asperity heights on the metal surface follows a Gaussian distribution, and the force pc generated by mechanical contact at the rough surface of the piston pin and connecting rod is calculated as
pc=KEGF(λ)
(37)
where constant K is given by
K=16152π(Naspβaspσ*)2σ*βasp
(38)
and constant EG is given by
EG=1(1ν(poi)12E1+1ν(poi)22E2)
(39)
In these equations, subscripts 1 and 2 indicate opposing metal surfaces, Nasp is the number of asperities per unit area, βasp is the radius of asperity curvature, σ* is calculated as [80,84]
σ*=0.7σ1+σ2
(40)
where σ1 and σ2 are standard deviations for asperity heights on the opposing metal surfaces. Assuming the height distribution of the asperities follows a Gaussian distribution, F(λ) in Eq. ( 37) can be calculated as [28,85]
F(λ)={4.4086×105(4λ)6.804(λ4)0(λ4)
(41)
where λ=hmt/σ*, and hmt is the distance between the two metal surfaces [85]. In other words, when the gap width becomes smaller and λ becomes 4 or less, we can judge that the asperities are in contact and the boundary lubrication state is reached. The current model does not allow analysis of full asperity contact and requires a thin lubricant film area. Numerical methods for conditions of complete asperity contact are currently under investigation.

3.4 Analysis Model and Numerical Analysis Conditions.

This section describes the analysis model and analysis conditions used for numerical analyses of EHL for piston pins and connecting rods. To compute the present system, we newly developed a custom solver named cavitatingFsiEHLFoam using the open source finite volume CFD code openfoam [86], which is in turn based on an extention of solids4Foam, an under-development openfoam solver [87,88]. This analysis requires solving for structural deformation, increasing computational loads over those for the analysis performed in the previous section, so we use a half model to reduce calculation times.

3.5 Fluid Analysis Model.

The fluid model used in this numerical analysis is very similar to the one described in the previous section, but as mentioned above, we use the half model (about 100,000 mesh cells) shown in Fig. 21 to reduce calculation times. Pressure and velocity boundary conditions are the same as in the analysis in Sec. 3. Also, this analysis considers temperatures, so we set the initial temperature condition as 350 K and apply the Neumann condition to the entire boundary (Table 3).

Fig. 21
Fluid mesh in the half model
Fig. 21
Fluid mesh in the half model
Close modal
Table 3

Boundary conditions for the fluid model

BoundaryPressureVelocityTemperature
InletConstant atmosphericContinuous in–out flowNeumann condition
OutletConstant atmosphericContinuous in–out flowNeumann condition
Piston pinNeumann conditionNo slipNeumann condition
Connecting rodNeumann conditionNo slipNeumann condition
Symmetry planeSymmetry conditionSymmetry conditionNeumann condition
BoundaryPressureVelocityTemperature
InletConstant atmosphericContinuous in–out flowNeumann condition
OutletConstant atmosphericContinuous in–out flowNeumann condition
Piston pinNeumann conditionNo slipNeumann condition
Connecting rodNeumann conditionNo slipNeumann condition
Symmetry planeSymmetry conditionSymmetry conditionNeumann condition

Table 4 shows physical property values for the lubricating oil used in the rheology model by Eyring–Houpert. SAE-Viscosity grade 10W–30 [65] was assumed for the physical properties of the oil.

Table 4

Physical characteristics of the lubricating oil [65]

LubricantVaporUnit
Viscosity1.0 × 10−58.97 × 10−6(m2/s)
Density8500.029(kg/m3)
Compressibility4.9 × 10−75.7 × 10−6(1/kPa)
Saturated vapor pressure5000(Pa)
Pressure–viscosity index0.689(–)
Thermo-viscosity index0.0476(–)
Eyring stress5.0 × 106(Pa)
Thermal conductivity0.150.025[W/(m · K)]
Heat capacity23001800[(J/(kg · K)]
LubricantVaporUnit
Viscosity1.0 × 10−58.97 × 10−6(m2/s)
Density8500.029(kg/m3)
Compressibility4.9 × 10−75.7 × 10−6(1/kPa)
Saturated vapor pressure5000(Pa)
Pressure–viscosity index0.689(–)
Thermo-viscosity index0.0476(–)
Eyring stress5.0 × 106(Pa)
Thermal conductivity0.150.025[W/(m · K)]
Heat capacity23001800[(J/(kg · K)]

3.6 Structural Analysis Model.

Figure 22 shows a numerical model of the piston pin–connecting rod structure, and Fig. 5 shows dimensions in the analysis model. There is a total of 224,000 mesh cells. Figure 23 shows the boundary names for both parts. The load applied to the piston pin in the combustion process described in Sec. 3 is applied to the wallTraction boundary, which connects the piston head and the piston pin. Table 5 shows structural displacement and the temperature boundary conditions. Both the piston pin and connecting rod are made of iron, and Table 6 shows material property values for the structure.

Fig. 22
Meshes for the piston pin structure (left) and the connecting rod (right)
Fig. 22
Meshes for the piston pin structure (left) and the connecting rod (right)
Close modal
Fig. 23
Structure boundary names
Fig. 23
Structure boundary names
Close modal
Table 5

Boundary conditions for the structure model

BoundaryDisplacementTemperature
ConnrodZero tractionNeumann condition
PinZero tractionNeumann condition
ConnrodWallsZero tractionNeumann condition
ConnrodWallBottomZero displacementNeumann condition
PinWallsZero tractionNeumann condition
ConnrodSymSymmetry conditionSymmetry condition
PinSymSymmetry conditionSymmetry condition
WallTractionDynamic tractionNeumann condition
BoundaryDisplacementTemperature
ConnrodZero tractionNeumann condition
PinZero tractionNeumann condition
ConnrodWallsZero tractionNeumann condition
ConnrodWallBottomZero displacementNeumann condition
PinWallsZero tractionNeumann condition
ConnrodSymSymmetry conditionSymmetry condition
PinSymSymmetry conditionSymmetry condition
WallTractionDynamic tractionNeumann condition
Table 6

Material property values for the structure

PropertyValueUnit
Density7850(kg/m3)
Young’s modulus209 × 109(Pa)
Poisson’s ratio0.3(–)
Thermal conductivity50(W/(m · K))
Thermal elasticity1.0 × 10−6(1/K)
Heat capacity500[(J/(kg · K)]
PropertyValueUnit
Density7850(kg/m3)
Young’s modulus209 × 109(Pa)
Poisson’s ratio0.3(–)
Thermal conductivity50(W/(m · K))
Thermal elasticity1.0 × 10−6(1/K)
Heat capacity500[(J/(kg · K)]

Table 7 shows information related to surface roughness in the Greenwood–Tripp model. Values in that table are based on measured data for the surface roughness of piston pin and conrod contact surfaces.

Table 7

Parameters related to structure surface roughness

PropertyValueUnit
Asperity density1.65 × 1010[1/m2]
Curvature radius2.81 × 10−5[m]
Standard deviation  2.0 × 10−6[m]
PropertyValueUnit
Asperity density1.65 × 1010[1/m2]
Curvature radius2.81 × 10−5[m]
Standard deviation  2.0 × 10−6[m]

3.7 Fluid-Side Analysis Results.

Figures 2426 show analysis results for the volume fraction of the vapor phase due to cavitation at each of the respective times above. Large negative pressure is generated in the upper part of the analysis region due to displacement of the piston pin in the negative direction of the y axis. In some cases, we can see that pressure in the lubricating oil was lower than the saturated vapor pressure, so cavitation occurred. Cavitation occurred only at the boundary on the piston pin side; no vapor phase was observed at the boundary on the connecting rod side. This vapor phase distribution was likely caused by sudden displacement of the piston pin due to a large load. In addition, cavitation can be seen at the connecting end of the oil supply hole at the upper part of the analysis region. Because cavitation occurs at the connecting end on the side of rotational direction of the oscillating motion, it is likely because pressure sharply drops when lubricating oil is drawn into the gap from the oil supply hole. As the cavitation region increases, viscosity of the lubricating oil decreases, and the compression rate increases, so there is a possibility that sufficient film thickness cannot be realized when the combustion process completes and an upward force is applied to the piston pin (Figs. 2729).

Fig. 24
Pressure distribution at t = 5.0 × 10−4 s (left) and pressure distribution at the z = 0 cross-section (right)
Fig. 24
Pressure distribution at t = 5.0 × 10−4 s (left) and pressure distribution at the z = 0 cross-section (right)
Close modal
Fig. 25
Pressure distribution at t = 2.0 × 10−3 s (left) and pressure distribution at the z = 0 cross-section (right)
Fig. 25
Pressure distribution at t = 2.0 × 10−3 s (left) and pressure distribution at the z = 0 cross-section (right)
Close modal
Fig. 26
Pressure distribution at t = 2.5 × 10−3 s (left) and pressure distribution at the z = 0 cross-section (right)
Fig. 26
Pressure distribution at t = 2.5 × 10−3 s (left) and pressure distribution at the z = 0 cross-section (right)
Close modal
Fig. 27
Cavitation distribution at t = 5.0 × 10−4 s
Fig. 27
Cavitation distribution at t = 5.0 × 10−4 s
Close modal
Fig. 28
Cavitation distribution at t = 2.0 × 10−3 s
Fig. 28
Cavitation distribution at t = 2.0 × 10−3 s
Close modal
Fig. 29
Cavitation distribution at t = 2.5 × 10−3 s
Fig. 29
Cavitation distribution at t = 2.5 × 10−3 s
Close modal

Figure 30 compares gap widths at the initial state and at t = 1.5 × 10−3 s, where load is maximized. In the rigid-body analysis, the minimum gap width narrowed to 1.8 μm and 9.95 μm as shown in Fig. 17, but in this analysis, the reduction was only 0.5 μm, for a minimum gap width of 11.25 μm. This is likely because elastic deformation of the structure suppresses the pressure increase, confirming that the gap width is overestimated by the rigid-body analysis alone, without considering elastic deformation.

Fig. 30
Piston pin–connecting rod separation at the bottom of the analysis model at the initial state (left) and at t = 1.5 × 10−3 s (right)
Fig. 30
Piston pin–connecting rod separation at the bottom of the analysis model at the initial state (left) and at t = 1.5 × 10−3 s (right)
Close modal

3.8 Structure-Side Analysis Results.

Figures 3133 show piston pin deformations and stress distributions. For ease of viewing, deformation amounts are magnified 100 times. The figures confirm that the piston pin exhibits a bow-like deformation as the load increases. The highest stress, approximately 500 MPa, is at the lower center of the piston pin. We can see that pressure in the lubricating oil was high at the lower center of the piston due to the squeeze effect, so the film thickness did not significantly change even when a large load was applied. Around the connecting rod edge, however, pushing the lubricating oil out suppresses pressure increase in the lubricating oil, and elastic deformation of the piston pin increases as the load increases. This is similar to the behavior predicted from the rigid-body analysis in Sec. 2.

Fig. 31
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 5.0 × 10−4 s (deformation amounts magnified 100×)
Fig. 31
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 5.0 × 10−4 s (deformation amounts magnified 100×)
Close modal
Fig. 32
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 2.0 × 10−3 s (deformation amounts magnified 100×)
Fig. 32
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 2.0 × 10−3 s (deformation amounts magnified 100×)
Close modal
Fig. 33
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 2.5 × 10−3 s (deformation amounts magnified 100×)
Fig. 33
Piston pin stress distribution (left) and x = 0 cross-section (right) at t = 2.5 × 10−3 s (deformation amounts magnified 100×)
Close modal

3.9 Discussion of the Piston Pin.

Figure 34 shows deformation of the piston pin and the connecting rod at t = 1.5 × 10−3 s. Deformation in both parts is magnified 100 times for ease of viewing.

Fig. 34
The piston pin and connecting rod at t = 1.5 × 10−3 s: (a) overhead view of the entire assembly and (b) cross-sectional view (deformation amounts magnified 100×)
Fig. 34
The piston pin and connecting rod at t = 1.5 × 10−3 s: (a) overhead view of the entire assembly and (b) cross-sectional view (deformation amounts magnified 100×)
Close modal

The film thickness is smallest at the lower end of the connecting rod when the piston pin exhibits a bow-like deformation. This part has the highest probability of seizure, and we can see that this is the same kind of seizure as the actual one shown in Fig. 1. We can thus say that this analysis qualitatively predicted the seizure location by calculation.

The amount of elastic deformation is much smaller in the connecting rod than in the piston pin, so deformation of the piston pin has the most significant effect on reducing the film thickness. In other words, minimizing bow-like deformation of the piston pin will reduce mechanical contact at the connecting rod end. An effective way of achieving this is to reduce the inner diameter of the piston pin. Piston pins generally have a hollow, cylindrical shape to reduce their weight, but while a larger inner radius makes the piston pin lighter, its rigidity against load is lessened. The analysis method developed in this research is expected to be very useful for numerical predictions of the minimal piston pin inner diameter that prevents seizure, allowing numerical analysis design using supercomputing without requiring high costs and long times for durability tests.

4 Conclusions

We conducted multiphase numerical analyses and investigations of thin-film fluid lubrication flow between a piston pin and connecting rod in an internal combustion engine, where there are particularly severe cavitating lubrication conditions at sliding parts with unsteady flow channel variation. We performed two analyses: a composite analysis with a rigid-body piston pin and a connecting rod, and a fluid–structure coupled analysis with a structure comprising both parts that can be elastically deformed. From these analyses, we obtained the following findings:

  • The analysis results for pressure distribution obtained from the composite rigid-body analysis were similar to those under the Sommerfeld condition, confirming the validity of this numerical analysis method and the analysis results. In the first half of the combustion process, when the load applied to the piston pin increases, pressure in the lubricating oil locally increases due to a squeeze effect resulting from rapid movement of the piston pin. At this time, movement of the piston pin caused large negative pressures on the side opposite of the high pressure, and cavitation was observed in some areas. Pressure increase in the lubricating oil was limited to the model’s central region, with no pressure increase observed around the edge.

  • We obtained the following findings from FSI analysis of elastic deformation of the piston pin and connecting rod. As in the rigid-body analysis, pressure in the lubricating oil increased in the loading direction due to the load applied to the piston pin, but the distribution was smoother than that in the rigid-body analysis. This is likely due to the effect of elastic deformation suppressing the pressure increase.

  • Fluid–structure coupled analysis showed bow-like deformation of the piston pin when a load is applied. We found that deformation of the connecting rod was smaller than that of the piston pin, and that the amount of piston pin deformation had a significant effect on the gap width between the two parts. We also showed that there is a high possibility of mechanical contact occurring at the connecting rod edge due to the bow-like deformation of the piston pin.

  • By comparison with a piston pin after actual use, we found that the fluid–structure coupled analysis qualitatively predicted the seizure location. One way to prevent piston pin seizure is to reduce its inner diameter to increase its rigidity. However, reducing the inner diameter leads to increased piston pin weight, so it is essential to use this numerical analysis method to select the optimal inner diameter. We showed that this numerical analysis method could predict seizure phenomena, which have been considered difficult to predict. This will allow some costly and time-consuming durability tests to be replaced with numerical analyses.

Acknowledgment

The authors would like to thank Dr. Zheng Zhang and Mr. Kiyoshi Yoshinaga (Softflow Co., Ltd., Japan) for their helpful discussions. This work was partly supported by the GCORE Program and the Collaborative Research Project of the IFS, Tohoku University. Part of the numerical results in this research was obtained using supercomputing resources at the Cyberscience Center, Tohoku University. The numerical simulations were partially performed using FUJITSU AFI-NITY scalar parallel computers at the Institute of Fluid Science, Tohoku University.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

Nomenclature

     
  • h =

    specific enthalpy, J/kg

  •  
  • k =

    thermal conductivity, W/(m · K)

  •  
  • n =

    unit vector normal to the surface

  •  
  • p =

    pressure, Pa

  •  
  • q =

    heat flux vector, W/m2

  •  
  • r =

    distance from the rigid-body, m

  •  
  • D =

    displacement, m

  •  
  • D =

    displacement vector, m

  •  
  • E =

    Young’s modulus, Pa

  •  
  • I =

    unit tensor

  •  
  • T =

    temperature, K

  •  
  • U =

    flow velocity vector, m/s

  •  
  • fs =

    body force, N/m3

  •  
  • Nasp =

    number of asperities per unit area

  •  
  • D/Dt =

    substantial derivative, 1/s

Greek Symbols

     
  • αpv =

    pressure–viscosity coefficient

  •  
  • αs =

    volume expansion rate 1/K

  •  
  • β =

    temperature–viscosity coefficient

  •  
  • βasp =

    radius of asperity curvature

  •  
  • ϵ =

    strain tensor

  •  
  • γ˙ =

    shear rate, Pa

  •  
  • η =

    viscosity coefficient, Pa · s

  •  
  • γ =

    displacement diffusion coefficient

  •  
  • μs, λs =

    Lame’s constants

  •  
  • μ =

    viscosity, Pa · s

  •  
  • ν(poi) =

    Poisson’s ratio

  •  
  • ρ =

    density, kg/m3

  •  
  • σLe =

    linear elastic stress, Pa

  •  
  • σ =

    stress tensor, Pa

  •  
  • σT =

    thermal stress, Pa

  •  
  • σ =

    normal stress, Pa

  •  
  • τ =

    viscous stress tensor, Pa

  •  
  • τ =

    shear stress, Pa

  •  
  • τ0 =

    Eyring stress, Pa

Subscripts

     
  • i =

    interface

  •  
  • l =

    liquid phase

  •  
  • s =

    structure side

  •  
  • v =

    vapor phase

  •  
  • 0 =

    reference value

  •  
  • sat =

    saturated vapor pressure

References

1.
Wellnitz
,
J.
,
Subic
,
A.
, and
Trufin
,
R.
,
2013
, “
Sustainable Automotive Technologies 2013
,”
Proceedings of the Fifth International Conference ICSAT 2013 (Lecture Notes in Mobility)
,
Ingolstadt, Germany
,
Sept. 25–27
,
Springer
,
New York
.
2.
Mihara
,
Y.
,
2016
, “
Research Trend and Latest Approach for Engine Tribology
,”
Engine Rev. Soc. Automotive Eng. Jpn.
,
6
(
3
), pp.
3
7
.
3.
Kurisu
,
T.
,
Kimura
,
S.
,
Shirai
,
H.
, and
Sugachika
,
N.
,
2015
, “
Tribology Analysis Technology Supports Engine Fuel Economy
,”
Mazda Tech. Rev.
,
32
, pp.
203
209
.
4.
Shi
,
F.
,
2011
, “
An Analysis of Floating Piston Pin
,”
SAE Int. J. Engines
,
4
(
1
), pp.
2100
2105
.
5.
Qiong
,
L.
,
Bin
,
M.
, and
Qinghua
,
Y.
,
2013
, “
Numerical Analysis and Study of Crankshaft System Hydrodynamic Lubrication and Dynamic Coupling System
,”
Appl. Math. Inform. Sci.
,
7
(
1L
), pp.
313
321
.
6.
Allmaier
,
H.
, and
Sander
,
D. E.
,
2020
,
Piston-Pin Rotation and Lubrication
,”
Lubricants
,
8
(
3
), p.
30
.
7.
Fang
,
C.
,
Meng
,
X.
,
Lu
,
Z.
,
Wu
,
G.
,
Tang
,
D.
, and
Zhao
,
B.
,
2019
, “
Modeling a Lubricated Full-Floating Pin Bearing in Planar Multibody Systems
,”
Tribol. Int.
,
131
, pp.
222
237
.
8.
Ferretti
,
A.
,
Giacopini
,
M.
,
Mastrandrea
,
L.
, and
Dini
,
D.
,
2018
, “
Investigation of the Influence of Different Asperity Contact Models on the Elastohydrodynamic Analysis of a Conrod Small-End/Piston Pin Coupling
,”
SAE Int. J. Engines
,
11
(
6
), pp.
919
934
.
9.
Rom
,
M.
, and
Müller
,
S.
,
2018
, “
An Effective Navier–Stokes Model for the Simulation of Textured Surface Lubrication
,”
Tribol. Int.
,
124
, pp.
247
258
.
10.
Rom
,
M.
, and
Müller
,
S.
,
2019
, “
A New Model for Textured Surface Lubrication Based on a Modified Reynolds Equation Including Inertia Effects
,”
Tribol. Int.
,
133
, pp.
55
66
.
11.
Wang
,
H.
,
Tang
,
L.
,
Zhou
,
C.
, and
Shi
,
Z.
,
2021
, “
Wear Life Prediction Method of Crowned Double Helical Gear Drive in Point Contact Mixed Elastohydrodynamic Lubrication
,”
Wear
,
484–485
, p.
204041
.
12.
Wang
,
H.
,
Zhou
,
C.
,
Lei
,
Y.
, and
Liu
,
Z.
,
2019
, “
An Adhesive Wear Model for Helical Gears in Line-Contact Mixed Elastohydrodynamic Lubrication
,”
Wear
,
426–427
, pp.
896
909
.
13.
Simon
,
V. V.
,
2019
, “
Improved Mixed Elastohydrodynamic Lubrication of Hypoid Gears by the Optimization of Manufacture Parameters
,”
Wear
,
438–439
, pp.
1
19
.
14.
Ghosh
,
M. K.
,
Majumdar
,
B. C.
, and
Sarangi
,
M.
,
2014
,
Fundamentals of Fluid Film Lubrication
,
McGraw-Hill Professional Pub
,
New York
.
15.
Habchi
,
W.
,
2018
,
Finite Element Modeling of Elastohydrodynamic Lubrication Problems
,
Wiley
,
London
.
16.
Peterson
,
W.
,
Singh
,
K.
, and
Sadeghi
,
F.
,
2022
,
Fluid–Solid Interaction Modeling of Elastohydrodynamic Lubrication Point Contacts
,”
ASME J. Tribol.
,
144
(
11
), p.
111601
.
17.
Hajishafiee
,
A.
,
Kadiric
,
A.
,
Ioannides
,
S.
, and
Dini
,
D.
,
2017
, “
A Coupled Finite-Volume CFD Solver for Two-Dimensional Elasto-Hydrodynamic Lubrication Problems With Particular Application to Rolling Element Bearings
,”
Tribol. Int.
,
109
, pp.
258
273
.
18.
Srirattayawong
,
S.
, and
Gao
,
S.
,
2012
, “
A CFD Study of the EHL Line Contact Problem With Consideration of the Surface Roughness Under Varied Loads
,”
HEFAT2012, Ninth International Conference on Heat Transfer, Fluid Mechanics and Thermodynamics
,
Malta
,
July 16–18
, pp.
723
729
.
19.
Ozasa
,
T.
,
Suzuki
,
S.
, and
Noda
,
T.
,
1996
, “
Elastohydrodynamic Lubrication Analysis of Con-Rod Big-End Bearings
,”
Toyota Central R&D Labs. R&D Rev.
,
31
(
4
), pp.
1
12
.
20.
Ba
,
L.
,
He
,
Z.-P.
,
Liu
,
Y.-H.
, and
Zhang
,
G.-C.
,
2015
, “
Analysis of Piston-Pin Lubrication Considering the Effects of Structure Deformation and Cavitation
,”
J. Zhejiang Univ. Sci. A
,
161
(
6
), pp.
443
463
.
21.
Dowson
,
D.
,
Godet
,
M.
, and
Taylor
,
C.
,
1974
, “
Cavitation and Related Phenomena in Lubrication
,”
Proceedings of the First Leeds-Lyon Symposium on Tribology held in the Institute of Tribology
,
Department of Mechanical Engineering, the University of Leeds, England, Institute of Tribology, Leeds University and the Institut National des Sciences Appliquées de Lyon
, pp.
1
238
.
22.
Elrod
,
H. G.
,
1981
, “
A Cavitation Algorithm
,”
ASME J. Lubr. Tech.
,
103
(
3
), pp.
350
354
.
23.
Liu
,
W.
,
Ni
,
H.
,
Wang
,
P.
, and
Chen
,
H.
,
2020
, “
Investigation on the Tribological Performance of Micro-Dimples Textured Surface Combined With Longitudinal or Transverse Vibration Under Hydrodynamic Lubrication
,”
Int. J. Mech. Sci.
,
174
, p.
105474
.
24.
Novotný
,
P.
,
Skara
,
P.
, and
Hliník
,
J.
,
2018
, “
The Effective Computational Model of the Hydrodynamics Journal Floating Ring Bearing for Simulations of Long Transient Regimes of Turbocharger Rotor Dynamics
,”
Int. J. Mech. Sci.
,
148
, pp.
611
619
.
25.
Taylor
,
C.
,
1975
, “
Gaseous Cavitation in Lightly Loaded Liquid Film Journal Bearings
,”
Int. J. Mech. Sci.
,
17
(
3
), pp.
177
185
.
26.
Meged
,
Y.
,
Venner
,
C.
, and
ten Napel
,
W.
,
1995
, “
Classification of Lubricants According to Cavitation Criteria
,”
Wear
,
186–187
(Part 2), pp.
444
453
.
27.
Choi
,
J.-K.
, and
Chahine
,
G. L.
,
2016
, “
Relationship Between Material Pitting and Cavitation Field Impulsive Pressures
,”
Wear
,
352–353
, pp.
42
53
.
28.
Zhao
,
B.
,
Hu
,
X.
,
Li
,
H.
,
Si
,
X.
,
Dong
,
Q.
,
Zhang
,
Z.
, and
Zhang
,
B.
,
2022
, “
A New Approach for Modeling and Analysis of the Lubricated Piston Skirt-Cylinder System With Multi-Physics Coupling
,”
Tribol. Int.
,
167
, p.
107381
.
29.
Tahmasebi
,
E.
,
Lucchini
,
T.
,
Errico
,
G. D.
,
Onorati
,
A.
, and
Hardy
,
G.
,
2017
, “
An Investigation of the Validity of a Homogeneous Equilibrium Model for Different Diesel Injector Nozzles and Flow Conditions
,”
Energy Convers. Manage.
,
154
, pp.
46
55
.
30.
Dhande
,
D.
, and
Pande
,
D.
,
2018
, “
Multiphase Flow Analysis of Hydrodynamic Journal Bearing Using CFD Coupled Fluid Structure Interaction Considering Cavitation
,”
J. King Saud Univ. – Eng. Sci.
,
30
(
4
), pp.
345
354
.
31.
Mithun
,
M.-G.
,
Koukouvinis
,
P.
, and
Gavaises
,
M.
,
2018
, “
Numerical Simulation of Cavitation and Atomization Using a Fully Compressible Three-Phase Model
,”
Phys. Rev. Fluids
,
3
, p.
064304
.
32.
Bicer
,
B.
, and
Sou
,
A.
,
2016
, “
Application of the Improved Cavitation Model to Turbulent Cavitating Flow in Fuel Injector Nozzle
,”
Appl. Math. Model.
,
40
(
7
), pp.
4712
4726
.
33.
Ishimoto
,
J.
,
Sato
,
F.
, and
Sato
,
G.
,
2010
, “
Computational Prediction of the Effect of Microcavitation on an Atomization Mechanism in a Gasoline Injector Nozzle
,”
ASME J. Eng. Gas Turbines Power
,
132
(
8
), p.
082801
.
34.
Schmidt
,
D. P.
,
Rutland
,
C. J.
, and
Corradini
,
M. L.
,
1999
, “
A Fully Compressible, Two-Dimensional Model of Small, High-Speed, Cavitating Nozzles
,”
Atomization Sprays
,
9
(
3
), pp.
255
276
.
35.
Feldermann
,
A.
,
Neumann
,
S.
, and
Jacobs
,
G.
,
2017
, “
CFD Simulation of Elastohydrodynamic Lubrication Problems With Reduced Order Models for Fluid–Structure Interaction
,”
Tribology–Mater. Surfaces Interfaces
,
11
(
1
), pp.
30
38
.
36.
Bair
,
S. S.
,
2019
,
High Pressure Rheology for Quantitative Elastohydrodynamics
, 2nd ed.,
Elsevier Science
,
Amsterdam
.
37.
Venner
,
C.
, and
Bos
,
J.
,
1994
, “
Effects of Lubricant Compressibility on the Film Thickness in Ehl Line and Circular Contacts
,”
Wear
,
173
(
1
), pp.
151
165
.
38.
Dowson
,
D.
,
1995
, “
Elastohydrodynamic and Micro-Elastohydrodynamic Lubrication
,”
Wear
,
190
(
2
), pp.
125
138
.
39.
Tsann-Rong
,
L.
, and
Jen-Fin
,
L.
,
1991
, “
Compressible Elastohydrodynamic Lubrication of Rolling and Sliding Contacts With a Power Law Fluid
,”
Wear
,
142
(
2
), pp.
315
330
.
40.
Dellacherie
,
S.
,
Faccanoni
,
G.
,
Grec
,
B.
, and
Penel
,
Y.
,
2019
, “
Accurate Steam-Water Equation of State for Two-Phase Flow Lmnc Model With Phase Transition
,”
Appl. Math. Model.
,
65
, pp.
207
233
.
41.
Jacobson
,
B. O.
,
2012
,
Rheology and Elastohydrodynamic Lubrication
(
Tribology Series, Vol. 19
),
Elsevier Science
,
Amsterdam
.
42.
Xiang
,
G.
,
Han
,
Y.
,
Wang
,
J.
,
Wang
,
J.
, and
Ni
,
X.
,
2019
, “
Coupling Transient Mixed Lubrication and Wear for Journal Bearing Modeling
,”
Tribol. Int.
,
138
, pp.
1
15
.
43.
Chnafa
,
C.
,
Mendez
,
S.
, and
Nicoud
,
F.
,
2014
, “
Image-Based Large-Eddy Simulation in a Realistic Left Heart
,”
Comput. Fluids
,
94
, pp.
173
187
.
44.
Nguyen-Schaefer
,
H.
,
2018
,
Computational Design of Rolling Bearings
,
Springer
,
Switzerland
.
45.
Biboulet
,
N.
, and
Lubrecht
,
A.
,
2018
, “
Efficient Solver Implementation for Reynolds Equation With Mass-Conserving Cavitation
,”
Tribol. Int.
,
118
, pp.
295
300
.
46.
Qiu
,
Y.
, and
Khonsari
,
M. M.
,
2009
,
On the Prediction of Cavitation in Dimples Using a Mass-Conservative Algorithm
,”
ASME J. Tribol.
,
131
(
4
), p.
041702
.
47.
Jakobsson
,
B.
, and
Floberg
,
L.
,
1957
, “
The Finite Journal Bearing, Considering Vaporization
,”
Trans. Chalmers Univ. Technol.
,
190
(
3
), pp.
1
116
.
48.
Floberg
,
L.
,
1957
, “
The Infinite Journal Bearing, Considering Vaporization
,”
Trans. Chalmers Univ. Technol.
,
189
(
2
), pp.
1
83
.
49.
Olsson
,
K.
,
1965
, “
Cavitation in Dynamically Loaded Bearings
,”
Trans. Chalmers Univ. Technol.
,
308
(
26
), pp.
1
60
.
50.
Etsion
,
I.
,
Kligerman
,
Y.
, and
Halperin
,
G.
,
1999
, “
Analytical and Experimental Investigation of Laser-Textured Mechanical Seal Faces
,”
Tribol. Trans.
,
42
(
3
), pp.
511
516
.
51.
Ausas
,
R.
,
Ragot
,
P.
,
Leiva
,
J.
,
Jai
,
M.
,
Bayada
,
G.
, and
Buscaglia
,
G. C.
,
2007
, “
The Impact of the Cavitation Model in the Analysis of Microtextured Lubricated Journal Bearings
,”
ASME J. Tribol.
,
129
(
4
), pp.
868
875
.
52.
Payvar
,
P.
, and
Salant
,
R. F.
,
1992
, “
A Computational Method for Cavitation in a Wavy Mechanical Seal
,”
ASME J. Tribol.
,
114
(
1
), pp.
199
204
.
53.
Bayada
,
G.
, and
Chupin
,
L.
,
2013
,
Compressible Fluid Model for Hydrodynamic Lubrication Cavitation
,”
ASME J. Tribol.
,
135
(
4
), p.
041702
.
54.
Wang
,
Y.
,
Wang
,
Q. J.
, and
Lin
,
C.
,
2003
, “
Mixed Lubrication of Coupled Journal-Thrust-Bearing Systems Including Mass Conserving Cavitation
,”
ASME J. Tribol.
,
125
(
4
), pp.
747
755
.
55.
Shahmohamadi
,
H.
,
Rahmani
,
R.
,
Rahnejat
,
H.
,
Garner
,
C.
, and
Dowson
,
D.
,
2015
,
Big End Bearing Losses With Thermal Cavitation Flow Under Cylinder Deactivation
,”
Tribol. Lett.
,
57
(
1
), pp.
1
17
.
56.
Schweizer
,
B.
,
2009
,
Numerical Approach for Solving Reynolds Equation With JFO Boundary Conditions Incorporating ALE Techniques
,”
ASME J. Tribol.
,
131
(
1
), p.
011702
.
57.
Kango
,
S.
,
Sharma
,
R.
, and
Pandey
,
R.
,
2014
, “
Thermal Analysis of Microtextured Journal Bearing Using Non-Newtonian Rheology of Lubricant and Jfo Boundary Conditions
,”
Tribol. Int.
,
69
, pp.
19
29
.
58.
Harp
,
S. R.
, and
Salant
,
R. F.
,
2000
, “
An Average Flow Model of Rough Surface Lubrication With Inter-Asperity Cavitation
,”
ASME J. Tribol.
,
123
(
1
), pp.
134
143
.
59.
Flores
,
P.
,
Ambrósio
,
J.
,
Claro
,
J. C. P.
,
Lankarani
,
H. M.
, and
Koshy
,
C. S.
,
2009
, “
Lubricated Revolute Joints in Rigid Multibody Systems
,”
Nonlinear Dyn.
,
56
, pp.
277
295
.
60.
Snyder
,
T. A.
,
Braun
,
M. J.
, and
Pierson
,
K. C.
,
2016
, “
Two-Way Coupled Reynolds and Rayleigh–Plesset Equations for a Fully Transient, Multiphysics Cavitation Model With Pseudo–Cavitation
,”
Tribol. Int.
,
93
, pp.
429
445
.
61.
Chasalevris
,
A.
, and
Sfyris
,
D.
,
2013
, “
Evaluation of the Finite Journal Bearing Characteristics, Using the Exact Analytical Solution of the Reynolds Equation
,”
Tribol. Int.
,
57
, pp.
216
234
.
62.
Weller
,
H. G.
, and
Tabora
,
G.
,
1998
, “
A Tensorial Approach to Computational Continuum Mechanics Using Object-Oriented Techniques
,”
Comput. Phys.
,
12
(
6
), pp.
620
631
.
63.
FunctionBay
,
2022
. https://functionbay.com/.
64.
Gummer
,
A.
, and
Sauer
,
B.
,
2014
, “
Modeling Planar Slider-Crank Mechanisms With Clearance Joints in RecurDyn
,”
Multibody Syst. Dyn.
,
31
, pp.
127
145
.
66.
Yasuna
,
J. A.
, and
Hughes
,
W. F.
,
1994
, “
Squeeze Film Dynamics of Two-Phase Seals: Part II-Turbulent Flow
,”
ASME J. Tribol.
,
116
(
3
), pp.
479
487
.
67.
Balcombe
,
R.
,
Fowell
,
M. T.
,
Olver
,
A. V.
,
Ioannides
,
S.
, and
Dini
,
D.
,
2011
, “
A Coupled Approach for Rolling Contact Fatigue Cracks in the Hydrodynamic Lubrication Regime: The Importance of Fluid/Solid Interactions
,”
Wear
,
271
(
5
), pp.
720
733
.
68.
Luo
,
L.
,
Jiang
,
Z.
,
Wei
,
D.
, and
Jia
,
F.
,
2021
, “
A Study of Influence of Hydraulic Pressure on Micro-Hydromechanical Deep Drawing Considering Size Effects and Surface Roughness
,”
Wear
,
477
, p.
203803
.
69.
Wang
,
H.
,
Yu
,
Y.
,
Yu
,
J.
,
Xu
,
W.
,
Li
,
X.
, and
Yu
,
S.
,
2019
, “
Numerical Simulation of the Erosion of Pipe Bends Considering Fluid-Induced Stress and Surface Scar Evolution
,”
Wear
,
440–441
, p.
203043
.
70.
Hartinger
,
M.
,
Dumont
,
M.-L.
,
Ioannides
,
S.
,
Gosman
,
D.
, and
Spikes
,
H.
,
2008
, “
CFD Modeling of a Thermal and Shear-Thinning Elastohydrodynamic Line Contact
,”
ASME J. Tribol.
,
130
(
4
), p.
041503
.
71.
Popp
,
A.
, and
Wriggers
,
P.
,
2018
,
Contact Modeling for Solids and Particles
(
CISM International Centre for Mechanical Sciences, Vol. 585
),
Springer International Publishing
,
Switzerland
.
72.
Wriggers
,
P.
,
2006
,
Computational Contact Mechanics
,
Springer-Verlag
,
Berlin, Heidelberg
.
73.
Eslami
,
M. R.
,
Hetnarski
,
R. B.
,
Ignaczak
,
J.
,
Noda
,
N.
,
Sumi
,
N.
, and
Tanigawa
,
Y.
,
2015
,
Theory of Elasticity and Thermal Stresses: Explanations, Problems and Solutions
(
Solid Mechanics and Its Applications
),
Springer Netherlands
,
New York
.
74.
Välikangas
,
T.
, “
Conjugate Heat Transfer in OpenFOAM
,”
Chalmers University of Technology
, http://www.tfd.chalmers.se/ hani/kurser/OS˙CFD˙2016.
75.
Houpert
,
L.
,
1985
, “
New Results of Traction Force Calculations in Elastohydrodynamic Contacts
,”
ASME J. Tribol.
,
107
(
2
), pp.
241
245
.
76.
Bruyere
,
V.
,
Fillot
,
N.
,
Morales-Espejel
,
G.
, and
Vergne
,
P.
,
2012
, “
Computational Fluid Dynamics and Full Elasticity Model for Sliding Line Thermal Elastohydrodynamic Contacts
,”
Tribol. Int.
,
46
(
1
), pp.
3
13
.
77.
Habchi
,
W.
, and
Bair
,
S.
,
2013
, “
Quantitative Compressibility Effects in Thermal Elastohydrodynamic Circular Contacts
,”
ASME J. Tribol.
,
135
(
1
), p.
011502
.
78.
Gamaniel
,
S.
,
Dini
,
D.
, and
Biancofiore
,
L.
,
2021
, “
The Effect of Fluid Viscoelasticity in Lubricated Contacts in the Presence of Cavitation
,”
Tribol. Int.
,
160
, p.
107011
.
79.
Profito
,
F.
,
Zachariadis
,
D.
, and
Dini
,
D.
,
2019
, “
Partitioned Fluid–Structure Interaction Techniques Applied to the Mixed-Elastohydrodynamic Solution of Dynamically Loaded Connecting-Rod Big-End Bearings
,”
Tribol. Int.
,
140
, p.
105767
.
80.
Greenwood
,
J.
, and
Tripp
,
J.
,
1970
, “
The Contact of Two Nominally Flat Rough Surfaces
,”
Proc. Inst. Mech. Eng.
,
185
(
48/71
), pp.
625
633
.
81.
Skuric
,
V.
,
Jaeger
,
P. D.
, and
Jasak
,
H.
,
2018
, “
Lubricated Elastoplastic Contact Model for Metal Forming Processes in Openfoam
,”
Comput. Fluids
,
172
, pp.
226
240
.
82.
Kudish
,
I. I.
,
2017
,
Elastohydrodynamic Lubrication for Line and Point Contacts: Asymptotic and Numerical Approaches
, 1st ed.,
CRC Press
,
Boca Raton, FL
.
83.
Lin
,
L.
,
2013
, “
Assessment of Effects of Surface Roughness and Oil Viscosity on Friction Coefficient Under Lubricated Rolling–Sliding Conditions
,” KOMATSU Technical Report, Vol.
55
, no.
166
, pp.
1
7
.
84.
Akalin
,
O.
, and
Newaz
,
G. M.
,
1999
, “
Piston Ring-Cylinder Bore Friction Modeling in Mixed Lubrication Regime: Part I—Analytical Results
,”
ASME J. Tribol.
,
123
(
1
), pp.
211
218
.
85.
Hu
,
Y.
,
Cheng
,
H. S.
,
Arai
,
T.
,
Kobayashi
,
Y.
, and
Aoyama
,
S.
,
1994
, “
Numerical Simulation of Piston Ring in Mixed Lubrication – A Nonaxisymmetrical Analysis
,”
ASME J. Tribol.
,
116
(
3
), pp.
470
478
.
86.
OpenFOAM
,
2021
. https://www.openfoam.com/.
87.
Cardiff
,
P.
,
Karac
,
A.
,
Jaeger
,
P. D.
,
Jasak
,
H.
,
Nagy
,
J.
,
Ivankovic
,
A.
, and
Tukovic
,
Z.
,
2018
, “
An Open-Source Finite Volume Toolbox for Solid Mechanics and Fluid–Solid Interaction Simulations
,” arXiv:1808.10736.
88.
Cardiff
,
P.
,
2021
, “
Solid Mechanics and Fluid–Solid Interaction Using the Solids4foam Tool Box
,” https://www.researchgate.net/profile/Philip-Cardiff.