Abstract
We consider the case of soft contacts in mixed lubrication conditions. We develop a novel, two scales contact algorithm in which the fluid- and asperity-asperity interactions are modeled within a deterministic or statistic scheme depending on the length scale at which those interactions are observed. In particular, the effects of large-scale roughness are deterministically calculated, whereas those of small-scale roughness are included by solving the corresponding homogenized problem. The contact scheme is then applied to the modeling of dynamic seals. The main advantage of the approach is the tunable compromise between the high-computing demanding characteristics of deterministic calculations and the much lower computing requirements of the homogenized solutions.
1. Introduction
Compliant contacts, most commonly known as soft-contacts, are very common in nature (e.g., cartilage lubrication, eye-eyelid contact) and technology (e.g., tires, rubber sealings, adhesives). It has long been stated that the friction and fluid leakage characteristics of wet soft-contacts are strongly related, among the other factors, to the local interactions occurring at the contact interface [1–5]. In the case of randomly rough surfaces, the basic understanding of the role played by the asperity-asperity and fluid-asperity interactions, occurring over a wide range of roughness length-scales, has been largely investigated and debated in the very recent scientific literature [6–12]. Given the (usual) fractal nature of random roughness, a number of interesting phenomena have been highlighted, as for example, the viscous-hydroplaning [6], the viscous flattening [9–15], the fluid-induced roughness anisotropic deformation [10, 11], the local [10, 11] and global [8, 16] fluid entrapment, and many others. The way to deal with random roughness contact mechanics, despite being nontrivial and suffering of a certain description fragmentation, is however well described in the current scientific literature.
On the other side, nowadays bio-inspired research [17, 18], together with the widely-spreading practice of surface engineering [19], is showing the many (mainly unexplored) opportunities offered by the physical-chemical ordered modification of surfaces in order to tailor targeted macroscopic contact characteristics, such as adhesion and friction. Bio-inspired adhesive research [20] is probably the best state of the art example of such research trend. However, investigating the combined effect of, let us say, quantized roughness and fluid action has not equally attracted the scientific community attention, apart from few experimental investigations [12, 21, 22] and basic theoretical investigations [23, 24]. This may be justified by the complexity of the numerical formulation of the problem, which is expected to not to present an analytical treatment. As a result, and to the best of authors’ knowledge, the combined effect of lubricant action and single-scale (contact splitter) roughness has not been practically investigated by tribologists’ community.
In this work we give our personal contribution to this research field and, in particular, we discuss a novel numerical scheme of soft mixed lubrication, which can be adopted to perform such investigations. In our model the roughness is split into two contributions. A threshold scale (where is the representative macroscopic size of the contact and is the threshold roughness length-scale) is identified. For length scales , that is , the system is investigated by using a deterministic approach, whereas for , that is, , the problem is treated by homogenizing the equations. This is of utmost help in performing numerical simulations, as for example in the case of microstructured or bio-inspired surfaces where the surface geometry roughness is characterized by a single scale texture (e.g., a pillars array) combined with random roughness at sub-micrometer scales. The main advantage of the proposed approach is, therefore, the strong reduction of numerical complexity. The paper is organized as follows. In Section 2 we describe the mixed lubrication model scheme (we refer to appendix for the details of numerics), whereas in Section 3 we report an example of application of the proposed model to the case of a lip sealing geometry operating in steady-sliding contact. For a detailed description on dynamic sealings modeling, the reader is referred to [25–29].
2. Problem Formulation
Here we consider a generic rough compliant solid in steady sliding contact with a rigid smooth counter surface, as shown in Figure 1(a). Due to the coupled action of asperity-fluid and asperity-asperity interactions, the compliant surface gains the actual (or deformed) configuration (Figure 1(b)). In soft-contacts, the shape difference between the deformed and initial configurations cannot be usually neglected, due to the occurrence of large deformations. Nevertheless, the small displacement assumption is usually adopted in the literature of soft-elastohydrodynamics [6, 30, 31] (soft-EHL), supported by a relevant experimental validation for ball-on-flat contact geometries. However, if one analyses the contact at much shorter length scales, that is, at the asperities length scales, then one finds out that the asperities may be completely flattened because of the high contact and fluid pressures. Hence, neglecting the influence of large deformation would lead to strong inaccuracy in describing the evolution of the system at the micro-scales. For this reason, in this work, the small displacement assumption has been relaxed. Moreover, in some cases, for example, for dynamic sealings modeling, the precise calculation of the large tangential displacement of the rubber at the interface of the sliding contact is a must in order to capture key phenomena as the well known reverse pumping effect [26].
(a) Initial configuration
(b) Deformed configuration
Consider now the schematic drawing of Figure 1(a). We assume that the whole texture belongs to the class of Reynolds roughness, that is, with (where is the surface height distribution). In such a case, the thin film lubrication formulation is sufficiently accurate in describing the fluid dynamics at the sliding interface at all roughness lengths. The (ensemble) average local separation, is calculated as: where and (the subscript or 2 or 3, see Figure 1(b)) are the three coordinates of the lip surface points at the initial () and at the deformed () state respectively. The quantity represents the macroscopic contact shape and the undeformed low-pass filtered (i.e., for ) surface roughness. is the generic component of the displacement vector describing the average local surface deformation of the compliant body. In the classical EHL approach the displacement is calculated within linear elasticity framework [31], by adopting a boundary element approach which requires operation where is the number of discretization points. In our case, however, we adopt a more general nonlinear rheology and the displacement components are calculated by employing a classical finite-element solver, which requires operations. Observe that the normal (locally averaged) stress is (see also [6]) where , is the locally averaged fluid pressure and the locally average solid contact pressure, the latter coming from the asperity-asperity interactions occurring at length-scales . We observe that [6]: where is the local normalized area of solid contact, where is the average fluid shear stress (see later for more details). We assume that the solid friction shear stress is constant (and directed along the sliding direction), as it happens in the case for a rubber-inert substrate contact [6]. Clearly, the model can be easily extended to include different boundary friction conditions.
The relation between the average solid contact pressure and the local average interfacial separation is obtained from Persson’s theory of contact mechanics [32, 33]. In particular, given the local elastic properties of the material and the power spectral density (PSD, see, e.g., Figure 2) of the surface ( is the surface roughness field, with average value , , where is the low frequency cut-off wavelength), the theory allows to calculate the local areal fraction of solid-solid contact as a function of the locally averaged solid-solid contact pressure as where the quantity is the root mean square gradient of the rough surface and is related to the PSD of the rough surface through the formula . The relation between the local interfacial separation and the solid-solid contact pressure is instead determined, in the adhesionless case, by requiring that the change of elastic energy per unit nominal contact area at the interface equals the work done by the contact pressure to deform the elastic body: In (6) the elastic energy per unit contact area is calculated as a function of the contact pressure (see [32, 33] for more details). The final formula linking the local separation and the solid-solid contact pressure is relatively complicated, but at large separation it simplifies and becomes [33] where the quantities and can be easily calculated once the PSD of the rough surface is known [33], and .
We consider now the case where a Newtonian fluid is sandwiched at the interface between the solids. We assume constant viscosity , constant density , and isothermal conditions. The adoption of a Reynolds roughness, and the assumption of a representative average interfacial separation value , suggests that the fluid velocity varies slowly with the coordinates and compared to the variation in the orthogonal direction . In such a case, combining the equilibrium with the continuity equation, the homogenized problem formulation reads [9, 34] where, as before, and are the locally averaged (over a length scale given by the longest wavelength surface roughness component, i.e., ) interfacial surface and fluid pressure, respectively. The fluid flow conductivity (tensor) and can be related, respectively, to the pressure flow factor (tensor) [7, 9, 35, 36] and to the shear flow factor (tensor) , which we assume to be a function of the average interfacial separation only. Here has been calculated on the basis of the Bruggeman’ effective medium theory and on the Persson’s contact mechanics [34]. Within this approach, the effect of solid contacts percolation, as well as of local fluid trapping, on the fluid flow can be taken approximately into account, as recently discussed in [16]. For simplicity, we have assumed , and , where is the identity tensor. Moreover, by considering the steady sliding condition, the homogenized lubricant equation simplifies into: where is the sliding velocity.
However, in order to include the effect of micro-cavitations occurring at large-scale roughness lengths, a JFO cavitation model (constant mixture pressure in the cavitation zone) is adopted. The Reynolds equation is then reformulated under a mass conservative equation valid throughout the cavitating/not cavitating domain: where we have adopted the cavitation index ( in the cavitation zones, and otherwise) and the dummy variable defined as where is the lubricant density, and is the representative solid shear elastic modulus. Note that (10) must be solved on the deformed configuration, which is an unknown of the problem: therefore, reformulating the fluid equation in the initial configuration may result numerically convenient. Thus (10) has been rephrased by using the mapping rule The Jacobian of the transformation can be calculated as the inverse of the Jacobian of the transformation Observe that the determinant of the Jacobian tensors must be necessarily larger than zero, that is, . The above transformation enables the generation of an adaptive mesh grid which follows the changing of the surface shape during the deformation process, thus without loosing spatial resolution. The detailed numerics derivation is reported in the Appendix.
The boundary conditions to be applied in the resolution of (10) depend, clearly, on the particular soft-contact problem under investigation. The dimensionless formulation of (10) is where and (note that has been made dimensionless with , whereas other lengths with ). corresponds in this case to the largest deterministic roughness wavelength.
3. Results
In this work we use the proposed model to investigate the lubricated contact of the lip seal schematically shown in Figure 3. Moreover, being much smaller than the shaft diameter, it is possible to reduce the computational domain to only a small angular fraction of the lip surface. Therefore, (10) is solved with the constant pressure boundary conditions at the low pressure side and at the high pressure side , and with the periodicity conditions on the circumferential -direction, with spatial period . The two-scale mixed lubrication model has been numerically solved, as described in the Appendix. We have opted for a fractal self-affine isotropic geometry. For any self-affine fractal surface the statistical properties are invariant under the transformation In such a case it can be shown that for isotropic surface the PSD is where , , is the Hurst exponent of the randomly rough profile, which is related to the fractal dimension . To apply our method we need to numerically generate the surface in the low frequencies spectrum, that is, for , see, for example, Figure 2. To this end we have utilized the spectral method described in [37] where the roughness is described by a periodic surface in the form of Fourier series where , . Since is real we must have . Moreover for randomly rough surfaces the following relation must be satisfied with , , where the symbol is the ensemble average operator. The PSD of Surface Equation (17) is from which it follows: For isotropic surfaces we have which simply gives and assuming self-affine fractal surface (see (16)) one obtains Hence the quantities can be determined once and the Hurst exponent of the fractal surface are known. However to completely characterize the rough profile we still need the probability distribution of the quantities . We first observe that the condition with , is satisfied if the phases of the complex quantities are random numbers uniformly distributed between and . We also recall the condition also implies that and that the quantities . So what we need now is only the probability distribution of . Of course there are several choices and the simplest one is to assume that the probability density function of is just a Dirac’s delta function centered at , that is, It can be shown that this choice guarantees also that the random profile has a Gaussian random distribution.
3.1. Sealing at Nearly Zero Sliding Velocity
Here we report on the model application to the case where the sliding velocity between the two mating surfaces is vanishing (), that is, for the case of static seals. The self-affine roughness of the lip surface presents a long-distance cut-off frequency , fractal dimension , , and small scale cut-off frequency . The threshold frequency adopted in the calculation is . In the following, the seal is assumed linear elastic. The calculation is performed at a constant rigid penetration of the lip. In Figure 4 we show the locally averaged interfacial separation field (note that is the circumferential direction, whereas is the axial direction). Interestingly the corrugation which can be observed in the figure is just a consequence of the deterministically included roughness (), at larger frequency the surface appears smooth since the high frequency () content of the roughness has been included through Persson’s statistic model, that is, by means of homogenization. The corresponding average solid contact pressure is shown in Figure 5(a). Observe first that, due to the presence of roughness, the contact is split into many contact patches in the whole lip apparent contact area and, correspondingly, normal stresses are concentrated into contact spots. Observe also that the average solid contact pressure is smooth at the contact borders, differently, from what instead expected in the case of smooth elastic contact (e.g., in the case of the Hertzian contact). This is due to the homogenized roughness contribution, which distributes, on a wider local contact area, the pressure acting on the single deterministic asperity. The average fluid pressure is shown in Figure 5(b). Observe that the fluid pressure presents steep variations (almost step-like) in those locations where the separation between the surfaces takes the smallest values. This is in perfect agreement with the critical junction theory of leakage in seals [8], which predicts that the hydraulic resistivities are only concentrated in few points where the fluid flow encounters strong restrictions. At these restrictions the pressure must present high gradients as clearly shown in Figures 5(b) and 6. In particular, in Figure 6 we report the average interfacial separation (black curve), the average solid contact pressure (green curve), the average fluid pressure (red curve), the total average normal stress (blue curve, slightly distinguished from the average solid contact pressure), and the solid contact pressure in the case of perfectly smooth contact (black dashed curve, useful for comparison), for (i.e., along the axial direction). Note that the largest pressure gradients occur in the very proximity of the largest values of , that is, where the local average separation takes its minimum value. It is also interesting to observe the normalized solid contact area for this static contact case, see Figure 7(a). Note that, due to the high local squeezing pressure, the asperities low-frequency features (i.e., the roughness asperities described by the deterministical model) have been squeezed so much to coalesce in larger contact patches. Figure 7(b) shows a sample of the leakage flux lines calculations (red lines) superposed to the locally averaged fluid velocity intensity (in gray scale). As expected [8], the presence of asperity-asperity contacts increases the hydraulic resistivity at the interface and, in particular, largest values of fluid velocity occur at the local minimum in the average interfacial separations.
| (a) Average solid contact pressure |
| (b) Average fluid pressure |
| (a) Normalized solid contact area |
(b) Sample of fluid flux lines
We have also carried out calculations assuming the seal obeys a Mooney-Rivlin model. The results are then compared with the elastic case. In Figure 8 we show the normalized area of real contact for the elastic (Figure 8(a)) and hyperelastic (Figure 8(b)) bulk rheology. Observe that for the linear case, the area of contact are slightly larger, as expected, then for the other case. This however, does not strongly affect the leakage flow for the current geometry, as shown in Figure 9. However, interestingly, the nonlinear rheology is characterized by a larger value of surface area change, evaluated as (see Figure 10), where we show in red colors that on part of the interaction domain there is a shrinking, that is, . Observe that in plane surface displacements are responsible for the frequency-shift of roughness PSD, for example, in the case of Figure 10(b) we expect a relevant modification of the homogenized power spectral density shape, which can therefore affect the local contact mechanics and flow factors. In the literature, surface displacements are not usually taken into account; however the present study shows that such a surface effect can be relevant for contact modeling.
(a) Linear elastic
(b) Mooney-Rivlin
(a) Linear elastic
(b) Mooney-Rivlin
(a) Linear elastic
(b) Mooney-Rivlin
3.2. Sealing at Non-Zero Sliding Velocity
In this section we show the case of non-zero sliding velocity. The macroscopic lip geometry is the same as before, as well as the roughness PSD. However, this time the part of roughness that is described deterministically is characterized by only one length scale, that is, just one frequency is included. The remaining roughness has been therefore included within the homogenized approach. In Figure 11 we show the cavitation areas (black areas), occurring over a circumferential portion of the lip-shaft macroscopic contact domain. Results are shown at different values of dimensionless sliding velocity. Interestingly, cavitation originates as expected in the low pressure side of the region (i.e., on the left side of the domain, whereas sliding velocity is directed from top to the bottom) and at the trailing edge of asperities. By increasing the sliding velocity, the cavitation extends from the low to the high pressure side. The cavitated areas may even coalesce at the largest sliding velocity, with the formation of cavitation fingers shown in Figure 11(e). The adoption of our two-scale approach allows then to capture the complex solid contact and fluid-dynamics characteristics of real contact geometry, which can help engineers to have useful insights into the mixed-EHL lubrication condition occurring at the interface of soft-contacts, especially in terms of friction and leakage, with much less computation effort. Figure 12 shows the axially component of fluid flow at the interface. Red areas corresponds to the counter-gradient flow (i.e., the flow directed in the opposite direction compared to the externally applied fluid pressure gradient). Observe that, on the high-pressure side, due to the asperities presence, a certain fluid recirculation is observed which tends to hamper the observed net leakage (blue fingers in the low-pressure side). The recirculation is due to the effect of fluid depressurization which occurs at the divergent part of the single asperity/substrate interfaces. This depressurization determines a relevant fluid suction from the high pressure side which is partially balanced by the flow induced by the counter-pressure gradient.
| (a) |
| (b) |
| (c) |
| (d) |
| (e) |
4. Conclusions
We have presented a novel two-scale approach for the description of the mixed lubrication regime for real soft-contacts. We modeled the asperity-asperity and asperity-fluid interactions with a deterministic or a statistical approach (DSA) depending on length scale at which the contact region is observed. The roughness at large length scales, which mainly determines the fluid flow at the interface, is deterministically included in the model while the roughness at short wavelengths, which strongly contributes only to the friction, is included by means of a homogenization process (recently developed in [9]). We have applied the DSA to lip sealings contact mechanics modeling, and we have analyzed the mixed lubrication characteristics at nearly-zero and non-zero shaft sliding velocity. In the case of nearly-zero sliding velocity the lip seal behaves as a classical static seal, and we showed that the fluid flow at the interface is determined only at the smallest constrictions along the leakage path, in agreement with recent developments of static seals theory. In the case of non-zero sliding velocity, we showed the occurrence of micro-cavitations and cavitation fingers, whereas leakage has been shown to be associated to asperity-induced fluid suction at the high pressure side (for the given geometry). Finally, we note that DSA-based mixed lubrication models, which belong to the class of multiscale contact mechanics models, provide a high-resolution description of very complex contact problems with a reduced or, at least, tunable computational effort, opening the perspective for its application in general purpose engineering software.
Appendix
Equation (10) has been discretized with the control volume approach. In particular, (10) can be integrated in a portion of the contact area considering that the elementary area gives where The forward difference for the Couette term has been used in order to get a stable scheme, in conjunction with the adoption of the Gauss-Seidel technique.
The average fluid shear stresses in the fluid zones can be calculated as: whereas in the cavitation zones:
The resolution scheme is summarized in Figure 13. The contact model is split into two, coupled problems, respectively, the deterministic and the homogenized problem. Given the macroscopic gap relation , the average surface stress is determined by solving the homogenized part, which consists of the Persson’s contact mechanics and of the homogenized fluid problem (previously discussed). The deterministic part follows, where the macroscopic deformation problem is solved as a function of the previously calculated average surface stress field. To do so, in this work we have used the Ansys finite element code; in particular, the lip geometry, similar to that described in [27], has been meshed with tetrahedral structural elements of type 92. The macroscopic gap relation is then finally updated, determining the loop restart. The solver iterates until certain convergence criteria are satisfied. In our case, convergence is checked on the average interfacial separation field . Under-relaxation, with relaxation factors in the range , is usually adopted to numerically damp the interfacial separation solution.
Nomenclature
| : | Fluid viscosity |
| : | Poisson’s ratio |
| : | Reduced root-mean-square roughness |
| : | Local normalized area of solid contact |
| : | Locally averaged interfacial separation |
| : | Parameter for the asymptotic interfacial separation law |
| : | Parameter for the asymptotic interfacial separation law |
| : | Threshold roughness wavelength |
| : | Full film density |
| : | Shear flow conductivity |
| : | Normal (locally averaged) stress |
| : | Pressure flow conductivity |
| : | Normal (locally averaged) fluid stress |
| : | Normal (locally averaged) solid contact stress |
| : | Tangential (locally averaged) fluid stress -component |
| : | Tangential (locally averaged) stress -component |
| : | Tangential (locally averaged) solid contact stress -component |
| : | Reduced sliding velocity |
| : | Roughness magnification or wavenumber |
| : | Threshold roughness wavenumber |
| : | Representative area of interaction |
| : | Area of solid contact in |
| : | Total roughness power spectral density |
| : | Power spectral density of the statistically-calculated roughness |
| : | Young’s modulus |
| : | Reduced elastic modulus |
| : | Cavitation index |
| : | Shear elastic modulus |
| : | Mapping rule Jacobian |
| : | Deterministic surface roughness height function |
| : | Inverse of the mapping rule Jacobian |
| Representative macroscopic size of the contact | |
| : | Shaft sliding velocity |
| : | Mean velocity |
| : | Macroscopic contact shape function |
| : | Locally stored elastic energy |
| : | Generic component of the average surface displacement vector |
| : | Initial state reference |
| : | Actual or deformed state reference. |
Acknowledgment
The authors acknowledge Regione Puglia for having supported part of this research activity through the constitution of the TRASFORMA Laboratory Network cod. 28.