The modeling and the analysis of frictional contact problems represent important topics both in Engineering Sciences and Applied Mathematics, see for instance the numerous works in these fields of research     and the references therein. The mathematical study of contact and friction models concerns the variational and the numerical analysis with also significant contributions related to numerical simulations. In the construction of the models, various contact conditions have been considered. We can cite, for instance, the so-called Signorini condition, introduced in , the so-called normal compliance contact condition, introduced in  and used in a large number of papers, see  . Recently, a more general contact condition, called the normal compliance condition restricted by unilateral constraint introduced in , models the contact with an elastic-rigid foundation. The mathematical analysis of models involving the frictionless contact condition with normal compliance and unilateral constraint can be found in . When friction is considered, the unique solvability of the variational problems can be proven by considering a smallness assumption of the friction coefficient, see for instance . Recently, the study of a model of frictional contact problem between a linearly elastic body and a moving rigid foundation, considered in , has shown both the obtention of multiple solutions when the coefficient of friction is larger than a critical value, and the existence and uniqueness of the solution at the low coefficient of friction.
The aim of this paper is to study the numerical analysis of the contact problem with a sliding version of Coulomb’s law of dry friction for rate-type viscoplastic materials within the framework of the Mathematical Theory of Contact Mechanics. We model the material’s behavior with a constitutive law of the form
where denotes the displacement field, represents the stress tensor and is the linearized strain. Here is the elasticity operator, assumed to be nonlinear, and represents the operator relaxation operator, assumed to be linear, with a fully discrete scheme for the numerical approximation of the problem, numerical simulations and its validity in the study of a two-dimensional numerical example. Therefore, the contact law with normal compliance and unilateral constraint was associated with a sliding version of Coulomb’s law of dry friction. Furthermore, both the material constitutive law of the body and the frictional contact model is characterized by memory terms in order to take into account physical relaxation behaviors. Such kind of history-dependent problems has been considered in .
Here, our goal is to provide the numerical analysis of the frictional contact model with normal compliance and unilateral constraint and to illustrate the error estimate of the discretization by numerical simulations. The mathematical model is based on a viscoelastic constitutive law with a long memory, contact conditions combining normal compliance, memory term, unilateral constraint and a frictional sliding version of Coulomb’s law. This nonstandard mathematical model can be formulated by a history-dependent quasi-variational inequality for the displacement field.
The rest of the paper is structured as follows. In Section 2 we present the mathematical model by describing its equations and boundary conditions. Some notation, assumptions and preliminary material have been introduced in order to derive a variational formulation of the problem. Section 3 is devoted to the numerical approximation and the numerical analysis of the variational problem considered in the previous section. Under certain solution regularity assumptions, we derive an optimal order error estimate that represents the main result of this work. Finally, in Section 4 we consider some numerical simulations in the study of a two-dimensional problem and we provide a numerical validation of the error estimate.
2. Frictional Contact Problem and Variational Formulation
First, we present some preliminary material use full for the setting of the problem. Let a regular domain of ( ) with its boundary that is partitioned into three disjoint measurable parts , and , such that . We use the notation for a typical point in and we denote by the outward unit normal defined almost everywhere (a.e.) on . We denote by , , and the displacement vector, the stress tensor, and the linearized strain tensor, respectively. Here and below the indices run between 1 and d and, unless stated otherwise, the summation convention over repeated indices is used. We note that sometimes we do not indicate the dependence of various vectors and tensors on the spatial variable and we recall that the components of the linearized
strain tensor are given by where the index that
follows a comma indicates a partial derivative with the corresponding component of the spatial variable , e.g. . We denote by t the time variable and a dot superscript represents the time derivative with respect to the time variable t, e.g. . Furthermore, we use the notation for the set of positive integers and will represent the set of non-negative real numbers, e.g. . Then, we denote by the space of second-order symmetric tensors on . The inner product and norm on and are defined respectively by
For all vector , we denote by and the normal and tangential components of on given by
We also recall that, if is a regular function, then the normal and tangential components of the stress field on are defined by
In this work, we consider a viscoelastic body that occupies the bounded domain with , its boundary. The body is clamped on and, therefore, the displacement field vanishes there. A volume force of density acts in , surface tractions of density act on . On , the body is in frictional contact with a moving obstacle, the so-called foundation. We suppose that the foundation is plane and moves steadily, i.e. its velocity is assumed to be larger than the tangential velocity on the surface contact (i.e.
), where denotes a given unitary vector in the tangential
plane and the value is also given. The contact between the body and the foundation is modeled by multivalued normal compliance and a unilateral constraint. The associated friction is based on a version of Coulomb’s law of dry friction, in which contact surfaces are assumed to be in relative slip status. The problem is studied in the interval of time .
Thereby, let us consider the formulation of our quasi-static frictional contact problem defined as follows.
Problem . Find a displacement field and a stress field such that, for all ,
and there exists a normal reaction that satisfies
Now, we shortly describe the physical meaning of relations (3)-(8). Equation (3) represents the viscoelastic constitutive law with long memory in which is the elasticity operator and is a relaxation tensor. Details and mechanical interpretation concerning such kinds of laws can be found in , for instance. Equation (4) represents the equation of equilibrium in which Div denotes the divergence operator for tensor valued functions ( ). Conditions (5) and (6) are the displacement boundary condition and the traction boundary condition, respectively. Finally, (7) and (8) represent the friction Coulomb’s law and the multivalued normal compliance contact condition with unilateral constraint and memory term, respectively. The friction condition (7) represents a regularized form of a version of Coulomb’s law in slip status where represents the coefficient of friction and is a regularization operator. Here, and represent the normal contact stress and the tangential friction stress on the contact surface , respectively. Condition (8) represents a version of the contact boundary conditions with normal compliance and unilateral constraint, in which the memory effects of the foundation are taken into account. Note that condition (8) models the contact with a foundation made of a rigid material and covered by a layer of soft material (asperities) of thickness g with a thin crust. Let us describe the different terms. g is the maximal penetration allowed in the foundation and represents the size of the soft material. is the normal compliance condition with memory term, where the normal compliance function p is a Lipschitz continuous increasing function which vanishes for a negative argument and describes the memory properties of the foundation. Here, F and are positive functions related to the history of the soft material.
We now turn to the variational formulation of Problem and, to this end, we consider standard notation for the Lebesgue and the Sobolev spaces associated to and . Also, we introduce the spaces
The spaces Q and Q1 are Hilbert spaces with the canonical inner product given by
and the associated norms and . Since , V is a real Hilbert space with the inner product
and the associated norm . By using the Sobolev trace theorem, a positive constant exists such that
In addition, for we denote by its normal component, in the sense of traces. Let be a linear continuous operator. Then, there exists a positive constant such that
Since is regular, the following Green’s formula holds:
Let us denote by the space of fourth-order tensor fields given by
Let us recall that is a real Banach space with the norm and, the following inequality holds,
Based on definitions of the norms and , proof of the inequality (12) is obtained by using simple calculations.
Furthermore, we need to consider the sets of admissible displacement and admissible constraints defined by
Let us consider the assumptions on the operators, tensors and data. To this end, we assume that the elasticity operator and the relaxation tensor satisfy the following conditions.
The densities of body forces and surface traction have the following regularities
The normal compliance function p and the surface yield function F satisfy
Finally, the surface memory function and the coefficient of friction verify
In what follows we assume that are sufficiently regular functions and satisfy (3)-(8). Here, let and be given. First, by using Green’s formula (11) and the equilibrium equation (4), we obtain that
To simplify the notation, we defined the operators and by
where represents the projection operator. With the inclusion , the operator R is well defined. Then we consider the space as well as operators , and , for all , defined respectively by
Finally, the functions and are defined by
Using the previous notation and substituting (3) in (22), we obtain the following variational formulation of Problem .
Problem . Find a displacement field and a stress field , for all , such that,
Under the assumptions (15)-(21) and an additional smallness assumption on the friction coefficient , a result of existence and uniqueness for the variational problem was provided in . Based on this previous variational formulation, our goal in the next section is to provide the numerical analysis of this specific and non-trivial contact problem.
3. Variational Approximation and Error Analysis
This section is devoted to the numerical discrete approximation of the Problem . In particular, an optimal error estimate is provided and represents the main result of the paper. More precisely, we are interested in solving the Problem over a finite time interval , with arbitrary but fixed. Thus, let N be a positive integer; we define the size of the time step and we consider the uniform temporal discretization characterized by the time instants for . We use the notations , . For simplicity, we assume that is a polygonal domain. Let consider a regular family of partitions of into triangles that are compatible with the partition of the boundary into , , and , in the sense that if the intersection of one side of an element with one of the three sets has a positive surface measure, then the side lies entirely in that set. Then we consider linear element spaces corresponding to ,
where stands for the space of polynomials of a degree less than or equal to 1 on T. We define a non-empty, convex, and closed set, approximating U by
where denotes the spatial discretization parameter. For the discretization of the integral term of Problem , we use a variant of the trapezoid method
where a prime indicates that the first and last terms in the summation are to be halved. Then, we introduce the approximations of operators S, M and N,
, for . For , we define .
From now on, c will represent various positive constants that do not depend on h and k and whose value may change from line to line. Then the fully discrete scheme of Problem is the following.
Problem . Find such that ,
Note that the existence and uniqueness of the solution of the Problem can be provided by using the same arguments as for the solvability of the variational problem .
We now focus on the error analysis between the solutions to Problems and . Taking and in (31), we have
Then, we consider the term
Furthermore, by using (15(c)), we obtain
We combine the inequalities (37), (38) and (39) to obtain
Now, we proceed to the estimation of each terms of (40). First, we use the definition (23) and the inequalities (10), (15(b)) and (18(c)) to obtain
For the term , we use arguments similar, to obtain
where and are constants, for more details about these constants.
By (34)-(36), we have
Moreover, we use the triangular inequality as well as the trapezoid method introduced in (33) to obtain
Under the following assumptions
And using the hypothesis (19)(b) of the F function, we have
Then we combine the inequalities (45), (46), (49), (50) and the trace inequality (10) to obtain
By similar reasoning applied for the term , we obtain that
Finally, we gather together the inequalities (40), (44), (51) and (52) to deduce that
Under the assumption of smallness for the uniqueness of the Problem , and by using the elementary inequality , , , we have
Now let us consider the following Gronwall inequalities.
Lemma 1. Let be given. For a positive integer N we define . Assume that and are two sequences of nonnegative numbers satisfying: for a positive constant independent of N or k. Then there exists a positive constant c, independent of N or k, such that: .
Applying Lemma 1, we obtain from (53) that if k is sufficiently small, then
To proceed further, we make the following solution regularity assumptions:
Then it can be shown that the relations (3)-(8) are satisfied pointwise a.e. and we have
So (54) reduces to
Let be the linear finite element interpolant of the solution . As the solution , i.e., , then . The standard finite element interpolation theory yields (cf.  ) and we have
In conclusion, we have shown the following result.
Theorem 2. Assume k is sufficiently small. Then under the assumptions (47), (55) and (56) concerning the regularity of the solution , we have the following optimal error estimate
with a positive constant c independent of k and h.
4. Numerical Simulations
This section provides computer simulation results on the contact Problem , including numerical evidence of the theoretical error estimates obtained in Section 3 for the discrete approximation of the variational Problem . The solution of Problem is based on numerical methods described in . For a large overview on numerical methods to solve contact problems, we can refer for instance to   .
Numerical example. The physical setting used for Problem is depicted in Figure 1. Here, we consider the frictional contact between a deformable body and a moving foundation. This specific foundation is composed of a rigid material covered by a thin crust and a deformable layer of asperities of depthg. Here g represents the maximum value of the allowed penetration in the foundation. When this value of penetration is reached, the contact follows a unilateral condition without any additional penetration. This kind of foundation is characterized by contact condition (8). Since the foundation is moving the friction condition is in a slip status within Coulomb’s form (7). The deformable body is
Figure 1. Reference configuration of the two-dimensional example.
a rectangle, , and its boundary is split as follows: , , .
The domain represents the cross-section of a three-dimensional linearly viscoelastic body subjected to the action of tractions in such a way that a plane stress hypothesis is assumed. On the part the body is clamped and, therefore, the displacement field vanishes there. Vertical compressions act on the part of the boundary and the part is traction-free. Constant vertical body forces are assumed to act on the viscoelastic body. The body is in frictional contact with an obstacle on the part of the boundary.
The compressible material’s behaviour of the domain is governed by a viscoelastic constitutive law of the form (3). In addition, we assume that the material is homogeneous and isotropic; then, the elasticity tensor and the relaxation tensor have the following forms
where the coefficients E and are Young’s modulus and the Poisson’s ratio of the material, respectively, and is a viscosity parameter. denotes the Kronecker symbol.
For the numerical simulation of Problem , the data concerning the material are the following: , , , , on , , , with . For the contact boundary conditions on , we recall that the friction is in a sliding status and the normal contact response follows a multivalued normal compliance condition with respect to the normal displacement and for which the maximal penetration is restricted by a unilateral constraint. Then the data concerning this specific frictional contact model (7) and (8) are the following: , , , , and .
Numerical solution of Problem . In Figure 2, we present the deformed configuration as well as the interface forces on .
On the right extremity of the boundary , we can see that a large number of contact nodes of are in several contact statuses with the foundation. We note that some of these nodes have first broken the crust (crust contact), others have crushed the asperities (normal compliance) and finally someone have reached the maximum value g of penetration (unilateral contact).
Errors and numerical convergence orders. The aim of this part is to illustrate the convergence of the discrete scheme and to provide numerical evidence of the optimal error estimate obtained in Section 3. To this end, we computed a sequence of numerical solutions by using uniform discretization of Problem according to the spatial discretization parameter h and the time step k, respectively. For instance, the deformed configuration and the interface forces plotted in Figure 2 are obtained for and and corresponds to a problem with 17028 degrees of freedom.
The numerical estimations of computed for several discretization parameters of h and k, have been estimated by using the energy norm
defined by . Since the exact solution can
not be calculated analytically, we consider as “reference” solution a numerical solution that corresponds to a fine approximation of Problem . For this procedure, the boundary of is divided into 1/h equal parts and the time interval is divided into 1/k time steps. We start with and which are successively halved. The numerical solution corresponding to and was taken as the “reference” solution. This fine discretization corresponds to a problem with 66820 degrees of freedom; the simulation runs in 129521 (expressed in seconds) CPU time on an IBM computer equipped with Intel Dual core processors (Model 5148, 2.33 GHz). The numerical results presented in Figure 3 provide good numerical evidence
Figure 2. Deformed mesh and interface forces on .
Figure 3. Numerical errors.
of the theoretically predicted order convergence of the numerical solution measured in the energy norm.
 Sofonea, M. and Souleiman, Y. (2015) Analysis of a Sliding Frictional Contact Problem with Unilateral Constraint. Mathematics and Mechanics of Solids, 22, 324-342.
 Sofonea, M. and Souleiman, Y. (2016) A Viscoelastic Sliding Contact Problem with Normal Compliance, Unilateral Constraint and Memory Term. Mediterranean Journal of Mathematics, 13, 2863-2886.
 Eck, C., Jarusek, J. and Stará, J. (2013) Normal Compliance Contact Models with Finite Interpenetration. Archive for Rational Mechanics and Analysis, 208, 25-57.
 Capatina, A. (2014) Variational Inequalities and Frictional Contact Problems. In: Gao, D.Y. and Ogden, R.W., Eds., Advances in Mechanics and Mathematics, Springer, Cham, 1-6.
 Oden, J.T. and Martins, J.A.C. (1985) Models and Computational Methods for Dynamic Friction Phenomena. Computer Methods in Applied Mechanics and Engineering, 52, 527-634.
 Souleiman, Y. (2021) Convergences and Numerical Analysis of a Contact Problem with Normal Compliance, Unilateral Constraint. African Journal of Mathematics and Computer Science Research, 14, 13-23.
 Jarusek, J. and Sofonea, M. (2008) On the Solvability of Dynamic Elastic-Visco-Plastic Contact Problems. ZAMM-Journal of Applied Mathematics and Mechanics, 88, 3-22.
 Barboteu, M., Cheng, X.L. and Sofonea, M. (2016) Analysis of a Contact Problem with Unilateral Constraint and Slip-Dependent Friction. Mathematics and Mechanics of Solids, 21, 791-811.
 Sofonea, M., Han, W. and Barboteu, M. (2013) Analysis of a Viscoelastic Contact Problem with Multivalued Normal Compliance and Unilateral Constraint. Computer Methods in Applied Mechanics and Engineering, 264, 12-22.
 Haslinger, J., Hlavácek, I. and Necas, J. (1996) Numerical Methods for Unilateral Problems in Solid Mechanics. Handbook of Numerical Analysis, 4, 313-485.