A frictional contact problem with wear diffusion

This paper constructs and analyzes a model for the dynamic frictional contact between a viscoelastic body and a moving foundation. The contact involves wear of the contacting surface and the diffusion of the wear debris. The relationships between the stresses and displacements on the contact boundary are modeled by the normal compliance law and a version of the Coulomb law of dry friction. The rate of wear of the contact surface is described by the differential form of the Archard law. The effects of the diffusion of the wear particles that cannot leave the contact surface on the surface are taken into account. The novelty of this work is that the contact surface is a manifold and, consequently, the diffusion of the debris takes place on a curved surface. The interest in the model is related to the wear of mechanical joints and orthopedic biomechanics where the wear debris are trapped, they diffuse and often cause the degradation of the properties of joint prosthesis and various implants. The model is in the form of a differential inclusion for the mechanical contact and the diffusion equation for the wear debris on the contacting surface. The existence of a weak solution is proved by using a truncation argument and the Kakutani--Ky Fan--Glicksberg fixed point theorem.


Introduction
This work studies a nonlinear dynamical model for the process of contact between a viscoelastic body and a reactive foundation when wear debris is generated and diffuses on the contact surface. The model includes subdifferential friction boundary condition, and considerably extends the model and the results in [20], which were announced in [19] and further developed in [8,10]. Additional information and details can be found in [21]. The research in [20] was motivated, in part, by biomechanical applications. Indeed, such problems arise in artificial joints after arthroplasty (knee, hip, shoulder, elbow, etc.) where debris is produced by articulating parts of the prosthesis and is transported to the bone-implant interface. The debris causes the deterioration of the interface, and is believed to be an important factor leading to prosthesis loosening (see, e.g., [17,18] and references therein). Thus, there is a considerable interest in modeling such complex contact problems arising in implanted joints. This pertains to both cement-less (the so-called "press-fit") and cemented implants.
We present a mathematical model for the dynamics of such problems. The contact process is assumed to include friction and wear between a viscoelastic body and a reactive foundation. Contact is described with a generalized compliance condition and friction with a general subdifferential law. We assume that the wear generation process takes place only on a part of the contact surface, and the wear rate is described by a generalized differential Archard condition that allows for the diffusion of the debris on the whole of the contact surface. This is the main novelty in the model. Such phenomena of wear diffusion can be found in many engineering settings, but in mathematical publications on contact and wear, it is tacitly assumed that the wear debris is removed from the surface once they are formed, which is the case some cases, such as car engines where the oil transports the debris away. The only mathematical works (that we are aware of) in which the wear debris remains on the surface and its diffusion is taken into account are [19,20], but there the contact surface was assumed to be planar. However, in most cases in applications, and those in joint replacements, the surfaces are curved. Moreover, in [19,20] the authors considered a quasistatic process and a moving foundation.
The novelty of this paper lies in that the process is assumed to be dynamic, the contact surface is a manifold and so we use of surface gradients and the Laplace-Beltrami operator instead of the linear diffusion equation. Also, we use a general nonmonotone subdifferential conditions to model friction, which is an extension of the classical formulation as a variational inequality with a subdifferential in the sense of convex analysis. In addition, the method of proof is new and very different from the usual one based on the use of results for variational inclusions.
The model for the processes consists of two coupled equations: the first is the dynamic equation of motion of a viscoelastic body and it contains a frictional multivalued term. The second one describes the diffusion of the wear debris on the contact surface of the body. Our key result is the theorem on the existence of a weak solution to the problem. In contrast to [8,10] (where the debris diffusion is modeled, but the contact surface is assumed to be flat) we do not use the Banach fixed point argument, but we base our approach on the Kakutani-Ky Fan-Glicksberg theorem that allows us to remove of the limitations on the constants present in the model at the cost of getting only existence, and not the uniqueness of a solution. In such a way we present a new way to obtain existence results for contact problems with friction and wear diffusion.
We remark here that we do not take into account adhesion effects in the model, and in many contact problems, one should also take into account the process of adhesion that is coupled with friction and wear diffusion. For instance, clinical practice shows that adhesion plays an important role at the bone-implant interface, and for further details we refer to [17,18] and the references therein.
The main mathematical difficulties of this paper lie in the formulation of the wear diffusion not on a subset of R 2 , but on a 2D manifold in R 3 . Also, due to the fact we do not impose any smallness condition on the constants in the model, we cannot use the Banach fixed point argument (such as it is done in [14]) that also asserts the solution uniqueness. In our approach we do not need any assumptions on the smallness of the data, but we obtain only the existence of a weak solution. Due to the rather general assumptions on the nonlinearities appearing in the problem, we are forced to use a truncation, and we first obtain the solutions to the truncated problem. We then obtain the necessary a priori estimates and remove this restriction by passing to the limit with truncation parameter.
The paper is organized as follows. Section 2 describes the 'classical model' for the process. We also describe shortly the equation for the wear diffusion on the contact manifold. Section 3 lists the assumptions on the problem data and derives its variational formulation, Problem P V . It is a system coupling an evolutionary differential inclusion for the displacements with a diffusion equation on the curved contact surface for the wear. Our main result, Theorem 3.2, states that under certain reasonable assumptions on the setting and problem data, there exists a solution of the variational problem, which is a weak solution for the 'classical' model. The proof to the main existence result is done in Section 4, and is based on the approach described above. Finally, Section 5 concludes with a short discussion and some open problems for further study.

