# Un principe mixte pour le problème de contact avec frottement de Coulomb

Anh Le Van, Tai Nguyen

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

## To cite this version:

Anh Le Van, Tai Nguyen. Un principe mixte pour le problème de contact avec frottement de Coulomb. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01408462

HAL Id: hal-01408462

https://hal.science/hal-01408462v1

Submitted on 5 Dec 2016

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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.

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

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Un principe mixte pour le problème de contact avec frottement de Coulomb

A. Le van¹, T.H.T. Nguyen¹

¹ GeM (Institut de Recherche en Génie Civil et Mécanique)
Université de Nantes, 2, rue de la Houssinière
BP 92208, 44322 Nantes Cedex 3, France
{anh.le-van,tai.nguyen}@univ-nantes.fr

**Résumé** — Nous proposons dans ce travail une relation mixte de type lagrangien augmenté en vue de résoudre le problème de contact en grandes transformations avec frottement de Coulomb. La relation est discrétisée par la méthode des éléments finis et quelques exemples numériques sont présentés permettant de montrer l'efficacité de la formulation proposée.

**Mots clés** — Contact, Frottement, Coulomb.

## 1 Formulation théorique

Considérons deux corps, repérés par les indices 1 et 2, subissant les transformations $\phi^{(1)}$ et $\phi^{(2)}$ dans l'espace 3D sur un certain intervalle de temps $[O, T]$. Nous supposons que les deux corps peuvent entrer en contact et nous formulons le problème de contact selon le point de vue lagrangien. Les notations, classiques, sont résumées comme suit. Dans la configuration de référence, le corps $i$ occupe le domaine $\Omega_o^{(i)}$. Les forces massiques dans le corps $i$ sont notées $f_i^{(i)}$. La frontière $S_o^{(i)}$ du corps $i$ est partitionnée en trois parties notées $S_{oU}^{(i)}$, $S_{oT}^{(i)}$ et $S_{oc}^{(i)}$, où $S_{oU}^{(i)}$ et $S_{oT}^{(i)}$ sont respectivement les parties où les déplacements et les contraintes sont imposés, et $S_{oc}^{(i)}$ est la partie où le contact a lieu. L'état des contraintes dans le corps $i$ est défini par le premier tenseur des contraintes de Piola-Kirchhoff $\Pi^{(i)}$. Le vecteur de contrainte nominale en tout point de $S_o^{(i)}$ et relatif à la normale $N^{(i)}$ est noté $T^{(i)} = \Pi^{(i)} \cdot N^{(i)}$. La position actuelle de la surface $S_o^{(i)}$ est notée $S^{(i)}$.

Pour chaque point $\mathbf{x} \in S_c^{(1)}$, on définit le point de contact (potentiel) $\mathbf{y} \in S_c^{(2)}$ comme étant le point le plus proche de $\mathbf{x}$ via la relation $\mathbf{y} = \arg\min_{\mathbf{x}^{(2)} \in S_c^{(2)}} \|\mathbf{x} - \mathbf{x}^{(2)}\|$. On définit aussi la proximité $g = -v(\mathbf{x} - \mathbf{y})$ où $v$ est la normale extérieure au point $\mathbf{y}$. Le point de référence $\mathbf{X} \in S_{oc}^{(1)}$ correspondant au point $\mathbf{x}$ en question vérifie $\mathbf{x} = \phi^{(1)}(\mathbf{X})$ et le point de référence $\mathbf{Y}(\mathbf{X}) \in S_{oc}^{(2)}$ correspondant au point $\mathbf{y}$ vérifie $\mathbf{y} = \phi^{(2)}(\mathbf{Y})$. La cinématique tangentielle est caractérisée par la vitesse de glissement $\mathbf{V}_T = \dot{\eta}^\alpha \mathbf{a}_\alpha$ qui est décomposée dans la base locale actuelle ($\mathbf{a}_1, \mathbf{a}_2, \mathbf{v}$) au point $\mathbf{y} \in S_c^{(2)}$. De même, le vecteur de contrainte nominale en tout point $\mathbf{X} \in S_{oc}^{(1)}$ est décomposé comme $\mathbf{T} = T_N\mathbf{v} - \mathbf{T}_T$. Le coefficient de frottement entre les deux corps est noté $\mu$.

