Effect of Non-Linear Properties of Intermediate Layer on Sliding Contact of Homogeneous and Coated Elastic Solids

Abstract

Lubrication in the contact area can be considered as a viscoelastic layer, especially in the presence of particles in it, or under conditions of low temperatures. The properties of this layer are often non-linear, in particular, they depend on local pressure. The paper presents a formulation and numerical-analytical method for solving the contact problem in the presence of viscoelastic layers, the compliance of which depends on the applied pressure and is included in the formulation as the corresponding operator. The layer is homogeneous or coated elastic half-space. For the selected type of operator, the influence of parameters, which characterize the nonlinearity of the model, on the distribution of contact pressure and the coefficient of friction due to hysteresis losses was analyzed. It is shown that for the nonlinear model, the maximum contact pressures are higher, and the friction coefficient is lower than for the linear model with constant compliance. The effect of non-linearity for a wide range of sliding velocities is considered. An analysis of principal shear and tensile-compressive stresses for a homogeneous elastic half-space and for a coating, in particular, for a coating-substrate interface, was also carried out.

1. Introduction

A contact problem with thin viscoelastic layers on the surface arises when the contact occurs in the presence of a ‘third body’. This could be a fretting process [1] or lubricated contact [2,3] under conditions of high pressure in the contact area [4,5,6,7,8], as well as relatively low temperatures [9], in which oil lubricants begin to behave like solid viscoelastic bodies.

In cases where the thickness and/or modulus of elasticity of the layer is relatively large, 3D models of material are used [10,11,12,13], in which the layer is considered as an isotropic viscoelastic body. When the layer is relatively soft and thin, it is possible to use simplified one-dimensional viscoelastic models, including the Maxwell model [2,3,14,15], the Kelvin–Voigt model [16,17,18,19], and the standard viscoelastic body model [20,21]. More information about the viscoelastic contact can be found in review [22].

In references [14,15], the rolling of cylinders was considered in the presence of a lubricant. The effect of thin surface layers of lubricant on contact characteristics and internal stresses was analyzed. To describe the properties of surface layers in the normal and tangential directions, the one-dimensional model of Maxwell viscoelastic layer was used [14,15], and the lubricant flow between the surfaces was described using the Reynolds equation [15].

In reference [2], the mutual sliding of two discs was modeled in the presence of high-pressured lubricant. To simulate the viscoelastic properties of an oil film in a sliding contact, a one-dimensional nonlinear (in tangential direction) Maxwell model was proposed. This model is in good agreement with the Eyring fluid flow theory. In reference [3], a similar model was used, differing only in the form of the nonlinear function.

The rolling of an elastic sphere over an elastic half-space covered with a viscoelastic layer was studied in references [17,18]. The viscoelastic layer was modeled using the one-dimensional Kelvin model in tangential direction and a modified Winkler model for the normal contact. The distribution of stick and slip zones in the contact, as well as stresses in the contacting bodies, was calculated and analyzed [18].

In ref [19], the possibility of modeling an elastohydrodynamic contact using a nonlinear one-dimensional Kelvin−Voigt model was considered. The model was developed based on an elastic one-dimensional model obtained in references [23,24,25] from the Navier-Cauchy equation, assuming that the lubricant thickness is small. The model allows for consideration of the local absence of lubrication; therefore, several contact conditions were considered including mixed lubrication.

The rolling contact of a rigid sphere over an elastic half-space covered with a viscoelastic layer (simulating a lubricant layer) was studied in reference [21]. Two approaches were proposed: a one-dimensional model of a standard viscoelastic body and a complete model based on the Navier-Cauchy equations using Papkovich−Neiber potentials. It was shown that for effective modeling of elastohydrodynamic contact, it is necessary to use pressure-dependent viscoelastic models.

