# Couplage “ Smoothed Particles Hydrodynamics - Eléments finis ” sous impact par la méthode Arlequin

Yann Chuzel-Marmot, Alain Combescure, Roland Ortiz

► To cite this version:

Yann Chuzel-Marmot, Alain Combescure, Roland Ortiz. Couplage “ Smoothed Particles Hydrodynamics - Eléments finis ” sous impact par la méthode Arlequin. 8e Colloque national en calcul des structures, CSMA, May 2007, Giens, France. hal-01504169

HAL Id: hal-01504169

https://hal.science/hal-01504169

Submitted on 8 Apr 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
# Couplage « SPH – EF »

## Couplage « Smoothed Particles Hydrodynamics – Eléments Finis » sous impact par la méthode Arlequin.

Yann Chuzel-Marmot* — Alain Combescure* — Roland Ortiz**

* LAboratoire de Mécanique des COntacts et des Structures, INSA de Lyon
UMR/CNRS 5514, Bâtiment Jean d'Alembert
18-20 rue des sciences
F-69621 Villeurbanne cedex

Yann.Chuzel-Marmot@insa-lyon.fr, Alain.Combescure@insa-lyon.fr

** Office National d'Etudes et de Recherches Aérospatiales, Centre de Lille
5 boulevard Paul Painlevé
F-59000 Lille cedex

Roland.Ortiz@onera.fr

**RÉSUMÉ.** La méthode Arlequin donne un cadre simple et efficace pour assembler des modèles utilisant des formulations différentes. Elle est développée ici, dans le cadre de la dynamique explicite afin d'assembler une zone maillée avec la méthode Smoothed Particles Hydrodynamics (notée SPH) permettant de simuler aisément des ruptures par fragmentation et une autre zone plus vaste, qui ne s'endommagera pas, maillée avec la méthode des éléments finis (notée EF). Ce document détaille donc les ingrédients de la méthode mise en place au sein du code de calcul EUROPLEXUS, sa validation sur des cas tests simples et une confrontation entre les simulations numériques et les résultats d'une campagne expérimentale sur la résistance de dalles en béton à l'impact de projectiles.

**ABSTRACT.** The Arlequin method gives a simple and effective framework to glue models using various formulations. It is developed here, as part of dynamic explicit in order to link a zone meshed with the Smoothed Particles Hydrodynamics formulation (noted SPH) showing ruptures by fragmentation and an other larger zone, without damage using the finite element method (noted EF). Thus, this paper details the method implemented in the EUROPLEXUS code, its validation on simple benchmarks and a confrontation between numerical simulations and results of an experimental study of concrete slabs resistance to projectile impact.

**MOTS-CLÉS :** Méthode arlequin, Couplage, Meshless, Sph, Dynamique explicite, Impact de projectile, Béton, Mazars.

**KEYWORDS:** Arlequin Method, Coupling, Meshless, Sph, Explicit dynamic, Projectile impact, Concrete, Mazars.
# 1. Introduction

La méthode Arlequin est appliquée ici au collage volumique d'une méthode éléments finis et d'une méthode SPH en dynamique explicite. Le formalisme présenté s'appuie largement sur les travaux de A. Ben Dhia (Ben Dhia, 1998) et de G. Rateau (Rateau, 2003) et se limite pour l'instant au cas de l'élasticité linéaire en formulation lagrangienne totale. Il peut cependant s'étendre à des problèmes plus généraux.

# 2. Présentation du problème continu mixte

Le domaine $\Omega$ est décomposé en deux sous domaines $\Omega_1$ et $\Omega_2$ qui se superposent sur une partie $V_g$ de volume non nul. Appelons $\Omega_1^H$ et $\Omega_2^H$ les parties des sous domaines qui se superposent. Les états mécaniques doivent alors être identiques dans les volumes $\Omega_1^H$ et $\Omega_2^H$. Pour ce faire, nous allons contraindre la restriction des déplacements des deux sous domaines $\Omega_1^H$ et $\Omega_2^H$ à être identiques. L'opérateur de couplage est alors défini par :

- un modèle de jonction : suivant soit un modèle rigide (approche lagrangienne) soit un modèle élastique (méthode de pénalisation),

- un espace médiateur : l'ensemble des éléments de comparaison (déplacements, déformations, ...) forme un espace pivot nommé médiateur,