La forme faible proposée ici se présente comme un principe mixte faisant intervenir à la fois les champs des déplacements $\mathbf{U}^{(i)}$, $i \in \{1, 2\}$, définis dans $\Omega_o^{(i)}$ et le champ des multiplicateurs $\lambda_N$ et $\lambda_T$ défini sur $S_{oc}^{(1)}$. Ainsi, les fonctions tests sont les déplacements virtuels $\mathbf{U}^{(1)*}$, $\mathbf{U}^{(2)*}$, et les multiplicateurs virtuels $\lambda_N^*$, $\lambda_T^*$. On suppose que les champs virtuels sont des champs arbitraires suffisamment réguliers pour qu'on puisse les intégrer ou les différentier à l'ordre souhaité. Ayant choisi deux scalaires positifs $\varepsilon_N, \varepsilon_T$ fixés, la forme faible du problème s'énonce [1] :
## PROPOSITION.

► Les champs solutions du problème fort – à savoir $\mathbf{U}^{(1)}, \mathbf{U}^{(2)}$ dans $\Omega_o^{(1)}, \Omega_o^{(2)}$ et $T_N, T_T, V_T$ sur $S_{oc}^{(1)}$ vérifient la relation suivante, où on fait $\lambda_N = T_N$ et $\lambda_T = T_T$ :

$$ \forall t, \forall \mathbf{U}^{(1)*}, \forall \mathbf{U}^{(2)*}, \forall \lambda_N^*, \forall \lambda_T^*, $$

$$ \begin{aligned} & \sum_{i=1}^{2} \left\{ -\int_{\Omega_o^{(i)}} \Pi^{(i)T} : \nabla_{\mathbf{X}^{(i)}} \mathbf{U}^{(i)*} d\Omega_o + \int_{\Omega_o^{(i)}} \rho_o^{(i)} \mathbf{f}^{(i)} \mathbf{U}^{(i)*} d\Omega_o + \int_{S_{oc}^{(i)} \cup S_{oT}^{(i)}} \mathbf{T}^{(i)} \mathbf{U}^{(i)*} dS_o \right\} \\ & + \int_{S_{oc}^{(i)}} \left[ \langle \hat{\lambda}_N \rangle \mathbf{v} - \left( 1 - \langle 1 - \frac{\mu \langle \hat{\lambda}_N \rangle}{\| \hat{\lambda}_T \|} \rangle \right) \hat{\lambda}_T \right] \left( \mathbf{U}^{(1)*}(\mathbf{X}) - \mathbf{U}^{(2)*}(\mathbf{Y}(\mathbf{X})) \right) dS_o \\ & + \int_{S_{oc}^{(i)}} \left\{ (\lambda_N - \langle \hat{\lambda}_N \rangle) \frac{\lambda_N^*}{\varepsilon_N} + \left[ \lambda_T - \left( 1 - \langle 1 - \frac{\mu \langle \hat{\lambda}_N \rangle}{\| \hat{\lambda}_T \|} \rangle \right) \hat{\lambda}_T \right] \frac{\lambda_T^*}{\varepsilon_T} \right\} dS_o = \sum_{i=1}^{2} \int_{\Omega_o^{(i)}} \rho_o^{(i)} \dot{\mathbf{U}}^{(i)} \mathbf{U}^{(i)*} d\Omega_o \end{aligned} \quad (1) $$

où on a noté pour abréger $\hat{\lambda}_N = \lambda_N + \varepsilon_N g$, $\hat{\lambda}_T = \lambda_T + \varepsilon_T V_T$, $\langle \cdot, \cdot \rangle$ désigne la fonction "partie positive" : $\langle a \rangle = a$ si $a \ge 0$, $= 0$ si $a < 0$, et où le terme $\langle 1 - \frac{\mu \langle \hat{\lambda}_N \rangle}{\| \hat{\lambda}_T \|} \rangle$ doit être remplacé par 0 en tout point de $S_{oc}^{(1)}$ où $\hat{\lambda}_T = \mathbf{0}$.