The model
We consider a viscoelastic body that occupies a bounded domain Ω ⊆ R d , d = 2, 3 that is acted upon by volume forces and surface tractions. Although the case R 2 is of interest mathematically, in this case the contact surface is a curve and there doesn't seem to be applied interest in such a case, so we have d = 3 in mind. As a result, the body may come in frictional contact with a foundation and, consequently, a part of the contacting surface may undergo wear. The wear particles or debris produced in this process remain on the contact surface and undergo diffusion. Thus, grooves and surface damage occur causing changes in the shape and properties of the contacting surface. We construct a mathematical model for the evolution of the mechanical state of the body during the time interval [0, T ], where 0 < T < +∞. The unknowns in the problem are the displacements and the surface wear function. We refer to [20] for a more thorough discussion and additional details of the process. The main novelty here is that the contact surface is curved, while there and in [8,10] the contact surface was assumed to be flat, and moreover, here the process is dynamic.
We let Γ denote the boundary of Ω that is assumed to be Lipschitz continuous. We assume that Γ consists of three pairwise disjoint sets: Γ D where the body is held fixed and µ d−1 (Γ D ) > 0; Γ N where surface tractions act; and Γ C that is the potential contact surface, where friction and wear take place. The set Γ C is assumed to be a C 2 manifold with smooth boundary ∂Γ C . We note here that the assumption µ d−1 (Γ D ) > 0 is not essential, but it allows to avoid certain technical difficulties, such as the lack of the Korn inequality. We use the notation Ω T = Ω × (0, T ), Γ T = Γ × (0, T ), and similarly for Γ DT , Γ N T and Γ CT .
The body is held clamped on Γ D and so the displacement field vanishes there. A volume force of density f 0 acts in Ω T and surface tractions of density f N are applied on Γ N T . An initial gap function g can exist between the potential contact surface Γ C and the foundation and it is measured along the outward normal ν.
We denote the displacement vector by u : Ω × [0, T ] → R d , the velocity vector by v = u ′ , where the prime represents the time derivative, the linearized strain tensor by so that ε ′ (u) = ε(v), and the stress tensor by σ = (σ ij ), all defined on Ω T . We write the normal components and tangential vectors on the boundary Γ C as and the normal and tangential stresses as We assume that the material is viscoelastic with linear constitutive relation are the viscosity and elasticity tensors, respectively. Thus, and summation over repeated indices is implied. The viscosity and elasticity tensors satisfy the following assumptions.
We remark that the presence of the viscosity tensor A that is coercive is crucial in the proof of our main result. Some problems with hyperbolic inclusions have been recently studied in e.g., [16,23], yet in our case they are not applicable and it remains an open problem to remove the viscosity term and consider a purely elastic material.
The displacement u satisfies the momentum law where f 0 : Ω × [0, T ] → R describes a volume force. The body is clamped on Γ D , hence, and the traction f 2 ise applied on Γ N , We turn to describe the wear process and note that in [10,19,20] the contact surface Γ C was divided into two sub-domains D d and D w and the wear took place only on the part D w , while the diffusion of particles took place on the whole of Γ C . In this work we assume that wear is generated and diffuses on Γ C , however, we note that it is straightforward to restrict wear generation to a part of Γ C by introducing the appropriate characteristic function, as was done in the articles above.
Before we continue, since we are interested in the diffusion of the wear debris on the surface, we need to introduce the concepts and notation related to diffusion on curved surfaces. We follow [7] (see also the references therein) and in particular, we refer the reader to [12, p.388], for the definition of hypersurfaces in R n and surface gradients on them. Let S be a smooth surface in R d , if G is a smooth function defined in a neighborhood of S, the surface or tangent gradient on S is defined as where G ν = ν · ∇G is the normal derivative of G on S, recalling that ν denotes the unit outer normal vector to S. Thus, the surface gradient at x ∈ S is the projection of the gradient at x onto the tangent plane to S at x. Note, that for the above definition of ∇ S G to make sense, we need to extend G from S to an open neighborhood in R d , however, such an extension always exists for smooth S and the value of ∇ S G does not depend on the choice of the extension (see, e.g., [5]). If we denote the components of the surface gradient by then the Laplace-Beltrami operator, which describes the spatial part of diffusion on the surface, is defined by the surface divergence of the surface gradient, i.e., where k = 1, . . . , d, and summation is implied. Next, we assume that the manifold S has a smooth boundary Γ S = ∂S and denote by ν S the unit outer normal to S on Γ S . Then, Green's formula on S is given by (see, e.g., [6]) and holds for each pair (ϕ, ψ) of smooth functions defined in a neighborhood of S. In our setting, S = Γ C and Γ S = ∂Γ C . For the sake of somewhat simplified notation we use ∇ Γ for the gradient and ∆ Γ for the Laplace-Beltrami operator on Γ C . We use the notion of the Sobolev space H 1 (S) of functions on the manifold S, i.e., functions in L 2 (S) that have their surface gradients belongs to L 2 (S) d , see [2] for the definition and properties of these functions on manifolds without boundary, and [1] for manifolds with smooth boundary, which is the case here.
To describe the wear process and its diffusion, we introduce the wear function θ that is defined on the contact surface Γ C × [0, T ], and its evolution is governed by a parabolic differential equation, and a zero flux boundary condition on ∂Γ C , since the debris cannot leave Γ C . We note that the rate form of the usual Archard's law of wear (see, e.g., [21]) states that the rate of surface wear is proportional to the frictional traction, and the relative velocity, i.e., the power of the friction resistance force, and is given by where η is the wear rate constant, µ is the friction coefficient, the function p ν describes the normal stress, and more details are below, and v τ (t) is the tangential velocity. As was done in [10,19,20], we extend the Archard law and allow diffusion of the wear debris on the surface, i.e., we generalize h w to a function satisfying It is straightforward to see that when the function p ν has at most linear growth, the wear source function satisfies this assumption. Then, the extended version of the Archard law for a pointwise wear process of growth and diffusion on Γ C is given by where κ is the wear diffusion constant. We note that the debris source h w depends on the wear and the surface speed, since the wear changes the surface geometry, the debris changes the friction resistance, and the friction coefficient is known to depend on the speed.
We turn to the contact conditions on Γ C . We describe the contact process on Γ C by a general condition of the form We impose the following hypotheses on h ν : An example of a law satisfying this conditions is the normal compliance condition (see, e.g., [21] and the references therein), where p ν is the normal compliance function that vanishes for negative arguments, since then there is no contact between the body and the foundation at the point of Γ C . In the literature it was typically assumed to be of the form where (·) + was the positive part, λ νc was assumed to be a large number and m ≥ 1 was the normal compliance exponent (see also [21]). We describe friction with a general subdifferential law where j is a locally Lipschitz function and ∂j stands for its Clarke subdifferential (see Section 3 below for details). We suppose that h τ and j τ satisfy x ∈ Γ C and moreover ζ · ξ 0 for ζ ∈ ∂j τ (x, ξ) for all ξ ∈ R n a.e x ∈ Γ C .
As an example of such a friction law, one may use a version of the Coulomb law, where µ is the coefficient of friction and µp ν is the friction bound, and That is, frictional resistance takes place only when there is relative motion and then it opposes it. We can write the condition in a condensed form as an inclusion where ∂|r| is the convex subdifferential of |r|, i.e., We use the formalism of Clarke subdifferentials in the friction law to account for possible nonmonotonicity in the relation between the tangential velocity and the friction force density. This represents the fact that kinetic friction can be less than static friction, i.e., a drop of the friction force can occur when motion starts. Finally, the initial conditions for the displacement, velocity and wear functions are,