In this paper, we study the sliding contact of homogeneous and inhomogeneous elastic bodies in the presence of a viscoelastic layer in the contact area. A one-dimensional nonlinear Maxwell model is proposed, in which the layer compliance depends on the local pressure. The aim of the study was to develop an adequate method for solving such problems, as well as to analyze the influence of the nonlinearity factor on contact and internal stresses, and on the friction force arising due to sliding resistance.

2. Problem Formulation

Let us consider the motion of a smooth slider of arbitrary shape over the boundary of homogeneous or covered elastic half-space, on the surface of which there is a viscoelastic layer of thickness $h_3$ (Figure 1). A normal force, Q, acts on the slider. The layer simulates an intermediate medium, which may consist of lubricant wear (or abrasive) particles. In references [14,15], the mechanical properties of the layer are described using the viscoelastic Maxwell model:

$${\dot{w}}_3(x,y)=C\left(\frac{p(x,y)}{T}+\dot{p}(x,y)\right),$$

(1)

where $w_3$ are the normal displacements of the layer boundary, $C$ is the constant compliance of the layer, $T$ is the relaxation time, and $p\left(x,y\right)$ is the contact pressure.

This model has unlimited creep, which in some cases describes well the properties of the intermediate medium. At the same time, it is logical to assume that the compliance of the layer should depend on the pressure, for example, due to the fact that the layer becomes thinner and, therefore, stiffer in the area of high pressures. In this regard, the following modification of (1) is proposed:

$${\dot{w}}_3(x,y)=C[p(x,y)]\left(\frac{p(x,y)}{T}+\dot{p}(x,y)\right),$$

(2)

Here, the dependence of compliance on pressure is determined by the specific features of the intermediate layer.

The sliding contact occurs at a constant velocity, $V$, so the problem is quasi-static. To provide sliding at the constant velocity over the viscoelastic layer, an a priori unknown constant tangential force, $Q_T$, is required, equal to the sliding resistance force due to the rheological properties of the layer. Let us position the Cartesian coordinate system, associated with the slider, in such a way that the XY plane coincides with the undeformed upper surface of the viscoelastic layer and the Z axis is normal to the layer surface. In this case, the contact condition will look like this:

$$\begin{array}{l}{w}_1\left(x,y\right)+w_3\left(x,y\right)=f\left(x,y\right)-D \hspace{1em} (x,y)\in \Omega \\ {\sigma }_z=0\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} \hspace{1em} \hspace{1em}\hspace{1em}\hspace{1em} (x,y)\notin \Omega \\ {\tau }_{xz}=0, {\tau }_{yx}=0 \end{array}$$

(3)

Here, $w_1\left(x,y\right)$ are the vertical displacements of the upper surface of a homogeneous or covered elastic half-space, $f\left(x,y\right)$ is a function that describes the shape of the slider, D is the penetration of the slider into the half-space, and $\Omega$ is the contact area. The material properties of a homogeneous half-space are determined by Young’s modulus $E_2$ and Poisson’s ratio ${\nu }_2$. In the case of coated half-space, we have $E_i$ and ${\nu }_i$, where $i=1, 2$ corresponds to the coating and the half-space, respectively. In the “spring” model under consideration, the normal stresses at the upper boundary of the layer ($z=0$) and at the boundary between the layer and the elastic half-space ($z=h_3$) have the same values, since the stresses do not change along the length of the spring.

It is assumed that the surface layer has the function of lubrication, so the friction force is not taken into account in the formulation of the contact problem (3).

An equilibrium condition is also satisfied:

$$Q={\displaystyle \underset{\Omega }{\iint }p\left(x,y\right)dxdy},$$

(4)

Since we are considering a smooth slider, the zero-pressure condition must be satisfied at the boundary of the contact area $\Omega$.

In the presence of a coating of thickness $h_1$, it is also necessary to formulate the conditions at the interface between the coating and the half-space:

$$\begin{array}{l} w^{\left(1\right)}=w^{\left(2\right)}, u_x^{\left(1\right)}=u_x^{\left(2\right)}, u_y^{\left(1\right)}=u_y^{\left(2\right)},\\ {\sigma }_z^{\left(1\right)}={\sigma }_z^{\left(2\right)}, {\tau }_{xz}^{\left(1\right)}={\tau }_{xz}^{\left(2\right)}, {\tau }_{yz}^{\left(1\right)}={\tau }_{yz}^{\left(2\right)}.\end{array}$$

(5)

Here, ${\sigma }_x^{\left(i\right)}$, ${\tau }_{xz}^{\left(i\right)}$, ${\tau }_{yz}^{\left(i\right)}$ are normal and tangential stresses; $w^{\left(i\right)}$, $u_x^{\left(i\right)}$, $u_y^{\left(i\right)}$ are normal and tangential displacements. Relations (5) correspond to the conditions of full adhesion at the coating-substrate interface, when the continuity of normal and tangential stresses and displacements takes place.

It is required that the contact area, the distribution of contact pressure and the slider penetration are found.

3. Method of Solution

The solution was obtained using the boundary element method and an iterative procedure. A rectangular area ${\Omega }^{\ast }\supset \Omega$ is selected, on which a mesh of n identical square elements of size $2a\times 2a$ is constructed. The shape of the slider is represented as a set of values $f_i\hspace{0.17em}(i=1\dots n)$ which are $f(x,y)$ for the center of each square element; the resulting pressure is obtained as a piecewise function $p_i\hspace{0.17em}(i=1\dots n)$. Knowing the dependences of the vertical displacements of the layer boundary on the uniformly distributed pressure inside a square element of the surface, we can transform Equations (3) and (4) into a system of $n+1$ linear algebraic equations:

$$\left[\begin{array}{cc}\left(\begin{array}{ccc}{A}_{11}& \cdots & A_{1n}\\ ⋮& \ddots & ⋮\\ A_{n1}& \cdots & A_{nn}\end{array}\right)& \begin{array}{c}1\\ ⋮\\ 1\end{array}\\ \left(\begin{array}{ccc}4a^2& \cdots & 4a^2\end{array}\right)& 0\end{array}\right]\cdot \left[\begin{array}{c}{p}_1\\ ⋮\\ p_n\\ D\end{array}\right]=\left[\begin{array}{c}{f}_1\\ ⋮\\ f_n\\ Q\end{array}\right],$$

(6)

The last equation here is the equilibrium condition (4) in terms of boundary elements: $A$ is a matrix of influence coefficients, which are boundary displacements related to constant pressure inside each square element. They are obtained by superposition of the displacements of a viscoelastic layer and a homogeneous (covered) half-space:

$$A_{ij}=A_{ij}^{(3)}+A_{ij}^{(1)},$$

(7)

where $A_{ij}^{\left(3\right)}$ corresponds to the normal displacements of the viscoelastic layer, and $A_{ij}^{\left(1\right)}$ corresponds to a homogeneous or covered half-space. Using direct integration of (2), and taking into account that non-zero displacements are only on the part of the layer surface that was or is under the influence of distributed pressure, we obtain the following expressions for $A_{ij}^{\left(3\right)}$:

$$\begin{array}{l}{A}_{ij}^{(3)}=0, \left(x_j-x_i\right)>a,-a>\left(y_j-y_i\right)>a\\ A_{ij}^{(3)}=-C_i\left(1+\frac{a-\left(x_j-x_i\right)}{T_{}V}\right),-a<\left(x_j-x_i\right)<a,-a<\left(y_j-y_i\right)<a\\ A_{ij}^{(3)}=-C_i\left(\frac{2a}{T_{}V}\right), \left(x_j-x_i\right)<-a,-a<\left(y_j-y_i\right)<a\end{array}$$

(8)

Here, $C_i$ are determined by the form of the operator C for a point with $(x,y)$ coordinates.

The matrix coefficients $A_{ij}^{\left(1\right)}$ are calculated using the relations for displacements of the elastic half-space boundary [26]:

$$A_{ij}^{(1)}=-\frac{\left(1-{\nu }_1{}^2\right)}{\pi E_1}{\displaystyle \underset{{\Omega }_i}{\iint }\frac{1}{\sqrt{{\left(x_j-\xi \right)}^2+{\left(y_j-\eta \right)}^2}}}d\xi d\eta ,$$

(9)

For a two-layer elastic half-space, the matrix coefficients $A_{ij}^{\left(1\right)}$ are obtained using a method based on double integral Fourier transforms [27]:

$$\begin{array}{l}{A}_{ij}^{(1)}=-\frac{a(1+{\nu }_1)}{\pi E_1}{\displaystyle \underset{0}{\overset{\pi /2}{\int }}{\displaystyle \underset{0}{\overset{\infty }{\int }}\Delta \left(\gamma ,\phi ,\lambda ,E_1,E_2\right)}}\cos \left(\left(x_i-x_j\right)\gamma \cos \left(\phi \right)/a\right)\\ \cos \left(\left(y_i-y_j\right)\gamma \sin \left(\phi \right)/a\right)d\gamma d\phi \end{array}$$

(10)

here $\gamma$, $\phi$ are coordinates in the space of double Fourier transforms; $\lambda =h_1/a$, $\Delta \left(\gamma ,\phi ,\lambda ,E_1,E_2\right)$ is an analytic function, which is obtained from boundary conditions (formulated for uniformly loaded surface elements) using the Fourier transforms.

Using an iterative procedure, at each step of which system (4) is solved for $A_{ij}=A_{ij}^{\left(1\right)}$ and the contact conditions (including zero pressure at the boundary of the contact zone) are satisfied, the unknown contact area, pressure distribution, and penetration D are obtained for an elastic body without the viscoelastic layer. Next, an iterative procedure is carried out in order to obtain a solution based on the following type of $C\left[p\left(x,y\right)\right]$:

$$C\left[p(x,y)\right]=C_0\left(1+\delta -\frac{2\delta }{1+e^{-p/p_0}}\right),$$

(11)

This type of operator is chosen from the following considerations. In the absence of pressure ($p=0$), the compliance of the layer is $C_0$. The action of pressure leads to a decrease in the thickness of the layer in the contact zone; it is assumed that the layer thickness at infinite pressure tends towards a finite value (including zero). Thus, the minimum compliance of the layer is defined as $C_0\left(1-\delta \right)$, where $\delta$ is a dimensionless parameter. The parameter $p_0$ has the dimension in Pascal and determines the sensitivity of the layer properties to changes in the contact pressure. Figure 2 illustrates the model for different values of $p_0$ and $\delta$.

The calculated contact pressure was used to find internal stresses in the coated elastic half-space. The method for calculating stresses in a two-layered elastic half-space is based on double integral Fourier transforms [10]. The expressions for such calculations were obtained for layered elastic half-space in reference [28] and used, for example, in reference [10].

The method was used to calculate tensile-compressive (${\sigma }_x$) stresses in the coated half-space.

The method for calculating stresses in a two-layered elastic half-space is based on double integral Fourier transforms. In the case of homogeneous half-space, the Love solution [29] was used for the constant pressure distributed inside a surface element.

The resulting stresses (both for the coated and the homogeneous half-space) are obtained by superposition.

4. Results and Discussion

The following dimensionless parameters are used to analyze the results for a spherical slider with radius R:

$$\begin{array}{l}\left(x^{\prime },y^{\prime },z^{\prime }\right)=\left(x,y,z\right)/R, C^{\prime }=C_0\cdot E_2/R, Q^{\prime }=Q/\left(R^2{E}_2\right),\\ p^{\prime }\left(x^{\prime },y^{\prime }\right)=p\left(x^{\prime },y^{\prime }\right)/E_2, V^{\prime }=VT/R, E_1^{\prime }=E_1/E_2, {p_0}^{\prime }=p_0/E_2,\\ {{\sigma }^{\prime }}_{x,y,z}={\sigma }_{x,y,z}/E_2, {\tau }_1^{\prime }={\tau }_1/E_2.\end{array}$$