► Inversement, les surfaces $S_c^{(1)}$ et $S_c^{(2)}$ étant paramétrisées par (un domaine de $\mathbb{R}^2) \ni \xi \mapsto \mathbf{x} \in S_c^{(1)} \subset \mathbb{R}^3$ et (un domaine de $\mathbb{R}^2) \ni \eta \mapsto \mathbf{y} \in S_c^{(2)} \subset \mathbb{R}^3$, supposons qu'il existe une bijection (au moins par morceaux) $S_c^{(1)} \ni \mathbf{x} \mapsto \mathbf{y} \in S_c^{(2)}$ entre les deux surfaces de contact. Alors, la relation (1) implique à tout instant $t$ les équations locales suivantes :

(1) L'équation locale de la dynamique dans les corps : $\forall i \in \{1, 2\}, \forall \mathbf{X}^{(i)} \in \Omega_o^{(i)}, \mathrm{div} \Pi^{(i)} + \rho_o^{(i)} \mathbf{f}^{(i)} = \rho_o^{(i)} \ddot{\mathbf{U}}^{(i)}$

(2) La relation suivante sur la portion de frontière $S_o^{(i)} \setminus S_{oc}^{(i)} = S_{oT}^{(i)} \cup S_{oU}^{(i)} : \Pi^{(i)} \mathbf{N}^{(i)} = \mathbf{T}^{(i)}$.

(3) Les égalités entre les contraintes nominales et les multiplicateurs sur la surface contactrice $S_{oc}^{(1)} \forall \mathbf{X} \in S_{oc}^{(1)}, \quad \mathbf{T}^{(1)} = \Pi^{(1)} \mathbf{N}^{(1)} = \lambda_N \mathbf{v} - \lambda_T \quad \Leftrightarrow \quad T_N = \mathbf{T}^{(1)} \mathbf{v} = \lambda_N \ge 0 \text{ et } \mathbf{T}_T = \lambda_T$

(4) Les lois de contact normal et tangentiel : $\forall t, \forall \mathbf{X} \in S_{oc}^{(1)}$,

(a) $g(\mathbf{X}, t) \le 0 \quad \begin{cases} \cdot \text{ si } g(\mathbf{X}, t) < 0, \text{ alors } T_N(\mathbf{X}, t) = 0 \\ \cdot \text{ si } g(\mathbf{X}, t) = 0, \text{ alors } T_N(\mathbf{X}, t) \ge 0 \end{cases}$

(b) $g(\mathbf{X}, t) < 0 \Rightarrow \mathbf{T}_T(\mathbf{X}, t) = \mathbf{0}$, pas de relation sur $\mathbf{V}_T(\mathbf{X}, t)$

(c) $g(\mathbf{X}, t) = 0 \Rightarrow \|\mathbf{T}_T(\mathbf{X}, t)\| \le \mu T_N(\mathbf{X}, t) \quad \begin{cases} \cdot \text{ si } \|\mathbf{T}_T\| < \mu T_N, \mathbf{V}_T = \mathbf{0} \\ \cdot \text{ si } \|\mathbf{T}_T\| = \mu T_N, \mathbf{V}_T \wedge \mathbf{T}_T = \mathbf{0} \text{ et } \mathbf{V}_T \cdot \mathbf{T}_T \ge 0 \end{cases}$

(5) La relation suivante qui traduit l'équilibre des contraintes de contact à l'interface de contact

$$ \forall \mathbf{X} \in S_{oc}^{(1)}, \quad \mathbf{T}^{(2)}(\mathbf{Y}(\mathbf{X},t))dS_{oc}^{(2)} = -(\lambda_N \mathbf{v} - \lambda_T) \frac{\|\mathbf{X}_{\xi1} \wedge \mathbf{X}_{\xi2}\|}{\|\mathbf{Y}_{\eta1} \wedge \mathbf{Y}_{\eta2}\|} \frac{D\xi}{D\eta} dS_{oc}^{(2)} = -\mathbf{T}^{(1)}(\mathbf{X})dS_{oc}^{(1)} \quad (2) $$

## 2 Equation semi-discrète

