<figure><img src="image_1.png" /></figure>

HAL
open science

# Solution of contact problems with friction

Pierre Joli, Zhi-Qiang Feng, Nadjet Talbi

► To cite this version:

Pierre Joli, Zhi-Qiang Feng, Nadjet Talbi. Solution of contact problems with friction: Pseudo-potential or Bi-potential? Uzawa or Newton?. 8e Colloque national en calcul des structures, CSMA, May 2007, Giens, France. hal-01493724

HAL Id: hal-01493724

https://hal.science/hal-01493724

Submitted on 22 Mar 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
# Solution of contact problems with friction

## Pseudo-potential or Bi-potential ? Uzawa or Newton ?

Pierre Joli* — Zhi-Qiang Feng** — Nadjet Talbi*

*Laboratoire IBISC, Université d'Évry-Val d'Essonne, 40 rue du Pelvoux, 91020 Évry

**Laboratoire de Mécanique d'Évry, Université d'Évry-Val d'Essonne, 40 rue du Pelvoux, 91020 Évry

{pjoli,feng,talbi}@iup.univ-evry.fr

**RÉSUMÉ.** Cet article traite de la modélisation numerique des problèmes de contact avec frottement en élasticité linéaire et en grands glissements. Deux formulations variationnelles appelées respectivement pseudo-potentiels et bi-potentiel sont présentées du point de vue algorithmique. Une comparaison quantitative a été menée pour la première fois entre une resolution numérique de type Uzawa et une autre de type Newton dans le cadre du bi-potentiel. Un exemple test est inclu afin d'évaluer les performances des algorithmes développés.

**ABSTRACT.** This paper is concerned with the numerical modeling of three-dimensional contact problems in elastostatics with Coulomb friction and large slips. Two variational formulations called pseudo-potential and bi-potential are presented from an algorithmic point of view. A quantitative comparison is made for the first time between the Newton algorithm and the Uzawa algorithm in the bi-potential framework. A test example is included to demonstrate the developed algorithms and to highlight their performance.

**MOTS-CLÉS:** Contact avec frottement, Uzawa, Newton, Pseudo-potentiel, Bi-potentiel

**KEYWORDS:** Contact and friction, Uzawa, Newton, Pseudo-potential, Bi-potential
# 1. Introduction

Problems involving contact and friction are among the most difficult ones in mechanics. In the last decade, substantial progress has been made in the analysis of contact problems using finite element procedures (Wriggers, 2002, Laursen, 2002). Heegaard and Curnier (Heegaard et al., 1993) presented a generalized Newton method for large slip frictionless contact problems. From a simple spring-wall contact problem, they have also provided some interesting discussion on convergence properties of the generalized Newton method compared to a global Uzawa algorithm. The bi-potential method proposed by De Saxcé and Feng provides a powerful tool to model dissipative constitutive laws such as Coulomb friction laws (de Saxcé et al., 1998). An iterative Uzawa algorithm can be used to solve the nonlinear implicit equations and this algorithm has been successfully applied by Feng et al. (Feng, 1995, Feng et al., 2003) to simulate large deformation contact problems.

The aim of the present paper is to develop a Newton like algorithm to solve the local contact nonlinear equations within the bi-potential framework. Characteristics of Uzawa and Newton algorithms are discussed. A test numerical example is performed in this study to show the validity of the developed algorithms. The pseudo-potential is only presented from an algorithmic point of view.

# 2. Governing equations

The finite element method is often used in computational mechanics. A typical solution procedure to solve nonlinear problems is the Newton-Raphson iterative procedure :

