{'fraction_non_alphanumeric': 0.03352550748510247, 'fraction_numerical': 0.011433554575844194, 'mean_word_length': 4.667764964305327, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
SIMULATION OF ADHESIVE FRICTIONAL CONTACT OF SOFT MATERIALS BY BI-POTENTIAL METHOD

Abstract. We propose a numerical formulation for the simulation of frictional contact with interfacial adhesion between soft bodies. In this work, the adhesive interface behaviour derives from a free surface energy and a pseudo-potential of the surface dissipation, which describe both the de-bonding process of the adhesive links due to tangential and normal interface deformation. Subsequently, a complete contact and friction law with the account for interface adhesion is formulated. Then the local contact nonlinear equations are solved using a Newton-like algorithm within the bi-potential framework. Numerical examples are performed to demonstrate the capacity of the proposed approach.


INTRODUCTION
Dynamic frictional contact involving soft materials and interface adhesion represents a frequent phenomenon of contact. In the area of numerical modelling, although many efforts have been exerted over the last decade to develop adhesive contact algorithms [1,10,12,18], literature on complete modelling schemes that can capture the entire bonding and de-bonding process of soft adhesives is still premature [9]. Such modelling schemes, in general, should involve proper definition of contact laws that reflect behaviours of the reversible interface adhesion, and development of robust and stable resolution algorithms, that can deal with the computational difficulties inherent to the problem non-linearities. Concerning the adhesive interface laws, a number of models have been developed over the past decades. The most prominent ones include Johnson-Kendall-Roberts (JKR) [8], Maugis-Dugdale (MD) [11], and Greenwood and Johnson models [5]. These models, proposed as early as the 1970s and became the reference in this area ever since, provided the theoretical basis for the contact and friction modelling of adhesive interfaces. Recent advances in this area include [6,9,13], which successfully represented reactions of adhesive interfaces under complex loads with tension and shear.
In this work, in order to simulate complex interface behaviours with adhesion, we adopt the RCC model [16,17,18] which is a cohesive interface model coupling adhesion, friction and unilateral contact. It is based on a adhesion intensity parameter β presented first in the work of Frémond [4]. Varying between 0 and 1, β characterizes the damage level of the interface adhesive bonds (0 refers to the state of complete de-bonding, 1 refers to complete bonding). Then the adhesive interface behaviour, evolving with β, derives from a free surface energy and a pseudo-potential of the surface dissipation [2]. A complete contact and friction law with extended Signorini condition and modified Coulomb friction law, accounting for adhesion is formulated.
However, the inherently severe non-linearities due to contact dynamics, and the non-smooth and multivalued nature of the adhesive interface law give rise to significant difficulties in numerical resolution, and makes it necessary to adopt robust, and stable algorithms to ensure iteration convergence, solution accuracy and balanced efficiency. Based on the augmented Lagrangian method, we propse a bipotential theory which couples the two variational inequalities of the unilateral contact and friction law into one single displacement based variational principle with one unique inequality. First introduced in the 1990s, the approach has recently been extended to problems involving elastic-plastic contact and friction [3,20], and interface wear [14,15].
The aim of this work is to propose a numerical formulation for the simulation of contact and friction between hyperelastic systems with adhesive interface. The article starts with a brief description of the contact kinematics. Then, we present the complete framework of the adhesive contact law, which is an extension of Signorini law and Coulomb rules. Then we present its implementation within the bipotential framework. In Section 3, the complete finite element formulation of the problem, including the resolution algorithm, is provided. To validate the framework, we present numerical examples in Section 4. In the end, a few conclucing remarks are drawn in Section 5.
2 Problem setting

Contact kinematics
We describe in this section the geometric definitions and notations related to the contact kinematics. Let's consider two deformable bodies B 1 and B 2 coming into contact with N c contact points. Each body is discretized with finite elements with nodal positions represented by X 1 (for B 1 ) and X 2 (for B 2 ) defined in the global coordinate system. Contact points belonging to B 1 are denoted by P α 1 (α = 1, 2, ...N c ), and accordingly P α 2 . We consider on each P α 1 a local orthogonal coordinate system, formed by T 1 , T 2 and N, representing respectively the tangential, and normal direction vectors defined with respect to the global coordinates. Therefore, P α 2 can be seen as the projection point of P α 1 on B 2 following N. We can build the relative position between P α 1 and P α 2 by with X(P α 1 ) and X(P α 2 ) the position vectors of P α 1 and P α 2 in the global coordinates. We can then introduce x α , the local relative position vector of the contact point α, by prjecting X α in the system (T 1 , T 2 , N) : We can thus express the local position vector x α as function of the global vector X : where H α is the transition matrix obtained by combining Eqs. (1,2). Similar relations can be determined with respect to contact forces. The local gap vector between two contact points can be derived from the incremental form of Eq.(3):  Then, let's denote the local and global contact force vectors on point α with respectively r α and R α . The relation between contact force vectors expressed in local and global coordinate systems: Here, due to the presence of adhesion on the contact interface, contact reaction r α is composed of the cumulative effects due to both dry contact and the interface adhesion, hence in which we user to denote dry contact reactions, andr contact forces due to interface adhesion.
We now assemble all the N c equations based on Eqs. (3,4,5), we obtain the following geometric and kinematic relations: with

RCC model for adhesion
The use in this work the RCC model to describe the effect of adhesion between contact surfaces. First introduced by Raous, Cangémi and Cocou [18], this model accounts for unilateral contact, friction and adhesion, based on an energy description of the contact interface, involving a free surface energy Ψ written as: and a pseudo-potential of the surface dissipation Φ: In these expressions, β is a scalar parameter measuring the intensity of adhesion. Introduced by Frémond [4], it varies within the range [0, 1]. Specifically, β = 0 represents no adhesion, β = 1 indicates perfect adhesion. Therefore, any β ∈ (0, 1) refers to partial adhesion between contact surfaces. Other parameters in Eq.(9,10) include: C t and C n are parameters characterizing the initial adhesive stiffness when adhesion is complete. w is decohesion energy threshold. I K and I P represent two indicator functions that assure unilateral contact (x n 0), and meaningful values of the degree of adhesion (β ∈ [0, 1]). The subscripts K and P indicate K = {θ | θ 0} and P = {η | 0 η 1}. µ is coefficient of friction, b is surface viscosity and r n is contact force.
Deriving the surface free energy Eq. (9), we obtain expression of the normal force of adhesion and the tangential force of adhesion : Both adhesion forces are dependant of the degree of adhesion β. Then deriving energy functions Eq. (9) and Eq.(10) with respect to β andβ yields the incremental expression of β which gives its evolution in time:

Modified Signorini law with adhesion
We recall the unilateral contact law, also called Signorini law, which for classical non adhesive contact is characterized by conditions of non-penetration and non adhesion. By usingr α n to denote local normal contact force on the point α due to dry contact, and x n the contact distance, we have The first relation eliminates geometric penetration between contact surfaces. The second inequality indicates the absence of adhesion forces between dry contact surfaces once they are separated. For adhesive contact, since contact forces result from both the effects of dry contact and adhesion, the classical conditions of unilateral contact should be modified by considering Eq. (6), hence Here, adhesive forcesr α n are zeros with surfaces in contact. They will appear when contact surfaces start to separate (the second relation), andr α n will tend to maintain the contact surfaces together. By considering Eq.(11), a modified Signorini condition with account for adhesion writes

Modified Coulomb friction rule with adhesion
Classically, friction problems are studied using Coulomb friction model which is characterized by a set of rate-independent slip rules. It describes tangential contact forces as function of normal forces in the context of dry friction: With consideration of adhesion, both tangential and normal contact forces are supplemented by contributions due to interface adhesion as shown in Eq.(6), the above rules become in whichr α t , the adhesive tangential force on contact point α can be calculated by considering Eq. (12): and in the normal direction, contact forces are With the consideration of interface adhesion, tangential friction is made from two contributions. The first follows the classical Coulomb rule and disappears once contact surfaces are separated. The second, r α t , the adhesive tangential force arises when slip occurs, and maintains even with the surface starting to separate.

Complete adhesive contact law within the bipotential method
By combining the modified Signorini law and Coulomb rule, we obtain the complete contact law with the account for interface adhesion as follows : x α n = 0 and x α t = 0, r α =r α Sliding : in whichr α n refers to the normal contact force on point α when surfaces are in contact. In the Sticking situation, since no relative motion occurs, adhesive forces are absent, contact force vectorr α lies in the classical Coulomb cone K µ , defined by The bipotential fuction and inequality of contact law is as follows: where I is indicator function. ℜ− and K u represent respectively the negative real numbers and Coulomb cone.
The indicator functions become null when the variables −x α and r α comply with the restraining conditions, α denotes any contact node (α ∈ 1, ..., N c ).
We multiply both sides of the inequality (24) a parameter ρ, which is used to ensure numerical convergence, and substitude (23) into (24): Taking into account the decomposition x = x t + x n n, the following inequality has to be satisfied: where the modified augmented contact force r α * is defined by: r α is the projection of r α * onto the closed convex Coulomb cone: According to the three different contact states, the projection procedure becomes: i f µ||r α * t || < −r α * n then r α = 0 separating elsei f ||r α * t || ≤ µr α * n then r α = r α * sticking else r α = r α * − ( ||r α * t ||−µr α * Since contact between soft bodies involves treatment of nonlinear kinematic relations, we formulate the nonlinear finite element problem within the framework of large deformations. In this work, we use Green-Lagrangian strain tensor E which comprises both linear and nonlinear terms, as function of nodal displacements u: where B L is the matrix relating the linear strain term to nodal displacements, and B NL (u), relates the nonlinear strain term to nodal displacements. From Eq.(30), the incremental form of the strain-displacement relationship can be written as : Using the principle of virtual displacement, we can write the virtual work δU of the problem as : where the second Piola-Kirchhoff stress tensor S. The vector of contact reaction force R is expressed in the global coordinate system. It is obtained by considering Eqs. (5,6,7) and includes in particular contributions due to adhesion : We can identify in Eq.(34) the vector of internal force: Since δu is arbitrary, a set of nonlinear equations can be obtained as It is noted that the stiffness effect is taken into account by the internal force vector F int . Eq.(36) can be transformed into with the initial conditions at t = 0u =u 0 and u = u 0 .
Taking the derivative of F int with respect to the nodal displacements u gives the tangent stiffness matrix as In addition, by considering Eqs.(31), the tangent stiffness matrix can be written as the sum of the elastic stiffness matrix K e , the geometric stiffness (or initial stress stiffness) matrix K σ and the initial displacement stiffness matrix K u : with

Numerical integration algorithm
Now we need to integrate Eq.(37) between consecutive time configuration t and t + ∆t. The Newmark method is the most common method which is based on a second order algorithm. However, higher order approximation does not necessarily mean better accuracy and may even be redundant in impact problems. When the contact conditions suddenly change (impact, release of contact), the velocity and acceleration are not continuous, and excessive regularity constraints may lead to serious errors. For this reason, Jean [7] has proposed a first order algorithm which is used in this work, Eq.(37) can be transformed into: This algorithm is based on the following approximations: In the iterative solution procedure, all the values at time t + ∆t are replaced by the values of the current iteration i + 1; for example, Finally, we obtain the recursive form of (42) in terms of displacements: where the so-called effective terms are given bȳ At the end of each time step, the velocity is updated bẏ By setting θ = 1 2 , this scheme is then called the implicit trapezoidal rule and it is equivalent to the Tamma -Namburu method in which the acceleration need not be computed [19].

Indentation on a hyperelastic material with adhesive surface
The adhesion effect is usually most significant in the normal direction. In order to clearly show the evolution of β during the complete process of bonding and de-bonding, the first example simulates the normal adhesive contact between an elastic semi-sphere and a hyperelastic block, shown in Figure 2(a). Material parameters of two bodies are: For the sphere indentor, Young's modulus E = 7 · 10 10 Pa, Poisson's ratio ν = 0.3, density ρ = 2500kg/m −3 ; For the block, Young's modulus E = 5 · 10 6 Pa, Poisson's ratio ν = 0.48, density ρ = 1000kg/m −3 . In this case, the sphere indentor behaves as a rigid body compared to the block.
A time dependent displacement is prescribed on the upper surface of the semi-sphere, so that a complete cycle of indentation is performed in 4 seconds. Figure 2 Figure 3(c). Cocou and al. in [2] investigated similar scenarios and obtained concordant results. Figure 3(a) shows the evolution of β on 7 contact nodes in case 1. On any contact point, its horizontal distance from the center point O determines the time length of the contact process on this point, involving bonding and de-bonding. The sequence of β evolution is thus distinctive on each point. Figure 3(b) shows the normal adhesion force of the 7 contact nodes in case 1. Similarly, the distance from the center point O determines the sequence of separation, which however does not influence the adhesion force at the moment of separation. Figure 3(c) shows evolution of β on the contact point A under 3 groups of different adhesive parameters. Figure 3(d) shows the normal adhesion force of the contact point A during the de-bonding process under 3 cases. We can see that the increase of C n makes the detachment difficult. The difference in adhesion force determines the rate of decrease of β as shown in Figure 3(c).

Adhesive friction between a hyperelastic plate and a deformable semi-cylinder
In this example, we investigate the adhesive friction between a hyperelastic plate that slides on top of a deformable semi-cylinder, as shown in Figure 4(a). Both the plate and the semi-cylinder are modelled by Blatz-Ko hyperelastic material, based on the same material property with shear modulus G = 10MPa. Radius of the cylinder is 5mm, and the plate thickness H = 2mm. The plate is sufficiently long so as to ensure contact between the plate and the cylinder during the simulation. While the plate is allowed to slide horizontally, the bottom surface of the semi-cylinder is fixed. The simulated scenario involves 2 L.B. Hu, Y. Cong and Z.-Q. Feng stages. On the first stage, the upper plate descends for 1mm to exert slight compression on the cylinder. Then on the second stage, lateral displacement is prescribed on the plate at a constant velocity. As a result of the combined effect of friction and adhesion, the cylinder is dragged to deform, and we investigate the interfacial behaviour during the process. In particular, by varying the descent velocity of the first stage, we modulate the total time of compression before sliding, during which bonding process takes place. This will have impact on the final adhesion level (characterized by β) before de-bonding starts at the onset of the sliding stage. In order to explore the influence of the adhesion level β on the subsequent adhesive friction behaviour, we set up 5 groups of cases with for each group a different descent velocity ( summarized by Figure 5(a)). Then for each group, we test on 5 different friction coefficients µ, so as to investigate the combined effect of friction and adhesion on the tangential interface behaviour. The adhesive parameters used in the simulations are: w = 20J · m −2 , C n = 2 · 10 9 N · m −3 and b = 0.1Ns · m −1 . Figure 4(b) and (c) present the Von Mises stress distributions of the sliding system, respectively at the onset of sliding, and during the sliding process. We post-process the frictional adhesive behaviour by isolating 2 nodes belonging to the system: as seen in Figure 4 the plate in contact with the cylinder, and one red node on top of the cylinder, in contact with the plate.
Here, we investigate the evolution of β during the first stage. By considering different descent velocities of the plate, varying from 1000mm/s to 200mm/s, we modulate for each case the time for the bonding process. As shown in Figure 5(c), the case with the plate slowly descending at 200mm/s (green curve) had sufficient time to achieve perfect bonding of adhesive links (β reached 1), whereas the most rapid descent (blue curve) did not allow enough time for the formation of complete bonding. In this case, debonding was already initiated after β reached 0.2. We then investigate the influence of friction coefficient µ on the de-bonding behaviour, by prescribing varying friction coefficients µ while considering the same plate descent velocity. We report in Figure 5(b) 5 simulations based on 5 values of µ ranging from 0 to 0.8. All the 5 cases consider the same plate descent velocity of 1000mm/s (Case 1 of Figure 5(a)). Results in Figure 5(b) indicate the formation of stronger bond (higher β) on rougher surfaces (greater µ). This can be interpreted by the fact that a rougher surface (higher µ) delays the onset of sliding motion, according to the Coulomb friction model, which results in longer time for better bonding of adhesive links. Therefore, we observe concordant trend on the curves reflecting tangential adhesion forces. With greater friction coefficient ( Figure 5(d)), the onset of de-bonding is delayed, creating increased level of adhesion force. Then, lower descent velocity on the first stage ( Figure 5(e)) also creates the effect of delaying the onset de de-bonding, permitting better bonding and more significant adhesion forces.

Conclusions
In this work, an interface adhesion model is proposed within the bipotential method to solve contact and friction problems between soft bodies. The model proposes a straightforward description of the interface adhesion based on a local scalar parameter, and allows coupling the effect of adhesion, friction and unilateral contact within a unified framework. Both normal and tangential effects are taken into account by the adhesive interface model. Interface behaviours can be tracked beyond the onset of tangential sliding or normal separation, with interface adhesion driven by the conditions of normal contact. Numerical examples have been performed to investigate the effects of friction and adhesion, including their combined effect, on the interface behaviour based on frictional contact scenarios in tangential and normal directions.