La relation mixte (1) est discrétisée en espace par la méthode des éléments finis. La surface contactrice $S_{oc}^{(1)}$ est maillée en éléments finis $e^{(1)}$ et la surface cible $S_{oc}^{(2)}$ en éléments finis $e^{(2)}$. Dans la suite, on utilisera la convention de notation où les accolades $\{ \}$ désignent un vecteur-colonne, les crochets $\langle \rangle$ un vecteur-ligne, et les crochets droits $[]$ une matrice quelconque.

On arrive en fin de compte au système semi-discrétisé non linéaire couplé suivant

$$ \begin{cases} [\mathbf{M}]\{\ddot{\mathbf{U}}\} + \{\Psi(\mathbf{U})\} = \{\Phi(\mathbf{U})\} + \{\Phi_{\mathrm{contact}}(\mathbf{U}, \Lambda)\} \\ \{\mathbf{R}_{\Lambda}(\mathbf{U}, \Lambda)\} = \{\mathbf{0}\} \end{cases} \quad (3) $$

où les vecteurs $\{\Phi_{\mathrm{contact}}(\mathbf{U}, \Lambda)\}$ et $\{\mathbf{R}_{\Lambda}(\mathbf{U}, \Lambda)\}$ s'obtiennent par assemblage des vecteurs élémentaires $\{\Phi_{\mathrm{contact}}\}^e j(\xi)$ et $\{\mathbf{R}_{\Lambda}\}^e j(\xi)$, $j(\xi)$ désignant le produit du jacobien d'aire et le poids du
point d'intégration $\xi \in \mathbb{R}^2$ antécédent de $X \in e^{(1)}$. Les vecteurs élémentaires $\{\Phi_{contact}\}^e$ et $\{\mathbf{R}_{\Lambda}\}^e$ sont quant à eux définis par

$$ \{\Phi_{contact}\}^e = \left\{ \begin{aligned} & [\mathcal{N}^{(1)}(\xi)]^{e(1)T} (\langle \hat{\lambda}_N \rangle \{v\} - (1 - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) \{\hat{\lambda}_T\}) \\ & - [\mathcal{N}^{(2)}(\eta)]^{e(2)T} (\langle \hat{\lambda}_N \rangle \{v\} - (1 - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) \{\hat{\lambda}_T\}) \end{aligned} \right\} \quad (4) \\ \{\mathbf{R}_{\Lambda}\}^e = [\mathcal{N}^{(1)}(\xi)]^{e(1)T} \left( \frac{1}{\varepsilon_N} (\langle \hat{\lambda}_N \rangle - \lambda_N) \{v\} + \frac{1}{\varepsilon_T} \{\lambda_T - (1 - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) \hat{\lambda}_T \} \right) $$

Les inconnues de (3) sont les déplacements nodaux $\mathbf{U} = (\mathbf{U}^{(1)}, \mathbf{U}^{(2)})$ dans les corps 1 et 2, ainsi que les multiplicateurs nodaux $\{\Lambda\}$ définis sur la surface de référence de contact $S_{oc}^{(1)}$. A part les vecteurs $\{\Phi_{contact}(\mathbf{U}, \Lambda)\}$ et $\{\mathbf{R}_{\Lambda}(\mathbf{U}, \Lambda)\}$ spécifiques au contact, les autres termes de (3) sont classiques : $[\mathbf{M}]$ est la matrice masse, $\{\Psi(\mathbf{U})\}$ le vecteur force interne dû à la déformation des corps, $\{\Phi(\mathbf{U})\}$ le vecteur force externe (force de contact excluse) qui peut dépendre du déplacement $\{\mathbf{U}\}$ s'il y a des forces suiveuses de type pression appliquées sur l'un des corps. Le système non linéaire semi-discrète (3) sera résolu par la méthode de Newton-Raphson afin d'obtenir simultanément les déplacements et les multiplicateurs nodaux, de manière analogue à ce qui est fait dans [2].