$$ \left\{ \begin{array}{l} \mathbf{K}_T^i \Delta \mathbf{U} = \mathbf{F}_{\text{int}}^i + \mathbf{F}_{\text{ext}} + \mathbf{R} \\ \mathbf{U}^{i+1} = \mathbf{U}^i + \Delta \mathbf{U} \end{array} \right. \qquad [1] $$

where $\mathbf{K}_T$ is the tangent stiffness matrix, $\mathbf{F}_{\text{int}}$ the vector of internal forces, $\mathbf{F}_{\text{ext}}$ the vector of external loads, $\mathbf{R}$ the vector of contact reaction forces and $\Delta \mathbf{U}$, the vector of nodal displacements correction.

It is noted that Eq.[1] can not be solved directly because $\Delta \mathbf{U}$ and $\mathbf{R}$ are both unknown. The key idea is to determine first the reaction vector $\mathbf{R}$ in a reduced system which only concerns the contact nodes. Then the displacement increments $\Delta \mathbf{U}$ can be computed in the whole structure, using contact reactions as external loading.

For each contact point, we define the reaction vector $\mathbf{r}$ and the gap vector $\mathbf{x}$ in its local coordinate frame ($\mathbf{n}, \mathbf{t}_1, \mathbf{t}_2$) as follows :

$$ \mathbf{r} = \mathbf{H}\mathbf{R}, \quad \mathbf{x} = \mathbf{H}\Delta\mathbf{U} + \mathbf{g} \qquad [2] $$

where $\mathbf{g}$ is the initial gap vector updated at the beginning of each load step and $\mathbf{H}$ is the rotation matrix. Solving a contact problem with $N_c$ contact points consists in having the $N_c$ vector solutions $\chi = (\mathbf{r} \ \mathbf{x})$ satisfy the Signorini conditions and the coulomb friction laws.
## 2.1. Signorini conditions and Coulomb friction laws

The unilateral contact law is characterized by a geometric condition of non-penetration, a static condition of no-adhesion and a mechanical complementary condition. These three conditions are known as Signorini conditions expressed, for each contact point, in terms of the signed contact distance $x_n$ and the normal contact force $r_n$ by

$$ \text{Signor}(x_n, r_n) \Leftrightarrow x_n \ge 0, r_n \ge 0 \quad \text{and} \quad x_n r_n = 0 \quad [3] $$

Classically, a rate independent dry friction law is characterized by a kinematic slip rule. In this work, the classic Coulomb friction rule is used and defined by

$$ \text{Coul}(\dot{\mathbf{x}}_t, \mathbf{r}_t) \Leftrightarrow \text{ if } \|\dot{\mathbf{x}}_t\| = 0 \text{ then } \|\mathbf{r}_t\| \le \mu r_n \text{ else } \mathbf{r}_t = -\mu r_n \frac{\dot{\mathbf{x}}_t}{\|\dot{\mathbf{x}}_t\|} \quad [4] $$

where $\mathbf{r}_t$ is the tangential contact forces, $\dot{\mathbf{x}}_t$ is the relative tangential velocity and $\mu$ is the coefficient of friction.

## 2.2. Variational approaches

The complete contact law (Signorini conditions + Coulomb friction laws) is a complex non smooth dissipative law including three statuses (no contact, contact with sticking and contact with sliding). In the classical approach (Duvaux et al., 1972), the gap vector and contact forces at each contact point are defined by means of two distinct pseudo-potentials: the first for unilateral contact and the second for the friction. In the bi-potential approach proposed by De Saxcé and Feng (de Saxcé et al., 1998), the unilateral contact and friction are coupled. In both cases, the augmented Lagrangian method can be used to regularize nondifferentiable potentials (Klarbring, 1992).

### 2.2.1. Pseudo-potential

We introduce a vector $\mathbf{Z} = (Z_n \quad \mathbf{Z}_t)$, the Signorini conditions are expressed by

$$ Z_n = r_n - \text{Proj}_{\mathbb{R}^+}(r_n^*) = 0 \quad \text{with} \quad r_n^* = r_n - \rho_n x_n \quad \forall \rho_n > 0 \quad [5] $$

The projection operator is defined by

$$ \text{Proj}_{\mathbb{R}^+}(\mathbf{y}) = \mathbf{y} \text{ if } \|\mathbf{y}\| < \mu r_n \text{ else } \text{Proj}_{\mathbb{R}^+}(\mathbf{y}) = 0 \quad [6] $$

Under the condition $r_n > 0$, the friction law is defined as follows

$$ \mathbf{Z}_t = \mathbf{r}_t - \text{Proj}_{\mathbf{D}_\mu}(\mathbf{r}_t^*) = \mathbf{0} \quad \text{with} \quad \mathbf{r}_t^* = \mathbf{r}_t - \rho_t \mathbf{x}_t \quad \forall \rho_t > 0 \quad [7] $$

The projection operator is defined by

$$ \text{Proj}_{\mathbf{D}_\mu}(\mathbf{y}) = \mathbf{y} \text{ if } \|\mathbf{y}\| < \mu r_n \text{ else } \text{Proj}_{\mathbf{D}_\mu}(\mathbf{y}) = -\mu r_n \frac{\mathbf{x}_t}{\|\mathbf{x}_t\|} \quad [8] $$
### 2.2.2. Bi-potential

In the bi-potential approach, only one projection operator is required to satisfy unilateral contact and friction laws :

$$ Z = \mathbf{r} - \mathrm{Proj}_{\mathbf{K}_{\mu}}(\mathbf{r}^{*}) = \mathbf{0} \quad \text{with} \quad \mathbf{r}^{*} = \mathbf{r} - \rho \mathbf{x}^{*} \quad \text{and} \quad \mathbf{x}^{*} = \mathbf{x} + \mu \|\mathbf{x}_{t}\| \mathbf{n} \quad [9] $$

The projection operation can be explicitly defined by

$$ \begin{align} \mathrm{Proj}_{\mathbf{K}_{\mu}}(\mathbf{r}^{*}) &= \mathbf{r}^{*} \quad \text{if} \quad \|\mathbf{r}_{t}^{*}\| < \mu r_{n}^{*}, \nonumber \\ \mathrm{Proj}_{\mathbf{K}_{\mu}}(\mathbf{r}^{*}) &= \mathbf{0} \quad \text{if} \quad \mu \|\mathbf{r}_{t}^{*}\| < -r_{n}^{*}, \tag{10} \\ \mathrm{Proj}_{\mathbf{K}_{\mu}}(\mathbf{r}^{*}) &= \mathbf{r}^{*} - \left( \frac{\|\mathbf{r}_{t}^{*}\| - \mu r_{n}^{*}}{1 + \mu^{2}} \right) \left( \frac{\mathbf{r}_{t}^{*}}{\|\mathbf{r}_{t}^{*}\|} - \mu \mathbf{n} \right) \quad \text{otherwise.} \nonumber \end{align} $$

## 3. Local solution methods

The algorithmic principle to solve a multi-contact problem consists in decomposing the global solution into $N_c$ successive local solutions of the following equations defined at each contact point :

$$ \mathbf{f}(\chi) = \begin{Bmatrix} \mathbf{x} - \mathbf{W}\mathbf{r} - \tilde{\mathbf{x}} \\ \mathbf{Z} \end{Bmatrix} = \mathbf{0} \qquad [11] $$

where **W** represents the local flexibility matrix (3 × 3) at the contact point and $\tilde{\mathbf{x}}$ is the part of the gap due to the initial gap, the external forces and the contribution of the other contact forces. This contribution is “frozen” when the local solution is computed. Then following the classical Gauss-Seidel approach, the contribution of the other contact forces is updated for the next contact point.

### 3.1. Uzawa algorithm

A numerical solution of the implicit equation [11] can be carried out by means of Uzawa's algorithm which leads to an iterative process involving one predictor-corrector step in case of bi-potential :

$$ \begin{align} \text{Predictor } \mathbf{r}^{*(k+1)} &= \mathbf{r}^{(k)} - \rho^{(k)} (\mathbf{x}^{*(k)}) \nonumber \\ \text{Corrector } \mathbf{r}^{(k+1)} &= \mathrm{Proj}_{\mathbf{K}_{\mu}} (\mathbf{r}^{*(k+1)}) \tag{12} \end{align} $$

where $k$ and $k+1$ are the iteration numbers at which the contact reactions are computed. In case of pseudo-potential, two predictor-corrector steps are necessary : the first one for computing normal contact force by Eq.[5] and the second one for computing tangential contact forces by Eq.[7].
## 3.2. Newton algorithm

The Newton algorithm applied to the implicit equation [11] leads to the following iterative numerical scheme :

$$ \left[ \frac{\partial f(\chi^i)}{\partial \chi} \right] \Delta \chi = -f(\chi^i) \\ \chi^{i+1} = \chi^i + \Delta \chi \qquad [13] $$

## 4. A test example

A test example with a linearly elastic law is considered in this study. The example concerns the contact between a 3D elastic block ABCDEFGH and a rigid surface Γ (Figure 1). The upper surface ABCD is given a rigid motion described by (a, b, θ)

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Contact problem under displacement control</figcaption></figure>

where **a** and **b** are respectively perpendicular and parallel to the rigid surface. The loading program is designed to apply first a vertical displacement following **a** and then a horizontal displacement following **b**. The lower surface EFGH comes into contact with the rigid surface whose normal vector is (0, 0, 1). Each side of the block has a length of 1 mm. The other characteristics of this example are : Young's modulus : $E = 210000 \text{ N/mm}^2$, Poisson's ratio : $\nu = 0.3$, coefficient of friction : $\mu = 0.3$, boundary conditions : $\|\mathbf{a}\| = 0.1 \text{ mm}$, $\|\mathbf{b}\| = 0.4 \text{ mm}$, $\theta = 60^{\circ}$. Fifty load steps are performed for this problem, so a horizontal or vertical displacement of 0.01 mm is applied to ABCD each step. Figure 2 shows the variation of normal contact forces of points F and H. The horizontal lines correspond to cases where the block slides in a stationary regime. The analysis was performed by means of both Uzawa and Newton algorithms developed above. The results are almost the same as shown in Figure 2. However, the number of iterations in the global contact solution is different, as expected. Figure 3 indicates the evolution of iterations with respect to cumulative iterations in the Newton-Raphson procedure. It shows that, globally, the Uzawa algorithm needs more iterations than the Newton algorithm.
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Variation of normal contact forces of points F and H</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Evolution of contact iterations</figcaption></figure>

**5. Conclusions**

In this work, three-dimensional contact problems with friction have been theoretically investigated and numerically implemented. Within the bi-potential framework, a Newton algorithm has been proposed. A comparative study has been made between the newly proposed Newton algorithm and the previously developed Uzawa algorithm. The characteristics of each algorithm are discussed. The numerical test indicates that both algorithms give good results. The numerical test shows also that the Uzawa algorithm needs more iterations than the Newton algorithm. But the Newton algorithm needs an additional local iterative procedure and, at each iteration, a small system of equations should be solved. Further investigation on the pseudo-potential is being undertaken.

**6. Bibliographie**

de Saxcé G., Feng Z.-Q., « The bi-potential method : a constructive approach to design the complete contact law with friction and improved numerical algorithms », *Mathematical and Computer Modeling*, vol. 28(4-8), p. 225-245, 1998.

Duvaux G., Lions J.-L., *Les inéquations en mécanique et en physique*, Dunod, Paris, 1972.

Feng Z.-Q., « 2D or 3D frictional contact algorithms and applications in a large deformation context », *Comm. Numer. Meth. Engng.*, vol. 11, p. 409-416, 1995.

Feng Z.-Q., Peyraut F., Labed N., « Solution of large deformation contact problems with friction between Blatz-Ko hyperelastic bodies », *Int. J. Engng. Science*, vol. 41, p. 2213-2225, 2003.

Heegaard J. H., Curnier A., « An augmented Lagrangian method for discrete large slip contact problems », *Int. J. Numer. Meth. Engng.*, vol. 36, p. 569-593, 1993.

Klarbring A., « Mathematical programming and augmented Lagrangian methods for frictional contact problems », in, A. Curnier (ed.), *Contact Mechanics Int. Symp.*, 1992. PPUR.

Laursen T. A., *Computational contact and impact mechanics : Fundamentals of Modeling Interfacial Phenomena in Nonlinear Finite Element Analysis*, Springer Verlag, 2002.

Wriggers P., *Computational contact mechanics*, John Wiley & Sons, 2002.