- un niveau de distribution d'énergie entre les modèles : chaque travail virtuel est pondéré par des fonctions de mélange $(\alpha, \beta)$ formant une partition de l'unité afin de ne pas doubler l'énergie liée à la zone de recouvrement.

Dans la formulation Arlequin, chaque sous domaine i respecte donc indépendamment les équations habituelles de la mécanique des milieux continus mais est contraint sur la zone de recouvrement par les efforts $C_i$:

$$ \frac{d}{dt} \int_{\Omega_i^{H0}} \rho_0 \underline{v}_i \cdot dV^0 = \int_{\Omega_i^{H0}} \alpha_i \nabla_X P \underline{f}_i \cdot dV^0 + \int_{\Omega_i^{H0}} \beta_i \rho_0 \underline{f}_i \cdot dV^0 + C_i \quad \text{pour } i = 1,2 $$

avec
* $\rho_0$, la masse volumique initiale,
* $\underline{v}$, le champ vitesse,
* $P$, la transposée du premier tenseur des contraintes de Piolat-Kirchhoff,
* $f$, les forces volumiques.
## 3. Formalisme du couplage « SPH – EF »

### 3.1. Choix du modèle de l'opérateur de couplage

Pour cette recherche nous avons choisi un couplage par multiplicateurs de Lagrange et un espace médiateur basé sur le seul champ de déplacements. Le terme de couplage se définit alors par le produit scalaire suivant :

$$C_i = C(\lambda, v_i) = \int_{V_g^0} \lambda \cdot v_i \quad \text{avec} \quad v_i = \underline{u}_{SPH} \quad \text{ou} \quad v_i = -\underline{u}_{EF}$$

avec $\lambda$, les multiplicateurs de Lagrange et $\underline{u}$, le champ déplacement.

### 3.2. Choix du niveau de distribution d'énergie

Dans le cas de couplage avec formulation SPH, le choix des fonctions de mélange $(\alpha, \beta)$ est assez restreint : les fonctions noyaux SPH habituelles ne supportent pas des variations brutales de propriétés mécaniques du matériau. Nous avons donc choisi des fonctions de mélange constantes sur la zone de recouvrement et égales à l'unité sur le reste de leur domaine, telles que :

$$\alpha_{SPH} = \beta_{SPH} = 0,999 \quad \text{sur} \quad V_g^0 \\ \alpha_{EF} = \beta_{EF} = 0,001$$

### 3.3. Interpolation des déplacements

L'interpolation des valeurs d'un champ pour la méthode SPH permet de déterminer la valeur approchée $u^\#$ des déplacements aux points de gauss éléments finis. Une méthode « Moving Least Square » (Rabczuk et al., 2004) à l'ordre 1 est utilisée :

$$\forall \underline{X}(x, y, z), u^\#(X, t) = \sum_{i=1}^{N_{\text{voisin}}} M_i(X) u_i(t), \left\{ \begin{array}{l} M_i = \sum_{j=1}^{T} \underline{p}(X) \underline{A}(X)^{-1} \underline{Q}(X_i) \cdot \underline{W}^0 \\ \underline{A}(X) = \sum_{i=1}^{N_{\text{voisin}}} \underline{Q}(X_i) \underline{W}^0 \underline{Q}(X_i) \\ \underline{p}(X) = [1 \ x \ y \ z] \\ \underline{Q}_{ij} = p_i(X_j) \quad \text{et} \quad \underline{W}_{ii}^0 = w^0(r_i, h) \end{array} \right.$$

avec $w^0$, la fonction noyau SPH lagrangienne totale définie sur la configuration initiale.
### 3.4. Ecriture de la matrice de couplage

La construction des matrices de couplage s'obtient par le biais de l'expression du travail de l'interface entre sous domaine :

$$E_{\text{interface}} = C_{\text{SPH}} + C_{\text{EF}} = \int_{V_g^0} \int_T \lambda (\underline{u}^{\text{SPH}} - \underline{u}^{\text{EF}}) dV_0$$

En discrétisant les multiplicateurs sur le même espace que les éléments finis et en se transposant sur le système de coordonnées parents, on obtient ainsi la relation discrétisée. En considérant alors un maillage compatible avec des éléments finis identiques, la nullité du travail d'interface permet d'identifier les matrices de couplage.

### 3.5. Schéma d'intégration explicite