Une quantité importante en contact qu'il faut discrétiser de manière adéquate est la vitesse de glissement $V_T = \dot{\eta}^\alpha \mathbf{a}_\alpha$ contenue dans (4). Notons $t$ l'instant actuel correspondant à l'incrément courant $n$ et $t_{n-1}$ l'instant correspondant à l'incrément précédent, et mettons-nous à une itération courante de l'incrément de temps courant $n$. Soit une particule $\mathbf{X}$ de la surface contactrice $S_{oc}^{(1)}$, sa position à l'instant $t_{n-1}$ est $\mathbf{x}_{n-1} = \phi^{(1)}(\mathbf{X}, t_{n-1})$ et sa position actuelle est $\mathbf{x} = \phi^{(1)}(\mathbf{X}, t)$. Définissons $\mathbf{y}_{n-1}$ comme la projection de $\mathbf{x}_{n-1}$ sur la surface cible $S_c^{(2)}(t_{n-1})$ à $t_{n-1}$, et $\mathbf{Y}_{n-1}$ l'antécédent de $\mathbf{y}_{n-1}$. En appliquant le schéma d'Euler implicite au pas de temps $\Delta t = t - t_{n-1}$, on arrive à l'expression discrète suivante de la vitesse de glissement

$$ \mathbf{V}_T \Delta t = \mathbf{y} - \Phi^{(2)}(\mathbf{Y}_{n-1}, t) \quad (5) $$

### 3 Matrice tangente de contact

On se bornera ici à la partie de la matrice tangente qui est relative au contact, étant donné que la partie restante due aux termes hors contact est bien connue. La matrice tangente de contact est obtenue par assemblage des matrices élémentaires $[\mathbf{K}_{contact}]^e$ de la forme

$$ [\mathbf{K}_{contact}]^e = \left[ \begin{array}{c|c} [\mathbf{K}_{UU}]^e & [\mathbf{K}_{UA}]^e \\ \hline [\mathbf{K}_{\Lambda U}]^e & [\mathbf{K}_{\Lambda \Lambda}]^e \end{array} \right] \quad (6) $$

où les blocs $[\mathbf{K}_{UU}]^e$, $[\mathbf{K}_{UA}]^e$, $[\mathbf{K}_{\Lambda U}]^e$ et $[\mathbf{K}_{\Lambda \Lambda}]^e$ sont listés dans la suite :

$$ [\mathbf{K}_{UU}]^e = H(\hat{\lambda}_N) \varepsilon_N [\mathbf{k}_{contact}^1]^e + \langle \hat{\lambda}_N \rangle [\mathbf{k}_{contact}^2]^e - H(1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|}) \mu H(\hat{\lambda}_N) \varepsilon_N [\mathbf{k}_{contact}^3]^e \\ + H(1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|}) \mu\langle\hat{\lambda}_N\rangle \varepsilon_T [\mathbf{k}_{contact}^4]^e + (1 - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) (\varepsilon_T [\mathbf{k}_{contact}^5]^e + [\mathbf{k}_{contact}^6]^e) \quad (7) $$

$$ [\mathbf{K}_{\Lambda U}]^e = \left[ \begin{array}{c} [\mathcal{N}^{(1)}]^e T \\ -[\mathcal{N}^{(2)}]^e T \end{array} \right] \left[ \begin{array}{c} -H(\hat{\lambda}_N) \{v\} \langle v \rangle + H(1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|}) \mu H(\hat{\lambda}_N) \frac{\{\hat{\lambda}_T\}}{\|\hat{\lambda}_T\|} \langle v \rangle \\ +H(1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|}) \mu\langle\hat{\lambda}_N\rangle \frac{\{\hat{\lambda}_T\}}{\|\hat{\lambda}_T\|} \langle v \rangle - (1 - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) [\mathbf{I} - v \otimes v] \end{array} \right] [\mathcal{N}^{(1)}]^e T \quad (8) $$