Let us present the results obtained for an elastic half-space with a layer for which the properties are described by relation (11). Here, and further, the sliding direction is positive (along the $0x^{\prime }$ axis). The sliding velocity was chosen for reasons of visual demonstration of the influence of viscoelastic properties on the characteristics of the contact. Figure 3 shows the distribution of the contact pressure and the compliance of the layer in the contact area. The contact pressure distribution is asymmetric, which is typical for the sliding contact of materials with rheological properties. The contact spot has a shape close to a semicircle. The compliance of the layer near the center of the contact area is close to minimal value.

The parameters $p_0$, $\delta$ describing the effect of local pressure on the properties of the layer in (11) also significantly affect the contact pressure distribution. The influence of $p_0$ is illustrated by the curves in Figure 4. Solid curves of contact pressure distribution lie between two dashed curves relating to linear models with the maximum and minimum compliance allowed by the non-linear model (11). The smaller the value of $p_0$, the greater the maximum contact pressure and the smaller the size of the contact area becomes.

A similar analysis was made for the parameter $\delta$ (Figure 5). The fact that the contact pressure increases with the increase in $\delta$ is predictable. Interestingly, the blue curve, calculated for the complete absence of the layer in the limiting case of infinite pressure, is more symmetric, i.e., the influence of the rheological properties of the layer decreases with increasing $\delta$.

The sliding resistance is directly related to the asymmetry of the contact pressure. The moment of forces and the coefficient associated of which are defined as:

$$M={\displaystyle \underset{\Omega }{\iint }x\cdot p\left(x,y\right)dxdy},\hspace{1em}\mu =\frac{{\displaystyle \underset{\Omega }{\iint }x\cdot p\left(x,y\right)dxdy}}{RQ},\hspace{1em}{Q}_T=-\mu Q,$$

(12)

Since the layer has lubricating properties, the sliding resistance due to viscosity can make a significant contribution to the total frictional force. Below is an analysis of the influence of the input parameters of the model on the value of $\mu$.

Figure 6 shows curves illustrating the influence of the relative compliance of the layer on the value. In this case, a decrease in the relative compliance of the layer is provided by increasing the rigidity of the half-space at a fixed load. The greater the rigidity of the half-space, the more significant the effect of the layer on the contact conditions. This explains the increase in the value of $\mu$, in which the asymptote probably tends to the value typical for the case of contact between the slider and the layer on a rigid base. It should be noted that the curves obtained for the linear and non-linear models are similar in shape and slope at each point, especially for relatively large values of $C^{\prime }$. The same can be said about the curves illustrating the effect of the thickness of relatively rigid coatings on $\mu$ (Figure 7), since an increase in thickness leads to an increase in the integral rigidity of a two-layered elastic body.

It is interesting to analyze the effect of the non-linearity factor on the friction coefficient as a function of sliding velocity (Figure 8). For a 3D model of the viscoelastic layer [10], as well as for the standard viscoelastic body model [20], this dependence usually is non-monotonic. For the Maxwell model the friction coefficient decreases with the velocity increase. Generally, for the non-linear model the friction coefficient is smaller than for linear one. The effect depends on the model parameters and the sliding velocity. The difference between $\mu \left(V\right)$ dependences for the linear and non-linear models of the layer (Figure 8b) demonstrates maximum at a fixed value of velocity. The value depends on the model parameters. The coefficient of friction here is significantly higher than for the results from the previous figures, since in this case a larger load value, Q′, was used for calculations.