Le schéma d'intégration temporelle est le schéma de dynamique explicite basé sur un schéma de type Newmark aux différences centrées. La condition de nullité du travail sur la zone de recouvrement permet d'exprimer les forces fictives de couplage entre les deux sous domaines :

$$\underline{F}_c^{n+1} = \underline{P}_c^T \cdot \lambda^{n+1}$$

$$\underline{\lambda}^{n+1} = [\underline{B}^{n+1}]^{-1} \cdot \underline{W}^{n+1} \quad \text{avec} \quad \underline{B}^{n+1} = \underline{P}_c \cdot \underline{m}^{-1} \cdot \underline{P}_c^T \cdot \underline{W}^{n+1} = \underline{P}_c \cdot \underline{a}^{n+1} - \underline{P}_c \cdot \underline{m}^{-1} \cdot (\underline{F}_{\text{ext}}^{n+1} - \underline{F}_{\text{int}}^{n+1})$$

avec $\underline{F}_c$, les forces de couplage.

## 4. Exemple d'application numérique

### 4.1. Cas simple d'une poutre

Plusieurs cas de chargement ont été testés pour une poutre de longueur 36 mm et de section 5×5 mm :

- une vitesse initiale imposée selon la forme $v = e^{-0,025(x-x_m)^2}$ comme dans la référence (Rabczuk et al., 2004) : les temps et l'amplitude des 2 premiers extremums (notés EXT- par la suite) de la réponse dynamique en vitesse des méthodes numériques EF et mixte (2/3 EF - 1/3 SPH) sont comparés,

- un essai de flexion simple à effort imposé (3000 N) : les temps et l'amplitude des 2 premiers extremums de la réponse dynamique en déplacement d'un maillage mixte 1:4 sont comparés à ceux d'un maillage purement EF.
<table><caption>**Tableau 1.** Comparaison analytique – simulations
lors de l'essai à vitesse imposée.</caption><thead><tr><th></th><th>Temps<br>EXT-1</th><th>EXT-2</th><th>Amplitude<br>EXT-1</th><th>EXT-2</th></tr></thead><tbody><tr><td>Valeur analytique</td><td>9,97 µs</td><td>19,94 µs</td><td>1 m.s<sup>-1</sup></td><td>1 m.s<sup>-1</sup></td></tr><tr><td>EF (Section 5x5)</td><td>10,03 µs</td><td>20,06 µs</td><td>0,999 m.s<sup>-1</sup></td><td>0,997 m.s<sup>-1</sup></td></tr><tr><td>Mixte (EF 5x5 – SPH 5x5)</td><td>10,39 µs</td><td>21,15 µs</td><td>0,894 m.s<sup>-1</sup></td><td>0,914 m.s<sup>-1</sup></td></tr><tr><td>Mixte (EF 5x5 – SPH 10x10)</td><td>10,25 µs</td><td>20,53 µs</td><td>0,961 m.s<sup>-1</sup></td><td>0,964 m.s<sup>-1</sup></td></tr><tr><td>Mixte (EF 5x5 – SPH 15x15)</td><td>10,21 µs</td><td>20,37 µs</td><td>0,976 m.s<sup>-1</sup></td><td>0,974 m.s<sup>-1</sup></td></tr><tr><td>Mixte (EF 5x5 – SPH 20x20)</td><td>10,20 µs</td><td>20,33 µs</td><td>0,981 m.s<sup>-1</sup></td><td>0,977 m.s<sup>-1</sup></td></tr><tr><td>Mixte (EF 5x5 – SPH 25x25)</td><td>10,19 µs</td><td>20,31 µs</td><td>0,985 m.s<sup>-1</sup></td><td>0,979 m.s<sup>-1</sup></td></tr></tbody></table>

<table><caption>**Tableau 2.** Comparaison modèle EF - modèle “SPH – EF couplé”
lors de l'essai de flexion simple.</caption><thead><tr><th></th><th>Temps<br>EXT-1</th><th>EXT-2</th><th>Amplitude<br>EXT-1</th><th>EXT-2</th></tr></thead><tbody><tr><td>EF (Section 5x5)</td><td>166,6 µs</td><td>489,4 µs</td><td>22,91 mm</td><td>23,14 mm</td></tr><tr><td>Mixte (EF 5x5 – SPH 20x20)</td><td>168,1 µs</td><td>496,9 µs</td><td>22,82 mm</td><td>22,83 mm</td></tr></tbody></table>