$$ [\mathbf{K}_{\Lambda U}]^e = [\mathcal{N}^{(1)}]^e T \\ \left[ \begin{aligned} & -H(\hat{\lambda}_N) \{v\} \langle v \rangle \langle \mathcal{N}^{(1)} \rangle^{e(1)}, -\langle \mathcal{N}^{(2)} \rangle^{e(2)} \\ & +H(1-\frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|})\mu H(\hat{\lambda}_N)\frac{\{\hat{\lambda}_T\}}{\|\hat{\lambda}_T\|} \langle v \rangle \langle \mathcal{N}^{(1)} \rangle^{e(1)}, -\langle \mathcal{N}^{(2)} \rangle^{e(2)} \\ & +H(1-\frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|})\mu(\hat{\lambda}_N)\frac{\{\hat{\lambda}_T\}}{\|\hat{\lambda}_T\|} \langle [\mathbf{W}] [\mathcal{N}^{(1)}]^{e(1)}, [\mathbf{I}-\mathbf{W}] [\mathcal{N}^{(2)}]^{e(2)} + [\mathbf{Z}^\alpha] [\mathcal{N}_\alpha^{(2)}]^{e(2)}, -\langle \mathcal{N}^{(2)}(\eta_{n-1}) \rangle^{e(2)} \\ & - (1-\langle 1-\frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle) [\mathbf{I}-\mathbf{v}\otimes\mathbf{v}] \langle [\mathbf{W}] [\mathcal{N}^{(1)}]^{e(1)}, [\mathbf{I}-\mathbf{W}] [\mathcal{N}^{(2)}]^{e(2)} + [\mathbf{Z}^\alpha] [\mathcal{N}_\alpha^{(2)}]^{e(2)}, -\langle \mathcal{N}^{(2)}(\eta_{n-1}) \rangle^{e(2)} \end{aligned} \right] \quad (9) $$
$$[\mathbf{K}_{\Lambda\Lambda}]^e = [\mathcal{N}^{(1)}]^{e^{(1)T}} \left[ \frac{1}{\epsilon_N} (H(\hat{\lambda}_N) - 1) \{v\} \langle v \rangle - H(1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|}) \mu H(\hat{\lambda}_N) \frac{1}{\epsilon_T} \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \langle v \rangle \right. \\ \left. - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \frac{\mu}{\epsilon_T} \langle\hat{\lambda}_N\rangle \frac{\{\hat{\lambda}_T\}\langle\hat{\lambda}_T\rangle}{\|\hat{\lambda}_T\|^3} - \langle 1 - \frac{\mu\langle\hat{\lambda}_N\rangle}{\|\hat{\lambda}_T\|} \rangle \frac{1}{\epsilon_T} [\mathbf{I} - v \otimes v] \right] [\mathcal{N}^{(1)}]^{e^{(1)}} \quad (10)$$

Dans la relation (7), les matrices $[\mathbf{k}_{contact}^1]^e$ à $[\mathbf{k}_{contact}^6]^e$ ont des expressions simples dépendant de la géométrie et des multiplicateurs, elles ne sont explicitées ici faute de place. Le symbole $H$ désigne la fonction de Heaviside. Dans la relation (9), $e_{n-1}^{(2)}$ désigne l'élément de la surface cible de référence $S_{oc}^{(2)}$ contenant le point $\mathbf{Y}_{n-1}$, qui est en général différent de l'élément $e^{(2)}$. Les tenseurs du second ordre $\mathbf{W}$ et $\mathbf{Z}^\alpha$, $\alpha \in \{1, 2\}$, non symétriques, sont définis par $\mathbf{W} = \bar{a}^{\alpha\beta}\mathbf{a}_\alpha \otimes \mathbf{a}_\beta$, $\mathbf{Z}^\alpha = g\bar{a}^{\alpha\beta}\mathbf{a}_\beta \otimes \mathbf{v}$, où la matrice $[\bar{a}^\cdot]$ est l'inverse de $[\bar{a}^\cdot]$, avec $\bar{a}^{\alpha\beta} = a^{\alpha\beta} + g\kappa^{\alpha\beta}$, $\kappa^{\alpha\beta} = v^\phi, \phi^{(2)}$. Enfin, malgré l'exposant $e$ rappelant l'élément, $[\mathbf{K}_{contact}]^e$ est en fait défini pour chaque point d'intégration de l'élément $e^{(1)} \subset S_{oc}^{(1)}$.

# 4 Exemples numériques