Variational formulation
We turn to the variational formulation of problem (2.1)-(2.9). To that end, we first introduce the concepts that are needed belowl, and then the variational formulation. In what follows, i, j = 1, . . . , d everywhere, the summation convention over repeated indices is used, and an index following a comma indicates a partial derivative. For a reflexive Banach space E, we denote by ·, · E * ×E the duality pairing between the dual space E * and E. If E is a Hilbert space, then the scalar product in E is denoted by (·, ·) E . Throughout this paper, we denote by C a generic positive constant that depends on the problem data and may change its value form line to line. By | · | we denote the Euclidean norm in R d or S d , the space of symmetric d × d matrices.
For a fixed and finite T > 0 we define the following standard time-dependent spaces: The Clarke subdifferential of a locally Lipschitz functional ϕ : R d → R is given by (see [3]) where N ϕ is a set of measure zero, outside of which ϕ is differentiable, and N is any set of measure zero. It is possible to generalize the notion of the Clarke subdifferential to functionals defined on Banach spaces, cf., [3,4,14], but for our purposed it is sufficient to consider this definition on R d . Now, we define the operators A, G : V → V * by We assume that f 0 (t) ∈ L 2 (Ω) d and f 2 (t) ∈ L 2 (Γ C ) d and this allows us to define f : (0, T ) → V * as Applying the Green formula and the usual manipulations, we are able to derive the following weak formulation of the problem governed by (2.1)-(2.9).
By ξ(t) ∈ S 2 ∂jτ (v τ (t)) we understand a L 2 -measurable selection out of the subdifferential ∂j at v τ . We also write, for the sake of simplicity, H(h) as a collection of the hypotheses H(h w ), H(h ν ) and H(h τ ).
We are now able to state the main theorem of this paper.
We conclude that the model (2.1)-(2.9) has a weak or variational solution. The uniqueness of the solution remains an unresolved question.

