# Utilisation des schémas d'intégration classiques pour traiter la plasticité en X-FEM

Alexandre Martin, Patrick Massin

► To cite this version:

Alexandre Martin, Patrick Massin. Utilisation des schémas d'intégration classiques pour traiter la plasticité en X-FEM. 12e Colloque national en calcul des structures, CSMA, May 2015, Giens, France. hal-01516632

HAL Id: hal-01516632

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

Submitted on 2 May 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.

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

Distributed under a Creative Commons CC0 1.0 - Universal - International License
# Utilisation des schémas d'intégration classiques pour traiter la plasticité en X-FEM

A. Martin<fnref n="1" />, P. Massin<fnref n="2" />

**Résumé** — La subdivision de l'élément parent en sous-cellules d'intégration utilisée dans le cadre de la méthode X-FEM convient difficilement à l'étude de la propagation d'une fissure dans un milieu plastifié. Il est alors notamment nécessaire de réaliser la projection des variables internes aux points d'intégration des sous-éléments. Cette contribution étend au cas de la plasticité une méthode d'intégration utilisant les formules de quadrature de Gauss standards, initialement développée dans le cadre de l'élasticité linéaire. La méthode proposée se limite aux cas des interfaces et des éléments linéaires.

**Mots clés** — Interfaces, Méthode des éléments finis étendue, Schémas d'intégration.

## 1 Introduction

Ventura a proposé une méthode d'intégration utilisant les formules de quadrature de Gauss standards pour le calcul de la matrice de rigidité et des forces internes pour un élément X-FEM traversé par une interface et dans le cadre de l'élasticité linéaire [2]. Cette contribution propose d'étendre cette méthode au cas de la plasticité.

## 2 Discrétisation des équations d'équilibre

Considérons un solide $\Omega$ de frontière $\partial\Omega$. La frontière $\partial\Omega$ est composée de $\Gamma_u = \Gamma_u^1 \cup \Gamma_u^2$, $\Gamma_t = \Gamma_t^1 \cup \Gamma_t^2$ et de l'interface $\Gamma$ de sorte que $\partial\Omega = \Gamma_u \cup \Gamma_t \cup \Gamma$. Une force volumique $f$ est appliquée à $\Omega$. Une force surfacique $T$ s'applique sur $\Gamma_t$ ainsi que des efforts surfaciques $t$ sur les lèvres $\Gamma_1$ et $\Gamma_2$ de l'interface $\Gamma$. Le déplacement $w$ est par ailleurs prescrit sur $\Gamma_u$. Les notations précédentes sont illustrées FIGURE 1. Les équations d'équilibre et les conditions aux limites du problème s'écrivent :

$$ 
\begin{align*}
& \nabla \cdot \sigma + f = 0, \text{ dans } \Omega, \\
& \sigma \cdot n = T, \text{ sur } \Gamma_t = \Gamma_t^1 \cup \Gamma_t^2, \\
& \sigma \cdot n_1 = t, \text{ sur } \Gamma_1, \\
& \sigma \cdot n_2 = -t, \text{ sur } \Gamma_2, \\
& u = w, \text{ sur } \Gamma_u = \Gamma_u^1 \cup \Gamma_u^2.
\end{align*}
$$

<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Notations du problème : a) configuration initiale et b) configuration déformée.</figcaption></figure>

<footnote n="1">¹ LaMSID, UMR EDF-CNRS-CEA-ENSTA 8193, 1 Avenue du Général de Gaulle, 92141 Clamart Cedex, alexandre-externe.martin@edf.fr</footnote>

<footnote n="2">² LaMSID, UMR EDF-CNRS-CEA-ENSTA 8193, 1 Avenue du Général de Gaulle, 92141 Clamart Cedex, patrick.massin@edf.fr</footnote>
où $n$, $n_1$ et $n_2$ sont les normales sortantes aux surfaces $\Gamma_l$, $\Gamma_1$ et $\Gamma_2$ respectivement, $\sigma$ le tenseur des contraintes de Cauchy et $u$ le champ de déplacements. Nous resterons dans la suite dans le cadre des petites déformations et petits déplacements.

Le champ de déplacement nodal $u_h$ a des composantes continues $u_r$ et des composantes discontinues $u_d$ de telle sorte qu'il peut être interpolé de la manière suivante :

$$u = \begin{bmatrix} N & HN \end{bmatrix} \begin{Bmatrix} u_r \\ u_d \end{Bmatrix}.$$

où $N$ désignent les fonctions de forme de l'élément fini et $H$ la fonction Heaviside généralisée [1] définie par :

$$H(x) = \begin{cases} +1 & \text{si } lsn(x) > 0 \\ -1 & \text{si } lsn(x) < 0 \end{cases}.$$

où $lsn(x)$ est la valeur de la level set normale au point $M(x)$, la level set normale représentant la distance signée du point $M(x)$ à l'interface $\Gamma$, l'orientation étant arbitraire si l'équation de la level set normale n'est pas donnée. Le tenseur des déformations est alors approximé de la manière suivante :