Les résultats numériques présentés dans la suite ont été obtenus à l'aide d'un code maison et comportent un minimum d'informations par manque de place.

## 4.1 Contact de deux poutres

Le premier exemple concerne le contact quasi-statique de deux poutres hyperélastiques, dont l'une est soumise à une force suiveuse $P$ à une extrémité. Les poutres ont la même géométrie de référence, longueur $L = 1m$ selon $x$, section droite rectangulaire de hauteur $0.05m$ selon $y$ et d'épaisseur $0.1m$ along $z$. Le matériau est hyperélastique de Saint-Venant Kirchhoff avec $E = 2.10^6 N/m^2$ et $\nu = 0.3$. Les poutres sont encastrées à $X = 0$. Dans la configuration de référence, la surface inférieure de la poutre supérieure est confondue avec celle supérieure de la poutre inférieure. Les résultats sont présentés dans la figure 1 pour deux cas : sans frottement et avec frottement $\mu = 0.3$. On remarquera que, dans le cas sans frottement et lorsque les poutres sont bien courbées, il y a perte de contact sur une grande partie de la surface de contact initiale.

<figure><img src="image_3.png" /><figcaption>Figure 1 – Contact entre deux poutres. Déformées et contraintes de contact au dernier incrément : (a) cas $\mu = 0$, (b) cas $\mu = 0.3$. (c) Flèche en bout $v$ en fonction de la force suiveuse $P$.</figcaption></figure>

## 4.2 Contact de deux membranes

Considérons deux membranes faites d'un matériau hyperélastique de type néo-hookéen, dont le potentiel s'écrit $w(\mathbf{C}) = \frac{\mu}{2}(I_1 - 3) - \frac{\mu}{2}\ln I_3 + \frac{\lambda}{8}(\ln I_3)^2$, dépendant des invariants $I_1 = \det \mathbf{C}$ et $I_3 = \det \mathbf{C}$ du tenseur de Cauchy droit $\mathbf{C}$. On prend $E = 2.10^6 N/m^2$, $\nu = 0.3$, où $E$ et $\nu$ sont reliés aux coefficients $\lambda$ et $\mu$ par les relations classiques de l'élasticité. Le problème est modélisé en déformation plane dans le plan $x-y$ et en contrainte plane dans l'épaisseur des membranes. A l'état
non déformé, les membranes ont une épaisseur de 1mm et occupent deux segments parallèles, chacun de 1m de longueur suivant x, distants de 1m suivant y et décentrés de 1m suivant x. Les bords parallèles à l'axe z sont fixes. Les membranes sont soumises à la même pression p qui les amène au contact. Les résultats sont présentés dans la figure 2.

<figure><img src="image_4.png" /><figcaption>Figure 2 – Contact de 2 membranes. Déformées et contraintes de contact : (a) $\mu = 0, p = 4023N/m^2$, (b) $\mu = 0.5, p = 5420N/m^2$. (c) Flèche d'un point en fonction de la pression.</figcaption></figure>

## 4.3 Impact de deux barres

On étudie ici l'impact de deux barres élastiques identiques et alignées, l'une lancée à la vitesse uniforme de 0.1m/s et l'autre immobile à l'instant initial. Les données géométriques et de matériau sont celles classiquement trouvées dans la littérature, à savoir : chaque barre est de longueur 10m et section droite 1m², le matériau est de Saint-Venant Kirchhoff avec $E = 100Pa$, $v = 0$ et la masse volumique est $\rho_o = 0.01kg/m^3$. Il n'y a pas de frottement. L'intégration temporelle est effectuée à l'aide du schéma de Moreau [3], avec un pas de temps $\Delta t = 0.005s$. Les résultats sont présentés dans la figure 3. Les déplacements trouvés concordent avec la solution analytique et l'énergie totale est bien conservée.

<figure><img src="image_5.png" /><figcaption>Figure 3 – Impact de 2 barres. (a) Déplacements des points de contact. (b) Vitesse. (c) Evolution des énergies.</figcaption></figure>

## 4.4 Disque impactant un plan