Further analysis concerns the stresses that arise in an elastic half-space during sliding. The tensile and principal shear stresses used in most strength criteria were considered. Figure 9 shows the distributions of tensile-compressive stresses on the surface of an elastic half-space in the presence of a viscoelastic layer characterized by different values of $\delta$. For a homogeneous half-space, this component of the stress tensor has maximum value on the surface. The compressive stresses in this case are negative, and the tensile stresses are positive. The tension maximum takes place at the boundary of the contact spot. As $\delta$ increases, the maximum values of both tensile and compressive stresses increase. It should be noted that, in contrast to the case of a viscoelastic half-space [30], the presence of a layer does not lead to an instantaneous change in the stress sign at the edge of the contact area, although the gradient of the curves near the frontal part of the contact area is quite large.

The principal shear stresses are calculated for a plane parallel to the sliding direction and passing through the point of initial contact (Figure 10). The shape of isolines is typical for frictional contact. An increase in the value of $\delta$, which characterizes the degree of nonlinearity of the layer properties, from 0.1 to 1.0 leads to an increase in the maximum stress values of almost 35 percent.

In the presence of coatings, it is important to evaluate the stresses not only at the surface but also at the interface, which often has poor strength properties. Distributions of tensile-compressive stresses were obtained for a plane parallel to the sliding direction and passing through the point of initial contact (Figure 11) for the coatings with different thickness (using parameters identical to those in the non-linear model). The non-linearity of the layer properties is manifested in an increase in the maximum contact pressure: for a thinner coating by 12, and for a thicker one by 16 percent. The left picture illustrates the case of relatively thin coating, in which the thickness is almost 40 times less than characteristic size of contact zone in $0x^{\prime }$ direction. The stresses inside the coating are almost the same as at the surface. For the case of the thicker coating (Figure 11b), the coating bend with compression at the surface and tension at the coating-substrate interface occurs. Figure 12 shows the distributions of tensile-compressive stresses for coatings of different thicknesses on the surface and at the interface. Since there is a jump of stresses at the interface, two curves are obtained for it (for coating and for substrate materials). The dotted lines are the distributions obtained for the linear model of the viscoelastic layer, and the solid lines are for the nonlinear model. The non-linearity factor significantly affects the distribution of surface stresses for the thicker coating but has almost no effect on the stresses at the interface. In the case of thin coating, the influence of this factor is approximately the same for the surface and the interface in percentage terms. For the relatively thick coating, maximal tension occurs at the coating-substrate interface under the contact spot. It means that brittle cracks can be initiated at the interface. For the case of relatively thin coatings, the maximal tension is at the surface. It should also be noted that for thin coatings there is a sharp change in tension-compression, which can be dangerous in terms of the formation of vertical cracks coming from the surface.

5. Conclusions

The developed model makes it possible to estimate the influence of the nonlinearity factor, which manifests itself in the dependence of the compliance of the viscoelastic Maxwell layer on pressure on the quasi-static sliding of a smooth indenter over the viscoelastic layer lying on a homogeneous or two-layered elastic body. A numerical-analytical method for solving the problem has been developed, which can be used for various types of operators relating compliance to local contact pressure. The results obtained for the selected operator, corresponding to the physical decrease in the layer thickness under the action of pressure, were analyzed. It is shown that in this case the nonlinearity factor leads to an increase in the maximum values of the contact pressure and some decrease in the size of the contact area. The study of the influence of the input parameters of the model on the coefficient of friction due to the resistance to sliding by the viscoelastic material showed that for the non-linear model this coefficient is smaller, with all other parameters being equal. The difference between the coefficient for linear and non-linear models depends on sliding velocity non-monotonically. The more rigid (relative to the viscoelastic layer) a homogeneous or two-layer half-space is, the greater the value of the coefficient. An analysis of the stresses in the elastic half-space under the layer showed an increase in their maximum values when the nonlinearity of the viscoelastic layer was taken into account. For bodies with rigid coatings, it is shown that the nonlinearity factor significantly affects surface stresses, and for relatively thin coatings, it also affects stresses at the interface.