Proof of Theorem 3.2
In this section we prove the existence theorem. The idea of the proof is as follows. First, we decouple the coupled Problem 3.1 by replacing the coupling terms with given functions and introduce truncation operators. We obtain the existence of solutions for the decoupled and truncated problems independently. Then, we apply the Kakutani-Ky Fan-Glicksberg fixed point theorem to show the existence result for the original problem. Finally, we pass to the limit with the truncation parameter. In the proof we always assume (H1)-(H8), and that u 0 ∈ V, v 0 ∈ H, θ 0 ∈ L 2 (Γ C ), f 0 ∈ L 2 (0, T ; V * ), f 2 ∈ L 2 (0, T ; L 2 (Γ N ) d ), and κ > 0, so we do not repeat these assumptions in the auxiliary lemmas below.
We start by recalling the fixed point theorem. Next, for l > 0, we define truncation operators N l : R d → R d and M l : R → R by x |x| l, |x| > l. The following lemma is straightforward to show, and we present the proof for the sake of completeness. Proof. We present the proof only for N l . Let x, y ∈ R d and we consider the three cases: |x|, |y| l, |x| > l, |y| l, and |x|, |y| > l. In the fist case, we immediately obtain the result. In the second case, we calculate the inner products in R d , In the last case, Now we fix l > 0 (large), choose the functions v ∈ W, ξ ∈ L 2 (0, T ; L 2 (Γ C ) d ) and θ ∈ W Γ , and let u given by (3.5) using v. Consider now the following two auxiliary problems.
Problem 4.4. Find the wear function θ ∈ W Γ such that We note that by using the given functions and the truncations, the two problems are uncoupled. Proof. For the proof of Lemma we refer to [14]. Proof. For the proof, we refer to classical results on parabolic problems, see, e.g., [13].
In the next step, we introduce the following coupled, but still truncated problem. In what follows, we check all the assumption of Theoremt 4.1, and summarize the steps in the lemmas. First, we derive the necessary a-priori estimates.
Proof. We choose w = v(t) in (4.1) and then it follows from (H3), (H5), (H6) and the Cauchy inequality with ε > 0 that for t ∈ (0, T ), Integrating (4.10) over (0, t) for t ∈ (0, T ) and choosing appropriate value of ε yields (4.10) Next, we choose η = θ(t) in (4.2), and then it follows from (H4) that Again, integrating (4.11) over (0, t) for t ∈ (0, T ) we obtain Using the Gronwall inequality we get for t ∈ (0, T ) Combining (4.12) and (4.13) it follows that hence from (4.13) and (4.14) we conclude that Now, we need to show the existence of ξ ∈ L 2 ((0, T ) × Γ C ) d such that ξ(t) ∈ S 2 ∂jτ (v τ (t)) for a.e. t ∈ (0, T ). To show this it is sufficient to prove the existence of a measurable selection, since the integrability, as well as the bound (4.9) follow from (H8). However, the existence of a measurable selection of the subdifferential follows from [4, Theorem 5.6.39], as the Clarke subdifferential of the locally Lipschitz integral functional is nonempty and its elements are measurable on the one-hand, and on the other-hand they are selections of the multifunction ∂j(v(x, t)), see also [11].
Next, we define the space Z = W × W Γ × L 2 (0, T ; L 2 (Γ C ) d ) and consider the solution operator Λ : Z → 2 Z , which assigns to a triple (v, θ, ξ) a triple (v, θ, ξ), where v and θ are the solutions of Problems 4.3 and 4.4, respectively, and ξ is a L 2 -measurable selection out of ∂j(v τ ). We have the following lemma. Lemma 4.9. There exist positive constants R 1 , R 2 and R 3 such that Λ(B) ⊂ 2 B , where the set B = B(R 1 , R 2 , R 3 ) is given by Proof. We use the estimates above and choose  Proof. The result follows from the convexity in the definition of the Clarke subdifferential, see eg. [3], the existence of ξ given in Lemma 4.8, and the existence and uniqueness of v and θ established in Lemmas 4.5 and 4.6.
Proof. We choose three sequences such that v n → v weakly in W, θ n → θ weakly in W Γ and ξ n → ξ weakly in L 2 (0, T ; L 2 (Γ C ) d ). Define v n , θ n and ξ n as, respectively, the solutions of Problems 4.3 and 4.4 corresponding to v n , θ n , ξ n , and the L 2 selection of ∂j τ (v nτ (x, t)). Assume that v n → v weakly in W, θ n → θ weakly in W Γ and ξ n → ξ weakly in L 2 (0, T ; L 2 (Γ C ) d ). We need to show that v and θ are the solutions of Problems 4.3 and 4.4 that correspond to v, θ and ξ, and that ξ is the L 2 selection of ∂j τ (v τ (x, t)).
First observe that the compactness of the embedding i : V → H 1−δ (Ω) d together with the Aubin-Lions lemma imply that the tangential components of the traces satisfy v nτ → v τ strongly in L 2 (0, T ; L 2 (Γ C ) d ).
To show that v and θ are the solutions of Problems 4.3 and 4.4 corresponding to v, θ and η, we need to write (4.1) and (4.2) for v n and θ n , and then pass to the limit n → ∞. It is clear that θ n (0) → θ(0) weakly in L 2 (Γ C ) and v n (0) → v(0) weakly in H, which implies that θ and v satisfy the same initial conditions as θ n , v n . Moreover, the following hold, To show the weak sequential closedness of Gr(Λ), we show the convergence of all the terms. The Aubin-Lions lemma implies v n → v strongly in L 2 (0, T ; H) and θ n → θ strongly in L 2 (0, T ; L 2 (Γ C )), (4.19) v n → v strongly in L 2 (0, T ; L 2 (Γ C ) d ), u n → u strongly in L 2 (0, T ; L 2 (Γ C ) d ), (4.20) θ n → θ strongly in L 2 (0, T ; L 2 (Γ C )). (4.21) Now, Lemma 4.9, (4.18) and (4.19), and the linearity of the duality pairings and the linearity and boundedness of operators A and G, we obtain that for every w ∈ L 2 (0, T ; V ), and every η ∈ L 2 (0, T ; H 1 (Ω)), as n → ∞. Next, we deal with the boundary integrals. To simplify the presentation, we omit the time dependence of the functions. We write The weak convergence ξ n → ξ in L 2 (0, T ; Moreover, by the continuity of h τ , the strong convergences (4.20) and (4.21) and the Lebesgue dominated convergence theorem, we find where the last term converges to zero as n → ∞. Hence, (4.22) and (4.23) yield that as n → ∞, By the direct application of the Lebesgue dominated convergence theorem and the continuity of h ν and h w , we find as n → ∞. This completees proof of the lemma.
The next step is essentially the last one.
Since the topology is weak, we need the following argument to show that this set is closed. But, , which is bounded, closed and convex in the reflexive space Z × Z,therefore, Gr(Λ| B(R 1 ,R 2 ,R 3 ) ) is sequentially compact, and so (w − Z) × (w − Z) is compact and (w − Z) × (w − Z) is closed. Taking into account Lemma 4.10, the assertion of the lemma follows now directly from Theorem 4.1.
We have shown that all the assumptions of the fixed-point theorem, Theorem 4.1, hold true and that establishes the following theorem, which guarantees the existence of a solution of the truncated problem. The last step in the proof of our main theorem is to show that we can remove the truncation operators from Problem 4.7.
Proof of Theorem 3.2. We need to obtain the relevant estimates on a solution v, θ of Problem 4.7 that are independent of the truncation parameter l. To that end, we choose η = v(t) in (4.3) and using again the Cauchy inequality with ε > 0, we obtain for t ∈ (0, T ). The hypotheses (H5)-(H7) and an appropriate choice of ε > 0 in (4.24) yields for t ∈ (0, T ). Straightforward manipulations that use the Cauchy inequality with ε again, the fact that u(t) = u 0 + t 0 v(t) dt and the inequality in [4], Lemma 8.4.12 show that Integrating (4.26) over(0, t), t ∈ (0, T ) and using (H1) and (H3), we get Choosing η = θ(t) in (4.5), and applying (H4) for t ∈ (0, T ), leads to the estimate Next, by using the continuous embedding H 1 (Ω) → L 4 (∂Ω), we find Integrating (4.31) over (0, t), t ∈ (0, T ), we find where C is independent of l. This estimate is crucial for the proof of the theorem. Now, let (v n , θ n , ξ n ) be a solution of Problem 4.7 with the truncation constant l = n. Then, (4.37) and the Aubin-Lions lemma imply that there is a subsequences such that v n → v strongly in L 2 (0, T ; H) and in L 2 (0, T ; L 2 (Γ C ) d ), and θ n → θ strongly in L 2 (0, T ; L 2 (Γ C )). Passing to the limit with the multivalued term follows exactly as in the proof of Lemma 4.11, so, we pass to the limit with all terms in Problem 4.7. To finally remove the truncations, we need to check that h ν (N n (u n )) → h ν (u), h τ (N n (u n ), N n (v n ), M n (θ n )) → h τ (u, v, θ), h w (N n (u n ), N n (v n )) → h w (u, v), strongly in L 2 (0, T ; L 2 (Γ C )). By continuity of h ν , h τ , h w it is enough to show that N n (v n ) → v strongly in L 2 (0, T ; L 2 (Γ C )), (4.38) N n (u n ) → u strongly in L 2 (0, T ; L 2 (Γ C )), (4.39) M n (θ n ) → θ strongly in L 2 (0, T ; L 2 (Γ C )). This proves (4.38). To show (4.39) and (4.40), we repeat similar calculations for u n and θ n . Hence, we can pass to the limit with truncation parameter l = n → ∞ in all terms. This completes the proof of Theorem 3.2.
Thus, the model has at least one solution. The question of uniqueness remains unresolved, but in view of the complexity of the system and its nonlinearities, it is unlikely. Indeed, the uniqueness of solution to Problem 3.1 does not follow from the presented argument, as it does in a case of the Banach fixed point theorem. Moreover, we suspect that uniqueness would require additional smallness assumptions on the data and stronger assumption on the functions h τ , h ν , h w .
As has been already mentioned, establishing an existence theorem for a purely elastic model is of considerable mathematical interest.