Normal contact between a rigid surface and a viscous body: verification of the method of reduction of dimensionality for viscous media Текст научной статьи по специальности «Математика»

y^K 539.62

Normal contact between a rigid surface and a viscous body: Verification of the method of reduction of dimensionality

for viscous media

S. Kurschner, A.E. Filippov1

Berlin University of Technology, Berlin, 1062, Germany 1 Donetsk Institute for Physics and Engineering of NASU, Donetsk, 83114, Ukraine

Analogies between elastic and viscous contact problems are used to apply boundary element algorithms and the method of reduction of dimensionality respectively to the problem of normal contact between a rigid body and a viscous half space. Some basic examples are examined in order to compare both methods and their results to each other.

Keywords: normal contact, viscous body, boundary element method, reduction of dimensionality

1. Introduction

Macroscopic tribological questions usually are multiscale problems. Therefore, direct numerical simulations of such processes are very time-consuming. Decreasing the number of relevant spatial dimensions by suitable transformations can reduce computation time significantly. Such kinds of transformations have been proposed by Popov and Psakhie [1]. Geike and Popov [2, 3] applied this method to the elastic contact problem. In [4], it is stated that the method of reduction of dimensionality can be extended straightforward to contacts of media with arbitrary linear rheology. However, a formal proof of this statement has not been provided till now. In the present paper we use some analogies between elastic and viscous problems published by Lee [5] and Radok [6] to prove the method of reduction of dimensionality in the case of viscous contact problem. Both analytic and numerical proofs are provided. The paper is organised as follows. Radok’s method is outlined and applied to the current problem in Sect. 2. In Sect. 3, these results are used to formulate a boundary element method for the viscous contact problem. In Sect. 4 the method of reduction of dimensionality is adapted to viscous contact problems by use of Radok’s method. In the following sections the boundary element method and the method of reduction of dimensionality are applied to some examples. In Sect. 5 a half space is indented by a smooth rigid rod. The same is done for a paraboloid of revolution in Sect. 6. In both cases the

validity of the method of reduction of dimensionality as formulated in [4] has been confirmed.

2. Radok’s method of functional equations

In this section we will reduce the viscous problem of normal contact to an equivalent elastic problem. In his paper on stress analysis in visco-elastic bodies Lee [5] published this concept first. It was extended and called the method of functional equations by Radok [6].

The stress-strain relation of an isotropic, linear elastic body can, according to [7], be written as

3 °ll 8ik + sik = K ell 8ik + 2Geik ■ (1)

As usual the strain tensor e is decomposed into pure compression ell8ik and pure shear eik ■ The stress tensor is split the same way. We use the Kronecker delta 8ik and the Einstein notation, K and G are the bulk and shear moduli and

eik = eik 3ell8ik ’ Sik ~ ^ik 3 8ik (2)

are the strain and stress deviators, respectively. A separated comparison shows

= 3Keu > (3)

Sj= 2Ger (4)

We now analyse an isotropic, linear viscous fluid. According to [8] the stress-strain relation reads

© Kurschner S., Filippov A.E., 2012

3 8ik + sik = K eii8 ik+Ze ii 8ik+2neik ■

(5)

The stress and strain tensors are decomposed as in the elastic case. A dot marks differentiation with respect to time, K is the fluid’s bulk modulus, Z and n are bulk and shear viscosity, respectively. Analogous to the elastic problem we get the relations

Vii = 3K ea + Ze ii, (6)

sij = 2n4 ■ (7)

We now assume incompressibility in both cases and hence ou = 0 follows from Eqs. (3) and (6) immediately.

According to Radok’s method we expand the right hand sides of (4) and (7) into series of linear time differential operators with constant coefficients

q dn „ q dn

Q = E qn tv > Q = 2 qn — ■

n=0 dtn n=0 dt

(8)

In the case of the elastic material the series consist of a term of order zero only

Sj = Qej, Q = 2G. (9)

We denote the viscous problem by a tilde. A single first order operator occurs in this case d

sy = Qej, Q = 2n^ •

(10)

Now we study the displacement of an incompressible, viscous fluid under the influence of normal force F acting on its surface. The fluid is treated as a viscous half space and the force is assumed to act on a small part of the surface. At this time we are not interested in tangential displacements or subsurface effects but in the normal direction of surface deflections only.

The solution of the comparable elastic problem is well known [7]:

F

u =-------, (11)

4nGr