Considérons un disque hyperélastique impactant un plan rigide (figure 4). Le disque est de rayon 0.1m et d'épaisseur 0.02m. Le matériau est de Saint-Venant Kirchhoff avec $E = 1.333 \times 10^6 Pa$, $v = 0.333$ et la masse volumique est $\rho_o = 6000kg/m^3$. Avant l'impact, le disque est en translation uniforme de vitesse $V = 200km/h$ à 45° et en rotation uniforme. Les calculs sont effectués avec trois valeurs de la vitesse de rotation : $\omega = 0, \omega = 200$ tours/min (effet slicé) et $\omega = -200$ tours/min (effet lifté). Le coefficient de frottement entre le disque et le plan est $\mu = 0.3$. L'intégration temporelle est effectuée à l'aide du schéma de Moreau comme dans l'exemple précédent, avec un pas de temps $\Delta t = 0.001s$. Les résultats sont présentés dans la figure 4. On constate que
dans le cas du slice, la déperdition d'énergie est la plus grande. Dans le même cas, signalons un
effet non visible ici, à savoir que le frottement pendant le contact a pour effet d'inverser le sens de
rotation du disque après le rebondissement.

<figure><img src="image_6.png" /><figcaption>Figure 4 – Disque impactant un plan rigide. Evolution des énergies : (a) $\omega = 0$, (b) $\omega < 0$, (c) $\omega > 0$. (d) Trajectoire du centre du disque.</figcaption></figure>

**5 Conclusion**

Les points marquants de la formulation proposée sont les suivants.

(i) La relation mixte (1) est une généralisation du principe des puissances virtuelles permettant de traiter le problème de contact en grandes transformations avec frottement de Coulomb. Ecrite sous la forme d'une égalité, elle implique les contraintes d'inégalité de contact. D'autre part, la relation (2) exprimant l'équilibre des contraintes de contact est aussi trouvée comme une conséquence de (1).

(ii) La forme faible (1), analogue à un principe de lagrangien augmenté, a l'avantage d'être un véritable principe des résidus pondérés faisant intervenir des champs virtuels arbitraires (suffisamment réguliers). Par conséquent, elle se prête facilement à la discrétisation par éléments finis.

(iii) L'expression discrétisée en temps (5) proposée pour la vitesse de glissement $V_T \Delta t$ ne s'exprime pas en termes de l'antécédent $\eta$ du point de projection $\mathbf{Y}$ et permet ainsi d'éviter la discontinuité indésirable en $\eta$ liée au découpage de $S_{oc}^{(2)}$ en éléments finis. Elle fait plutôt intervenir les variables $\mathbf{y}$ et $\mathbf{Y}_{n-1}$ qui possèdent des significations géométriques intrinsèques indépendantes du maillage de $S_{oc}^{(2)}$.

(iv) Concernant la régularité des fonctions cinématiques, on remarquera que la formulation proposée contient des dérivées d'ordre 2 au plus en espace (voir la définition de la courbure $\kappa_{\alpha\beta}$ après (10)), contrairement à d'autres formulations qui nécessitent des dérivées d'ordre supérieur.

(v) Comme la forme faible (1) inclut les lois de contact, on est dispensé de l'intégration locale de la loi de frottement lors de la résolution numérique. Ce point, analogue à une propriété implicite présentée dans [2], diffère de la plupart des travaux dans la littérature tel que [4].

**Références**

[1] A. Le van, T.H.T. Nguyen. Une formulation variationnelle du problème de contact avec frottement de Coulomb. *Comptes Rendus Mécanique*, Académie des Sciences, Paris, 33, 606-611, 2008.

[2] G. Pietrzak, A. Curnier. Large deformation frictional contact mechanics : continuum formulation and augmented Lagrangian treatment. *Computer Methods in Applied Mechanics and Engineering*, 177, 351–381, 1999.

[3] D. Vola, E. Pratt, M. Jean, M. Raous. Consistent time discretization for a dynamical frictional contact problem and complementary techniques. *Revue Européenne des Eléments Finis*, 7, 1-2-3, 149–162, 1998.

[4] T.A. Laursen. *Computational Contact and Impact Mechanics - Fundamentals of Modeling Interfacial Phenomena in Nonlinear Finite Element Analysis*, Springer, 2003.