## 4.2. Seconde application : Impact d'un projectile sur une dalle béton Description

M. H. Zhang et al. (Zhang et al., 2000) ont réalisé de multiples expériences afin de caractériser la résistance du béton haute performance sous impact. Un test a été simulé : l'impact d'un projectile de 15 g à tête ogivale de 12,6 mm de diamètre sur une dalle en béton de 300×170×150 mm à une vitesse de 664 m.s⁻¹ (contraintes limites : traction – 19,3 MPa ; compression – 115 MPa ; module de Young : 30 GPa ; coefficient de poisson : 0,22 ; masse volumique : 2300 kg.m⁻³).

La loi de comportement utilisée est basée sur les travaux de J. Mazars (Mazars, 1984) auxquels s'ajoutent des modifications inhérentes à l'utilisation du modèle en dynamique rapide. En effet, la déformation équivalente proposée par J. Mazars n'est pas adaptée à la rupture en dynamique. Sous un chargement de compression, la déformation équivalente est déterminée à partir des déformations positives. Or, sous chargement rapide, des ruptures en compression sont observées qui ne peuvent pas être prédites avec ce choix de déformation équivalente. La mesure de la déformation équivalente a donc été modifiée afin d'endommager le matériau sous un chargement de compression rapide. La formule suivante a été retenue :

$$ \tilde{\epsilon} = \sqrt{\sum_{i} (H_{i} \epsilon_{i})^{2}} + \sqrt{2\nu} \sqrt{\sum_{i} ((1 - H_{i}) \epsilon_{i})^{2}} \text{ avec } H_{i} \text{ la fonction de Heaviside.} $$
L'utilisation d'une limitation du taux de croissance évite les localisations artificielles obtenues avec les lois adoucissantes et rend le calcul indépendant du maillage. (Suffis *et al.*, 2002)

$$ \begin{cases} \dot{D} = 1/\tau_c \exp[-a(D_{nc} - D)] \text{ si } D \le D^c \\ D = D^c \text{ sinon} \end{cases} $$

avec a et $\tau_c$, 2 paramètres matériaux - nos valeurs : a=1 et $\tau_c=8\mu$s, D, $D_{nc}$ et $D^c$ les endommagements courant, non corrigé et critique.

<table><caption>**Tableau 3. Comparaison expérience - simulation.**</caption><thead><tr><td></td><td>Profondeur</td><td>Diamètre du cratère</td></tr></thead><tbody><tr><td>Expérience (moyenne sur 2 tests)</td><td>31,75 mm</td><td>63,5 mm</td></tr><tr><td>Simulation</td><td>33 mm</td><td>61,6 mm</td></tr></tbody></table>

## 6. Conclusion

Le formalisme Arlequin a été appliqué au cas du couplage « SPH – EF » et implémenté dans le code de calcul dynamique explicite EUROPLEXUS. Les résultats sur des cas tests simples et la comparaison avec une expérience d'impact attestent de la validité de l'approche. D'autres cas tests et comparaisons expérimentales sont en cours de réalisation.

## 7. Bibliographie

H. BenDhia, « Problèmes mécaniques multi-échelles : la méthode arlequin », *Comptes Rendus de l'Académie des Sciences*, Série IIb, vol. 36, 1998, p. 899-904.

J. Mazars, « Application de la mécanique de l'endommagement au comportement non linéaire et à la rupture du béton de structure », Thèse de doctorat d'état, Université Paris VI, 1984.

T. Rabczuk, T. Belytschko, S.P. Xiao, « Stable particle methods on Lagrangian kernels », *Comput. Methods Appl. Engrg.*, vol. 193, 2004, p. 1035-1063.

G. Rateau, « Méthode Arlequin pour les problèmes mécaniques multi-échelles. Application à des problèmes de jonction et de fissuration de structures élancées », Thèse Ecole Centrale Paris, 2003.

A. Suffis, A. Combescure, « Modèle d'endommagement à effet retard, analyse numérique et analytique de l'évolution de la longueur caractéristique », *Revue Européenne des Eléments Finis*, vol. 11, 2002, p. 593-620.

M.H. Zhang, V.P.W. Shim, G. Lu, C.W.Chew, « Resistance of high-strength concrete to projectile impact », *Int. Journal of Impact Engineering*, vol. 31, 2005, n°7, p. 825-841.