$$\varepsilon(u) = \begin{bmatrix} B & HB \end{bmatrix} \begin{Bmatrix} u_r \\ u_d \end{Bmatrix}.$$

où $B = \nabla_s N$ avec $\nabla_s$ l'opérateur gradient symétrique. Le champ de déplacement $u_h$ est solution du système d'équations algébriques :

$$f^{\text{int}}(u_h) - f^{\text{ext}} = 0, \quad (1)$$

où $f^{\text{int}}$ est le vecteur des forces internes et $f^{\text{ext}}$ est le vecteur des forces externes.

## 3 Élasticité linéaire

Dans le cadre de l'élasticité linéaire le système (1) est linéaire et peut être réécrit sous la forme :

$$K u_h = f^{\text{ext}}.$$

où $K$ est la matrice de rigidité. Soit $e$ un élément X-FEM occupant le domaine $\Omega_e$. La matrice de rigidité élémentaire $K_e$ s'écrit :

$$K_e = \int_{\Omega_e} \begin{bmatrix} B^T \\ HB^T \end{bmatrix} E \begin{bmatrix} B & HB \end{bmatrix} d\Omega.$$

où $E$ est le tenseur de Hooke. Classiquement le calcul de $K_e$ repose sur une partition de $\Omega_e$ conforme à l'interface $\Gamma$ :

$$K_e = \int_{\Omega_1^e} \begin{bmatrix} B^T EB & (+1)B^T EB \\ (+1)B^T EB & B^T EB \end{bmatrix} d\Omega + \int_{\Omega_2^e} \begin{bmatrix} B^T EB & (-1)B^T EB \\ (-1)B^T EB & B^T EB \end{bmatrix} d\Omega,$$

où $\Omega_e = \Omega_1^e \cup \Omega_2^e$. L'élément parent $e$ est divisé en sous-éléments conformes à l'interface $\Gamma$ pour lesquels les formules de quadrature de Gauss peuvent être utilisées.

Ventura [2] a introduit une fonction $\tilde{H}$, continue sur $\Omega_e$, non uniforme et telle que :

$$\int_{\Omega_e} \tilde{H} B^T E B d\Omega = \int_{\Omega_1^e} (+1) B^T E B d\Omega + \int_{\Omega_2^e} (-1) B^T E B d\Omega.$$

La matrice de rigidité élémentaire devient alors :

$$K_e = \int_{\Omega_e} \begin{bmatrix} B^T E B & \tilde{H} B^T E B \\ \tilde{H} B^T E B & B^T E B \end{bmatrix} d\Omega.$$

Le calcul de $K_e$ s'appuie alors sur les formules de quadrature de Gauss dans l'élément entier.
Les composantes de $B^T E B$ appartiennent à l'espace $\mathcal{P}$ des fonctions polynomiales de degré inférieur ou égal à $d$, où $d$ dépend du type d'élément. Ventura choisit de chercher $\tilde{H} \in \mathcal{P}$ tel que :

$$ \forall f \in \mathcal{P}, \int_{\Omega_e} \tilde{H} f d\Omega = \int_{\Omega_1^e} (+1) f d\Omega + \int_{\Omega_2^e} (-1) f d\Omega. $$

Soit $\mathcal{P}$ une base de $\mathcal{P}$ tel que $\tilde{H}$ se décompose sous la forme :

$$ \tilde{H} = \mathcal{P}h, $$

et $h$ est solution du système linéaire :

$$ \underbrace{\left( \int_{\Omega_e} P^T P d\Omega \right)}_{=V} h = \underbrace{\int_{\Omega_1^e} (+1) P^T d\Omega + \int_{\Omega_2^e} (-1) P^T d\Omega}_{=b} \quad (2) $$

Ventura a calculé $\mathcal{P}$ pour les éléments triangles, quadrangles et tétrahèdres ainsi que résolu analytiquement le système (2) pour une position de l'interface donnée.

Cette contribution propose une stratégie générique pour résoudre le système (2) où l'intégration de la matrice de masse $V$ utilise les formules de quadratures de Gauss dans l'élément parent, tandis que l'intégration du second membre $b$ s'appuie sur le découpage de l'élément parent en sous-éléments. Cette méthode permet d'étendre l'approche aux éléments hexaèdres et n'exige aucune hypothèse supplémentaire sur la position de l'interface [3].

## 4 Extension à la plasticité

Dans le cas d'un comportement non linéaire, l'algorithme de Newton-Raphson est utilisé pour trouver la solution du système algébrique non linéaire (1). Le système linéaire à résoudre pour l'intégration $n$ et l'incrément de chargement $i$ s'écrit :

$$ K_{T_i}^{n-1} \delta u_{h_i}^n = f_i^{\text{ext}} - f_i^{\text{int},n-1}. $$

La matrice de rigidité tangente élémentaire $K_{T_i}^{n-1}$ s'écrit :

$$ K_{e_i}^{n-1} = \int_{\Omega_1^e \cup \Omega_2^e} \left[ \begin{bmatrix} B^T \\ HB^T \end{bmatrix} \right] E_{T_i}^{n-1} \begin{bmatrix} B & HB \end{bmatrix} d\Omega, $$