where r denotes the distance between the point of interest and the point the force is acting on. Application of Radok’s method to the solution of the elastic problem provides a solution of the viscous problem. The shear modulus is substituted according to (9) and the Laplace transform is applied to the elastic solution. Therefore some additional assumptions have to be made to the problems history. For simplicity we consider F(t) = 0 if t < 0 and F(t) = F0 for all t > 0. Within Laplace domain the operator series Q describing the elastic material is replaced by the viscous one. After that the equation is taken back to time domain by using the inverse Laplace transform. The solution of the viscous problem reads F

(12)

where u denotes the velocity field of the surface deflection. The total surface deflection field u is provided by integration with respect to time.

However, equation (11) is a Green’s function describing the normal deflection of an elastic half space at the surface under a normal point force. In the same manner one can identify Eq. (12) to be a Green’s function describing the deflection velocity of a viscous half space. It is easily seen that Eqs. (11) and (12) agree with common models.

3. Boundary element algorithm

We now study the problem of dipping a rigid body into a viscous fluid. As before, we assume the fluid to meet the half space assumption. The process is driven by a constant force F, pressing the rigid into the viscous body. Inertia forces are neglected.

We study this problem by means of the boundary element method. As described in the previous section we can take use of the comparable elastic problem and methods to solve it. Therefore we firstly examine the elastic normal contact closely following Geike’s dissertation [9].

One major challenge in solving this problem is that almost nothing is known a priori, neither deformations nor stresses or even position or size of the contact area. Thus, one has to guess and iterate these quantities. A common algorithm is drafted in the left hand side of Fig. 1. The calculations of displacement fields (outside contact area) and stress distributions are based on the Green’s function (11) directly. In Fig. 1 they are drawn as rounded-cornered boxes.

In the previous section it was shown that Eq. (11) is transformed into the Green’s function of the viscous problem by replacing shear modulus by shear viscosity and displacement by deflection rate as a result of Radok’s method. Hence the described algorithm can be used to calculate stress distributions and velocity fields for viscous contact problems

Start configuration

Start configuration

Displacement field (inside contact area only)

I I

| Increase of contact area |