avec :

$$ E_{T_i}^{n-1} = E_T(\varepsilon(u_{e_i}^{n-1})). $$

où $E_T$ est le comportement tangent. Le vecteur des forces internes élémentaire $f_i^{\text{int},n-1}$ s'écrit :

$$ f_{e_i}^{\text{int},n-1} = \int_{\Omega_1^e \cup \Omega_2^e} \left[ \begin{bmatrix} B^T \\ HB^T \end{bmatrix} \right] \sigma_i^{n-1} d\Omega, $$

où

$$ \sigma_i^{n-1} = \sigma(\varepsilon(u_{e_i}^{n-1})). $$

Soit $\varepsilon|_{H=+1}$ le tenseur de déformation obtenu en supposant que $H = +1$ dans l'élément entier et soit $\varepsilon|_{H=-1}$ le tenseur de déformation obtenu en supposant que $H = -1$ dans l'élément entier. Le tenseur des contraintes peut être décomposé sous la forme :

$$ \sigma_i^{n-1} = \bar{\sigma} + H \Delta \sigma, \quad (3) $$

avec

$$ \bar{\sigma} = \frac{1}{2} [\sigma(\varepsilon|_{H=+1}(u_{e_i}^{n-1})) + \sigma(\varepsilon|_{H=-1}(u_{e_i}^{n-1}))], \\ \Delta\sigma = \frac{1}{2} [\sigma(\varepsilon|_{H=+1}(u_{e_i}^{n-1})) - \sigma(\varepsilon|_{H=-1}(u_{e_i}^{n-1}))]. $$
<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Sphère creuse en plasticité parfaite : a) géométrie et conditions aux limites b) erreur relative en termes de semi-norme H¹ en fonction du raffinement du maillage.</figcaption></figure>

De même $E_{T_i}^{n-1}$ peut être décomposé sous la forme :

$$E_{T_i}^{n-1} = \bar{E} + H \Delta E, \quad (4)$$

avec

$$\bar{E} = \frac{1}{2} [E_T(\varepsilon|_{H=+1}(\mathbf{u}_{e_i}^{n-1})) + E_T(\varepsilon|_{H=-1}(\mathbf{u}_{e_i}^{n-1}))], \\ \Delta E = \frac{1}{2} [E_T(\varepsilon|_{H=+1}(\mathbf{u}_{e_i}^{n-1})) - E_T(\varepsilon|_{H=-1}(\mathbf{u}_{e_i}^{n-1}))].$$

Il est à noter que cette décomposition nécessite d'évaluer le comportement en chaque point d'intégration deux fois : une fois en supposant $H = -1$ et une fois en supposant $H = +1$.

En utilisant (3), (4) et en remarquant que $H^2 = 1$, on obtient :

$$K_{e_i}^{n-1} = \int_{\Omega_1^e \cup \Omega_2^e} \left[ \begin{array}{ll} B^T \bar{E} B + H B^T \Delta E B & H B^T \bar{E} B + B^T \Delta E B \\ H B^T \bar{E} B + B^T \Delta E B & B^T \bar{E} B + H B^T \Delta E B \end{array} \right] d\Omega, \quad (5)$$

et

$$f_{e_i}^{\text{int},n-1} = \int_{\Omega_1^e \cup \Omega_2^e} \left[ \begin{array}{l} B^T \bar{\sigma} + H B^T \Delta \sigma \\ H B^T \bar{\sigma} + B^T \Delta \sigma \end{array} \right] d\Omega. \quad (6)$$

Finalement, la fonction $\tilde{H}$ peut être substituée à $H$ dans (5) et (6). Les calculs de $K_{e_i}^{n-1}$ et de $f_{e_i}^{\text{int},n-1}$ s'appuient alors sur les formules de quadrature de Gauss dans l'élément entier.

# 5 Tests numériques

La méthode présentée dans cette contribution a été validée par des analyses de convergence en termes de norme en énergie (cas élastiques) et de semi-norme $H^1$ (cas élasto-plastiques) [3]. Les résultats obtenus dans le cas d'une modélisation axisymétrique d'une sphère creuse en plasticité parfaite, soumise à une pression interne, sont présentés FIGURE 2. Les taux de convergence obtenus avec la méthode standard et la méthode proposée sont identiques et en accord avec la valeur théorique attendue.

# Références

[1] N. Moës, J. Dolbow, and T. Belytschko. A finite element method for crack growth without remeshing. *International Journal for Numerical Methods in Engineering*, 46(1):131–150, 1999.

[2] G. Ventura. On the elimination of quadrature subcells for discontinuous functions in the eXtended Finite-Element Method. *International Journal for Numerical Methods in Engineering*, 66(5):761–795, 2006.

[3] A. Martin, J.-B. Esnault, and P. Massin. About the use of standard integration schemes for X-FEM in solid mechanics plasticity. *Computer Methods in Applied Mechanics and Engineering*, 283:551–572, 2015.