[ Stress distribution |

| Decrease of contact area |

[ Displacement field Iteration

Velocity field (inside contact area only)

Stress distribution I _

Decrease of contact area

Iteration

Velocity field

ix:

Time step length + ~

Displacement field I -----------------

Increase of contact area

Time step

Fig. 1. Boundary element algorithms to solve elastic (a) or viscous (b) normal contact problems

too. Only some minor changes have to be made, most of them connected to geometric conditions.

The displacement field is provided by integration with respect to time. This is the main drawback in this topic. In the course of time, the indentation depth changes and in general the contact area changes too. Hence one has to solve the contact problem in every single time step or at least if the change in contact area is larger than any chosen threshold. This algorithm is shown in the right hand side of

Fig. 1.

4. Method of reduction of dimensionality

A contact of two bodies is generally a three-dimensional problem. Geike and Popov [2, 3] published an idea how to map such problems onto one-dimensional ones. In [10] it was applied to a study of frictional forces between a rigid rough surface and a simple linear viscous body. In [11] this approach has been extended to an arbitrary linear rheology. Such simulations became only possible due to reduction of the dimensionality of the corresponding problems, as the method of reduction of dimensionality decreases computation time drastically. Geike and Popov [2] have shown in the case of normal contacts that essential contact properties are preserved under this mapping. We will now apply Radok’s method to the reduction method.

We again start with the elastic problem. One can imagine the one-dimensional elastic model as a row of spring elements independent of one another (elastic foundation). As before, we restrict ourselves to an incompressible material. According to [2], the elements stiffness per unit length has to be chosen as

c = 4G.

(13)

The radii of curvature R3D have to be reduced by half under the mapping

(14)

Remembering the results of Sect. 2, it is clear how the viscous one-dimensional model has to be constructed. The spring elements have to be replaced by linear dashpots with viscosity

d = 4n (15)

per unit length. Equation (14) remains unchanged. This viscous one-dimensional model is drafted in Fig. 2. In the next sections we will examine these models with some simple examples.

5. Example: Indentation by a smooth rod

Consider a smooth cylindrical rod (radius a) which is pressed into an incompressible, elastic body (shear modulus G) by a constant force F. This case is examined in detail in [4]. The radial distribution of normal stress accords to

P(r) = Po

, r 2 TV2

1 - a2

V

r < a■

The central stress p0 is governed by

P0 = W

and the indentation depth reads F

8aG

(16)

(17)

(18)

Two important contact properties which are closely connected to each other are contact length l and contact stiffness k. The latter one can be defined by the quotient of acting force and resulting indentation depth

(19)

The contact length can be interpreted as the sum of diameters of all contact regions. In this case, it is the diameter of the rod. In a more general manner the contact length can be calculated by

l = — ■ (20) 4Gu

Both properties have been studied recently in relation to contact problems of elastic bodies with fractal surfaces [12].

Leaving this preliminary considerations behind us, we are now interested in the comparable viscous contact problem. Thus we press the same rod into an incompressible viscous body (shear viscosity n). By the results of Sect. 2 one can see: stress distribution (16) and central stress (17) remain unchanged but (18) is mapped onto

u =F ■ (21)

8an

Fig. 2. Dipping a rigid body into a viscous half space (a) and the one-dimensional model (b)

Fig. 3. Indentation of a viscous half space by a smooth cylindrical rod, through the centre of contact area

r/a

displacement rate (a) and normal stress distribution (b) along an axis

But this remains a more general result. Obviously the viscous problem does not depend on the displacement itself. Hence dipping of any rigid body with a circular cross-section area (as cones, spheres or paraboloids of revolution) is governed by Eqs. (16), (17) and (21). Although the size of contact area may differ in time, it can be treated similarly to a rod at any time. We will return to this subject in the next section.

By the same mapping contact stiffness (19) and contact length (20) are transformed into F

k = — = 8 an, (22)

u

Consider a rod with radius a = 1 m is pressed by a normal force F = 1 N into a viscous half space with shear viscosity n = 1 Pa • s. The central stress and the indentation rate can be calculated analytically by Eqs. (17) and (21):

pA = 0.1592Pa, uA = 0.125ms. (24)

For comparison the face plane of the rod was discretised into 31417 quadratic elements and the problem was solved by means of the boundary element method described in Sect. 3. The following values were obtained:

p0B = 0.1608 Pa, uB = 0.1252 ms. (25)

As one can see, the analytical (24) and boundary element result (25) match rather well. Figure 3(a) shows the displacement rate of the free surface along an axis through the centre of the contact area. Inside the contact area (| r| < a) this velocity is equal to uA in case of the analytical calculation (marked by dots) respectively to u B by means of boundary element calculation (solid line). In addition the latter provides the velocities decrease outside the contact area.

Figure 3(b) compares the stress distributions by means of analytical (dots) and boundary element (solid line) cal-

culations. These do fit as well, even near the singularities at the boundary of the contact area.

We now treat the problem as a one-dimensional one by means of the method of reduction of dimensionality. The rod is modelled as a constant line (length — 2 m) and pressed into a viscous foundation according to Sect. 4. In this special case, the number of elements is of minor importance. One obtains, due to modelling, the exact result for the displacement velocity

uR = 0.125 ms. (26)

But the model does not provide any information about the stress distribution or surface displacement outside the contact area.

6. Example: Indentation by paraboloids of revolution

Now we indent a viscous body by a rigid paraboloid of revolution. Thus the indenters shape can be represented by

h(x, y) = c(x2 + y2), (27)

where a parameter c describes its level of curvature. As before, we are interested in the normal contact problem in connection with an incompressible half space only. We restrict ourselves to the case that the indenter axis of symmetry is perpendicular to the free surface. Thus, the contact area is circular. In addition the whole contact area is displaced with the same velocity, due to the indenter’s rigidity. However, the size of the contact area changes gradually. But according to the last section this process is driven by Eqs. (16), (17) and (21) at any point in time.

We use this example to compare the results provided by means of boundary element calculations and by means of the method of reduction of dimensionality. We choose indentation depth and contact length as functions of time for comparison of these two calculation methods.

The boundary element model consists of a square, rigid indenter with a symmetric surface described by Eq. (27). It

Fig. 4. Paraboloids of revolution with different levels of curvature indent a viscous half space. Indentation depth (a) and contact length (b) are shown, calculated by means of the boundary element method (symbols) and the method of reduction of dimensionality (solid lines)

is discretised with 257x257 quadratic elements. The total edge length is 2.56 m in x- and j-direction each. Five inden-ters with different levels of curvature (c = 0.05, 0.2, 1, 5, 20) are used. Each of them is pressed into a half space of viscosity n = 1 Pa • s due to a constant normal force F = = 1 N. Each time five seconds of indentation were simulated, starting from the first touch. One indentation took up to 1.5 hours of computation time. The algorithm includes multilevel multisummation techniques developed by Venner and Lubrecht [13] and the iteration itself is done by the conjugate gradient method by Polonsky and Keer [14]. In addition stress distribution and velocity field were recalculated only, if the change in contact area became larger than 1 % compared to the contact area. This speeds up the calculation enormously but may cause a minor loss of accuracy. The indentation depth is provided by numerical integration with respect to time as described in Sect. 3. The contact length is calculated according to Eq. (23).

In accordance to Sect. 4, the one-dimensional model consists of a parabolic indenter with a one-dimensional height distribution represented by

h( x) = 2cx2, (28)

where c is the same parameter as before. Its total length is 2/Vn-2.56 m according to [4] and it is discretised with 66 049 elements too. The other parameters are left unchanged in comparison to the boundary element method. One indentation lasts around a second without use of any speedup techniques. This is a huge advantage over the boundary element method.

Indentation depth and contact length are shown in Fig. 4. Solid lines are related to the one-dimensional model. Single time steps by means of the boundary element method are marked by symbols as described in the legends of the diagrams. The parameter c refers to the level of curvature used in Eqs. (27) and (28). The results of both methods match

nearly perfectly. Just when the contact length comes close to saturation one can see slight differences.

7. Conclusions

We used some analogies between elastic and viscous contact problems to adapt the boundary element method and the method of reduction of dimensionality to the latter one. Both methods were applied to simple contact problems. The independently obtained results match each to another rather well. Hence, the method of reduction of dimensionality is successfully applicable to viscous contact problems, at least in simple cases. This means a huge reduction of computation time. However, some further studies are necessary.

The authors are grateful to V.L. Popov for the idea of this work and many valuable discussions and suggestions to the manuscript. We are grateful to R. Pohrt for providing algorithms and useful discussions. A.E. Filippov acknowledges financial support of the Deutsche Forschungsgemein-schaft.

References

1. Popov V.L., Psakhie S.G. Numerical simulation methods in tribology // Tribol. Int. - 2007. - V. 40. - No. 6. - P. 916-923.

2. Geike T., Popov V.L. Mapping of three-dimensional contact problems into one dimension // Phys. Rev. E. - 2007. - V. 76. - No. 3. - P. 036710 (5 pp.).

3. Geike T., Popov V.L. Reduction of three-dimensional contact problems to one-dimensional ones // Tribol. Int. - 2007. - V. 40. - No. 6. -P. 924-929.

4. Popov V.L. Contact Mechanics and Friction, Physical Principles and Applications. - Berlin: Springer-Verlag, 2010. - 362 p.

5. Lee E.H. Stress analysis in visco-elastic bodies // Quart. Appl. Math. -

1955. - V. 13. - No. 2. - P. 183-190.

6. Radok J.R.M. Visco-elastic stress analysis // Quart. Appl. Math. -1957. - V. 15. - No. 2. - P. 198-202.

7. Landau L.D., Lifschitz E.M. Lehrbuch der Theoretischen Physik, Band 7: Elastizitatstheorie. - Berlin: Akademie-Verlag, 1965. - 183 p.

8. Landau L.D., Lifschitz E.M. Lehrbuch der Theoretischen Physik, Band 6: Hydrodynamik. - Berlin: Akademie-Verlag, 1991. - 683 p.

9. Geike T Theoretische Grundlagen eines schnellen Berechnungsver-fahrens fur den Kontakt rauer Oberflachen. Dissertation. - Berlin: TU Berlin, 2008.

10. Popov VL., Filippov A.E. Force of friction between fractal rough surface and elastomer // Tech. Phys. Lett. - 2010. - V. 36. - No. 6. -P. 525-527.

11. Popov VL., Dimaki A.V Using hierarchical memory to calculate friction force between fractal rough solid surface and elastomer with

arbitrary linear rheological properties // Tech. Phys. Lett. - 2011. -V. 37. - No. 1. - P. 8-11.

12. Pohrt R., Popov VL. Normal contact stiffness of elastic solids with fractal rough surfaces // Phys. Rev. Lett. - 2012. - V. 108. - P. 104301 (4 pp.).

13. Venner C.H., Lubrecht A.A. Multilevel Methods in Lubrication. -Amsterdam: Elsevier, 2000. - 379 p.

14. Polonsky I.A., Keer L.M. A numerical method for solving rough contact problems based on the multi-level multi-summation and conjugate gradient techniques // Wear. - 1999. - V. 231. - No. 2. - P. 206-219.

nocTynnna b pe^aKunro 10.06.2012 r.

CeedeHun 06 aem0pax

Kurschner Silvio, Dipl.-Ing., Berlin University of Technology, Germany, [email protected]

Filippov Alexander E., Prof., Donetsk Institute for Physics and Engineering of NASU, Ukraine, [email protected]