# Méthodes sans maillage de type éléments naturels pour la simulation des procédés de mise en forme

Iciar Alfaro, Julien Yvonnet, Elías Cueto, Francisco Chinesta, Pierre Villon, Manuel Doblaré

► To cite this version:

Iciar Alfaro, Julien Yvonnet, Elías Cueto, Francisco Chinesta, Pierre Villon, et al. Méthodes sans maillage de type éléments naturels pour la simulation des procédés de mise en forme. Revue Européenne de Mécanique Numérique/European Journal of Computational Mechanics, 2006, 15 (1-2-3), pp.29-40. 10.3166/remn.15.29-40. hal-01007032

HAL Id: hal-01007032

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

Submitted on 9 Oct 2016

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 CC BY 4.0 - Attribution - International License
# Méthodes sans maillage de type éléments naturels pour la simulation des procédés de mise en forme

Iciar Alfaro\* — Julien Yvonnet\*\* — Elias Cueto\*  
Francesco Chinesta\*\* — Pierre Villon\*\*\* — Manuel Doblaré\*

\* Aragon Institute of Engineering Research. University of Zaragoza
Edificio Betancourt. María de Luna, 5, E-50018. Zaragoza, Spain
iciar@unizar.es ecueto@unizar.es mdoblare@unizar.es

\*\* Laboratoire de Mécanique des Systèmes et des Procédés
UMR 8106 CNRS-ENSAM-ESEM, 151 boulevard de l'Hôpital, F-75013 Paris
julien.yvonnet@paris.ensam.fr francisco.chinesta@paris.ensam.fr

\*\*\* Laboratoire Roberval, Université de Technologie de Compiègne
BP 20529, F-60205 Compiègne
pierre.villon@utc.fr

RÉSUMÉ. Dans cet article, des exemples de simulations de différents procédés sont présentés. Une méthode de type éléments naturels est utilisée. L'une des difficultés de la méthode est le traitement des domaines non convexes. Pour contourner le problème, deux approches sont utilisées, la méthode des formes alpha, et la méthode C-NEM. Une autre difficulté est l'amélioration de la qualité en convergence, la NEM ne possédant qu'une consistance de premier ordre. Une nouvelle interpolation de type Hermite NEM est proposée pour construire une approximation NEM de second ordre en convergence.

ABSTRACT. In this paper, numerical simulations of forming processes are presented. A natural neighbor Galerkin method is used. One issue in this technique is the treatment of non-convex domains. To circumvent this problem, two approaches are used, the alpha-NEM and the C-NEM. Another difficulty in the NEM is to achieve p-adaptivity. For this purpose, a new Hermite-NEM scheme is proposed to construct higher-order approximation schemes in the NEM.

MOTS-CLÉS : méthodes sans maillage, éléments naturels, simulation des procédés.

KEYWORDS: meshless methods, natural elements, forming processes simulation
# 1. Introduction

La maîtrise des procédés de mise en forme passe par l'élaboration d'un modèle mécanique, prenant en compte tous les aspects multiphysiques, et d'un modèle numérique associé, bien adapté au caractère 3D, multi-échelle (localisation en espace-temps des champs de déformation, des contraintes, des températures, de l'endommagement, etc. Il est nécessaire de prendre également en compte la localisation précise des éventuelles fissures, mais aussi les grandes transformations qui caractérisent les différents procédés. Les limitations de l'approche numérique pour la simulation des procédés de mise en forme des matériaux sont encore actuellement nombreuses. L'une des limitations caractéristiques des approches de type éléments finis dans la simulation des procédés est la nécessité de remailler régulièrement la pièce. Ceci est la conséquence des déformations importantes subies par la matière. Le remaillage est une opération complexe et coûteuse, et qui nécessite une projection des champs de variables internes entre l'ancien et le nouveau maillage. Pour alléger ces difficultés, une option est l'utilisation des méthodes dites sans maillage (meshless methods), alternatives à la méthode des éléments finis. De par leur développement récent, de nombreuses difficultés numériques restent encore à aborder. Une future analyse complète de la simulation des procédés dans le cadre sans maillage implique l'étude des points suivants :

- les très grandes transformations géométriques et la problématique du remaillage, ainsi que la projection des champs (Martínez et al. 2004);

- la description précise des bords du domaine matériel : la problématique liée aux surfaces libres ; la gestion des surfaces mobiles (Yvonnet et al. 2005a);

- la localisation et le raffinage adaptatif, l'estimation de l'erreur (Yvonnet et al. 2005b);

- la gestion des fissures et la problématique de la dépendance de sa direction de propagation vis-à-vis de la méthode de discrétisation (recherche en cours);

- les formulations mixtes dans le cadre des modélisations multi-physiques couplées : prise en compte de l'incompressibilité et des conditions inf-sup (González et al. 2004, Yvonnet et al. n.d.);

- les formulations multi-domaines capables de profiter du calcul parallèle (recherche en cours).

# 2. La méthode des éléments naturels (NEM) : fondements et avancées récentes

La notion de voisins naturels et d'interpolation naturelle fut introduite par Sibson en 1980 pour construire des interpolations sur des nuages de nœuds très irréguliers (déstructurés). Cet interpolant est basé sur les constructions géométriques connues sous le nom de diagramme de Voronoi ou de triangulation de Delaunay. Le diagramme de Voronoi est une construction géométrique réalisée à partir d'un ensemble de nœuds distincts $n_1,..., n_N$ dont les cellules sont définies par :

$$ T_I = \{ \mathbf{x} \in \mathbb{R}^n : d(\mathbf{x}, n_I) < d(\mathbf{x}, n_J) \quad \forall J \neq I \} \quad [1] $$
<figure><img src="image_2.png" /><figcaption>**Figure 1.** (a) Définition des voisins naturels et ses coordonnées ; (b) Fonction de forme NEM (permission N. Sukumar)</figcaption></figure>

où d est une distance et $T_I$ est la cellule de Voronoi associée au nœud $n_I$. La cellule de Voronoi $T_I$ est une région de l'espace telle que tout point x à l'intérieur de cette cellule est plus proche du nœud $n_I$ que de n'importe quel autre nœud $n_J$. Deux nœuds partageant une face de leurs cellules de Voronoi respectives sont appelés voisins naturels. On peut définir une cellule de second ordre $T_{IJ}$ de la façon suivante

$$ T_{IJ} = \{ \mathbf{x} \in \mathbb{R}^n : d(\mathbf{x}, n_I) < d(\mathbf{x}, n_J) < d(\mathbf{x}, n_K) \quad \forall K \neq J; \forall K \neq I \}. \quad [2] $$

Un point quelconque se trouvant dans $T_{IJ}$ a comme voisin le plus proche $n_I$ et comme second voisin le plus proche $n_J$.

Si l'on considère le cas bidimensionnel de la figure 1(a), la fonction de forme associée au nœud $n_1$ calculée en $\mathbf{x}$ est donnée par le rapport entre les aires de la cellule de second ordre définie par $n_1$ et $\mathbf{x}$, et celle de premier ordre définie par $\mathbf{x}$ :

$$ \phi_1^{sib}(\mathbf{x}) = \frac{\text{Area}(T_{x1})}{\text{Area}(T_x)} = \frac{\text{Area}(abfe)}{\text{Area}(abcd)}. \quad [3] $$

La figure 1(b) montre l'allure d'une fonction de forme. Le schéma d'interpolation pour un champ $\mathbf{u}$, scalaire, vectoriel ou tensoriel, est de la forme :

$$ \mathbf{u}(\mathbf{x}) = \sum_{i=1}^{N} \phi_I(\mathbf{x}) \mathbf{u}_I \quad [4] $$

$\mathbf{u}_I$ représentant les valeurs nodales du champ considéré.

Les propriétés fondamentales des fonctions de forme de la méthode des éléments naturels sont : (i) l'interpolation ; (ii) la partition de l'unité ; (iii) la consistance linéaire (iv) la linéarité stricte de l'approximation sur les bords des domaines convexes et (v) l'adaptation automatique de la connectivite pour n'importe quelle discrétisation.
La problématique liée aux domaines non convexes est liée aux influences parasites entre ces bords. Pour éviter ces influences il existe deux alternatives principales :

– l'introduction des formes alpha. Cette approche est bien adaptée aux simulations avec bords libres comme les simulations de fluides, car elle permet d'extraire le bord du domaine à partir du nuage de nœuds sans avoir besoin d'une description explicite du bord (Cueto et al. 2003);

– l'introduction d'un critère de visibilité afin de restreindre les voisins naturels d'un point, lorsque les voisins ne sont pas visibles depuis le point, en considérant les bords comme opaques (Yvonnet et al. 2004). Cette approche est plus adaptée à l'étude des solides, ou de fortes concavités telles que les fissures ne pouvant être décrites efficacement par les formes alpha. Par contre, le suivi des bords libres avec cette technique est délicat.

La simulation des procédés peut être caractérisée par l'apparition de localisation de déformations ou d'autres champs. Il faut alors adapter la densité nodale pour bien représenter ces localisations. Nous avons dans ce but proposé différents indicateurs d'erreur, basés sur la différence entre une contrainte évaluée dans les nœuds et interpolée partout dans le domaine, et une contrainte équivalente définie dans chaque cellule de Voronoi. On peut alors, par une stratégie basée sur le diagramme de Voronoi, raffiner le nuage même si celui-ci est non structuré. Un exemple de raffinement dans le cadre de la simulation du cisaillage adiabatique est présenté sur la figure 7 (Yvonnet et al. 2005b).

Le concept de diagramme de Voronoi contraint (Yvonnet et al. 2004) permet de modéliser des milieux hétérogènes à l'échelle macroscopique ainsi que des interfaces mobiles. En construisant un diagramme de Voronoi contraint dans chaque sous-domaine (on peut imposer différentes conditions de transmission, avec continuité où discontinuité de la fonction interpolée) (Yvonnet et al. 2005a), la fonction interpolée devient linéaire sur les bords de chaque sous-domaine, assurant une continuité stricte entre les domaines. Cette propriété a été exploitée pour simuler la propagation des interfaces dans les problèmes de changement de phase (Yvonnet et al. 2005a). De même, si l'on considère deux degrés de liberté dans chaque nœud placé sur une interface, nous pouvons représenter également une transmission avec discontinuité. Les domaines hétérogènes ont été également traités dans le cadre des α-NEM (Cueto et al. 2003).

Les simulations concernant les écoulements incompressibles ont été traitées par une formulation mixte vitesse-pression (le premier champ étant interpolé dans le cadre NEM, le deuxième étant discontinu – constant à l'intérieur de chaque cellule de Voronoi –). Bien que ce schéma ne vérifie pas la condition inf-sup (LBB), nous n'avons jusqu'à présent rencontré aucun blocage numérique. Pour s'affranchir de la vérification de la condition LBB González et al. (2004) proposent un enrichissement de l'interpolation des vitesses dans le cadre de la partition de l'unité (PUM). L'approximation retenue est similaire à celle du mini-élément (élément fini bulle). L'introduction des fonctions bulle (comme fonctions d'enrichissement de l'interpolation des vitesses)
construites à partir du produit des fonctions de forme C-NEM des trois sommets de chaque triangle de Delaunay permet également de vérifier la condition LBB (Yvonnet et al. n.d.).

Enfin, la supériorité de la méthode NEM (par rapport aux autres techniques sans maillage) pour décrire avec précision les bords et interfaces, ouvre une perspective nouvelle pour une modélisation multidomaine (FETI par exemple) ainsi que pour un possible couplage avec les éléments finis.

### 3. Formulation de type Hermite

Une voie d'amélioration de la NEM est l'amélioration de la vitesse de convergence, seulement de premier ordre. Une possibilité pour obtenir une consistance d'ordre supérieur est d'introduire les fonctions de forme NEM comme des fonctions poids dans un schéma de type moindres carrés mobiles, en utilisant une base polynomiale d'ordre supérieur. Cependant, il a été montré dans (Liu et al. 1995) qu'un point d'intégration doit être influencé par un nombre suffisant de nœuds pour reproduire la base, c'est-à-dire que le nombre de nœuds influents doit être supérieur ou égal au nombre de monômes contenus dans la base. Dans la méthode des éléments naturels, le nombre de nœuds influents ne peut être choisi par l'utilisateur. Dans le pire des cas, un point d'intégration est contenu dans un triangle de Delaunay unique. Il n'y a dans ce cas que trois nœuds influents, ce qui ne permet de reproduire qu'une base linéaire. L'idée de l'approche présentée ici est d'introduire des degrés de liberté supplémentaires aux nœuds afin de pouvoir reproduire des polynômes d'ordre supérieur. Dans ce but nous considérons un schéma d'approximation de la forme :

$$u^h(\mathbf{x}) = \sum_{i=1}^{n} \psi_i(\mathbf{x})u_i + \sum_{i=1}^{n} \psi_i^x(\mathbf{x})\frac{\partial u_i}{\partial x} + \sum_{i=1}^{n} \psi_i^y(\mathbf{x})\frac{\partial u_i}{\partial y} \quad [5]$$

où $\psi_i(\mathbf{x})$ sont les fonctions de forme associées avec les inconnues principales $u_i$, $\psi_i^x(\mathbf{x})$ et $\psi_i^y(\mathbf{x})$ sont les fonctions de forme associées avec les dérivées spatiales de $u_i$ par rapport à $x$ et $y$, respectivement. Dans ce contexte, $u_i$, $\frac{\partial u_i}{\partial x}$ et $\frac{\partial u_i}{\partial y}$ sont des inconnues du problème (degrés de liberté). Afin de construire les fonctions de forme, nous considérons le schéma d'approximation suivant :

$$u^h(\mathbf{x}) = \mathbf{p}(\mathbf{x})^T \mathbf{a} \quad [6]$$

où $\mathbf{p}(\mathbf{x})$ est une base polynomiale, par exemple $\mathbf{p}(\mathbf{x}) = \{1, x, y, xy, x^2, y^2\}$ et $\mathbf{a}$ est un vecteur de coefficients inconnus. Afin de déterminer $\mathbf{a}$, nous considérons la fonctionnelle suivante :

$$J = \frac{1}{2} \sum_{i=1}^{n} w_i(\mathbf{x}) \left\{ \left[ \mathbf{p}^T(\mathbf{x})\mathbf{a} - u_i \right]^2 + \alpha \left[ \frac{\partial \mathbf{p}^T(\mathbf{x})\mathbf{a}}{\partial x} - \frac{\partial u_i}{\partial x} \right]^2 \right. \\ \left. + \alpha \left[ \frac{\partial \mathbf{p}^T(\mathbf{x})\mathbf{a}}{\partial y} - \frac{\partial u_i}{\partial y} \right]^2 \right\} \quad [7]$$
où $n$ est le nombre de voisins naturels du point $x$, $w_i(\mathbf{x})$ sont les fonctions de forme NEM calculées au point $x$, $\frac{\partial \mathbf{p}^T}{\partial x}(\mathbf{x})$ et $\frac{\partial \mathbf{p}^T}{\partial y}(\mathbf{x})$ représentent les dérivées de la base polynomiale par rapport à $x$ et $y$, respectivement, $\alpha$ étant un paramètre dimensionnel que nous fixons à 1 dans nos simulations. En minimisant $J$ par rapport à $\mathbf{a}$, ($\frac{\partial J}{\partial \mathbf{a}} = 0$), on obtient le système d'équations suivant :

$$A\mathbf{a}(\mathbf{x}) = B\mathbf{q} \qquad [8]$$

avec $\mathbf{q} = \{u_1, \frac{\partial u_1}{\partial x}, \frac{\partial u_1}{\partial y}, u_2, \frac{\partial u_2}{\partial x}, \frac{\partial u_2}{\partial y}, \dots, u_N, \frac{\partial u_N}{\partial x}, \frac{\partial u_N}{\partial y}\}$. Les matrices $A$ et $B$ sont exprimées par :

$$A_{ij}(\mathbf{x}) = \sum_{k=1}^{n} w_k(\mathbf{x}) \left\{ p_i(\mathbf{x}_k) p_j(\mathbf{x}_k) + \alpha \frac{\partial p_i(\mathbf{x}_k)}{\partial x} \frac{\partial p_j(\mathbf{x}_k)}{\partial x} + \alpha \frac{\partial p_i(\mathbf{x}_k)}{\partial y} \frac{\partial p_j(\mathbf{x}_k)}{\partial y} \right\} \qquad [9]$$

$$B_{i(3j-2)} = w_j(\mathbf{x})p_i(\mathbf{x}_j) \\ B_{i(3j-1)} = \alpha w_j(\mathbf{x}) \frac{\partial p_j(\mathbf{x})}{\partial x} \\ B_{i(3j)} = \alpha w_j(\mathbf{x}) \frac{\partial p_j(\mathbf{x})}{\partial y}$$

Les dérivées des fonctions de forme sont obtenues de manière standard (Belytschko et al. 1994) et sont fonction des dérivées des poids $w_i(\mathbf{x})$. Une expression analytique des dérivées des fonctions de forme de Sibson peut être trouvée dans (Piper 1995). Si l'on réecrit l'expression des fonctions de forme comme :

$$\phi_i(\mathbf{x}) = \frac{\omega_i}{\sum_{j=1}^{j=n} \omega_j} \qquad [10]$$

où $\omega_i$ est la mesure de la cellule du second ordre $T_{ix}$ dont l'expression du gradient $\nabla\omega_i$ est donnée par

$$\nabla\omega_i = \frac{f_i}{d_i}(\mathbf{c}_i - \mathbf{x}) \qquad [11]$$

avec $f_i$ la longueur (aire en 3D) de la face de Voronoi entre $\mathbf{x}$ et $n_i$, $d_i$ la distance entre $\mathbf{x}$ et $n_i$, et $\mathbf{c}_i$ les coordonnées du centre de gravité de la face de Voronoi entre $\mathbf{x}$ et $n_i$. Les dérivées des fonctions de forme peuvent être obtenues à partir des équations [10] et [11]. Les fonctions de forme obtenues sont représentées sur la figure 2.

### 3.1. Exemple numérique dans le cadre de l'approximation H-NEM

Afin d'évaluer l'approximation H-NEM, un problème simple de Poisson est présenté. Le problème aux conditions aux limites est donné par :

$$\left\{ \begin{array}{ll} -\Delta u = f & \text{dans } \Omega = ]0, 1[ \times ]0, 1[ \\ u = \bar{u} & \text{dans } \Gamma_u, \end{array} \right. \qquad [12]$$

Nous choisissons :

$$\left\{ \begin{array}{l} \bar{u} = 0, \\ f = 4\pi^2 \{2\cos(2\pi x)\cos(2\pi y) - \cos(2\pi x) - \cos(2\pi y)\} \end{array} \right. \qquad [13]$$
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Fonctions de forme Hermite-NEM</figcaption></figure>

qui conduit à la solution exacte suivante :

$$u^h(\mathbf{x}) = \{1 - \cos(2\pi x)\} \{1 - \cos(2\pi y)\} \quad [14]$$

La forme faible associée à l'équation [12] est donnée par :

Trouver $u \in H_0^1(\Omega)$ tel que :

Trouver $u \in H_0^1(\Omega)$ tel que :

$$\int_{\Omega} \nabla u \cdot \nabla \delta u d\Omega = \int_{\Gamma_u} f \delta u d\Gamma, \quad \forall \delta u \in H_0^1(\Omega) \quad [15]$$

où $H_0^1(\Omega)$ est l'espace fonctionnel de Sobolev usuel. L'interpolation Hermite-NEM que nous venons de construire est utilisée pour approximer $\mathbf{u}$ et $\delta\mathbf{u}$, avec la seule contribution des nœuds internes.

L'erreur en norme énergie est donnée par :

$$\|\mathbf{u} - \mathbf{u}^h\|_{E(\Omega)} = \left( \frac{1}{2} \int_{\Omega} (\nabla u^{ex} - \nabla u^h)^T (\nabla u^{ex} - \nabla u^h) \right)^{1/2} \quad [16]$$

Pour l'évaluation des équations [15] et [16], les cellules de Voronoi sont triangulées et un schéma d'intégration de Gauss est appliqué dans chaque triangle, avec 3, 6 et 12 points. La figure 3 compare la convergence de la méthode Hermite-NEM (H-NEM) et de la méthode NEM standard. Un nuage régulier de nœuds a été utilisé et puis raffine par homothétie. Si trois points d'intégration de Gauss sont utilisés, la qualité de la méthode H-NEM est supérieure à celle de la méthode NEM, mais la différence de vitesse de convergence n'est pas significative par rapport à la NEM. Si un schéma d'intégration suffisamment fin est utilisé (6 points), la convergence dans la méthode H-NEM devient quadratique.
<figure><img src="image_4.png" /><figcaption>**Figure 3.** *Convergence en norme énergie pour le problème 2D de Poisson*</figcaption></figure>

**4. Exemples numériques : simulations de procédés de mise en forme**

Dans cette section nous allons réunir quelques exemples d'application des tech-
niques C-NEM et alpha-NEM dans le cadre de la simulation des procédés. Les simu-
lations d'extrusion (figure 4) et d'usinage (figure 5) d'un matériau viscoplastique ont
été réalisées en utilisant la méthode α-NEM, avec le même nuage de nœuds pendant
toute la simulation. A chaque pas de temps les nœuds sont transportés avec la vitesse
matérielle. Le diagramme de Voronoi est alors reconstruit, permettant le calcul des
fonctions de forme dans les points d'intégration qui servent à évaluer la formulation
variationnelle de puissances virtuelles. Comme nous l'avons déjà précisé, la méthode
α-NEM ne nécessite aucune description explicite de la frontière du domaine au cours
du temps. La frontière associée à chaque nuage de nœuds est extraite automatiquement
en faisant appel aux formes alpha. Un avantage supplémentaire de la méthode α-NEM
est la possibilité de gérer l'interaction entre différents fronts fluides qui se rencontrent
comme on peut l'apprécier dans la figure 6. Par contre le problème de l'auto-contact
ou de la gestion de la propagation d'une fissure reste plus délicat.

La simulation du cisaillage adiabatique a été réalisée dans le cadre de la méthode C-NEM couplée avec un indicateur d'erreur de type Z-Z et un algorithme de raffinage adaptatif (Yvonnet et al. 2005b) (voir figure 7). La méthode C-NEM a été choisie par sa capacité à décrire avec précision des interfaces telles que les fissures, rendant donc possible la gestion précise de la séparation de la matière. Par rapport à la méthode α-NEM, la méthode C-NEM nécessite la description explicite du bord du domaine ainsi
que de toutes les interfaces. Cette description est faite à l'aide d'une triangulation du bord du domaine, qui doit être raffinée quand elle subit des déformations importantes.

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Extrusion d'un matériau thermoviscoplastique 3D. (a) Filière ; (b) Distribution de la température [K] ; (c) distribution de la vitesse de déformation équivalente [$s^{-1}$]</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** *Simulation viscoplastique de la coupe (α-NEM)*</figcaption></figure>

Dans tous les cas (α-NEM ou C-NEM) la qualité de la solution est assez insensible vis-à-vis de la qualité des triangles (tétraèdres) de Delaunay, ce qui justifie le caractère « sans maillage » de la méthode des éléments naturels. Pour d'autres applications concernant l'application de la méthode des éléments naturels en mise en forme, le lecteur peut se reporter à (Martínez et al. 2004).
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Extrusion d'un tube mettant en jeu plusieurs fronts de matière</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 7.** Localisation d'une bande de cisaillement (thermo-élasto-plasticité en grandes transformations) (Yvonnet et al. 2005b)</figcaption></figure>
# 5. Conclusions

Dans cet article, quelques exemples de simulation de différents procédés ont été présentés : extrusion, coupe, découpe par cisaillage adiabatique... obtenus en utilisant une méthode de type éléments naturels. Deux techniques permettant de gérer les domaines non convexes ont été présentées. La méthode des formes alpha permet d'extraire automatiquement la frontière du domaine avec la seule information provenant du nuage de nœuds. La méthode C-NEM se base sur l'introduction d'un critère de visibilité permettant de bien représenter les fissures et les interfaces.

Dans le cadre de la simulation des écoulements de fluides visqueux, où des formulations mixtes sont nécessaires pour éviter les problèmes liés à l'incompressibilité, différentes solutions ont été proposées. Une formulation de type Hermite a été également proposée pour obtenir une approximation à consistance d'ordre supérieur, afin d'améliorer la qualité dans la NEM. La validation de cette approche ainsi que la gestion de la propagation des fissures sont, entre autres, des recherches en cours.

# 6. Bibliographie

Belytschko T., Lu Y. Y., Gu L., « Element-free Galerkin methods », *International Journal for Numerical Methods in Engineering*, vol. 37, p. 229-256, 1994.

Cueto E., Sukumar N., Calvo B., Martínez M. A., Cegoñino J., Doblaré M., « Overview and recent developments in Natural Neighbour Galerkin methods », *Archives of Computational Methods in Engineering*, vol. 10 (4), p. 307-384, 2003.

González D., Cueto E., Doblaré M., « Volumetric locking in Natural Neighbour Galerkin methods », *International Journal For Numerical Methods In Engineering*, vol. 61 (4), p. 611-632, 2004.

Liu W. K., Jun S., Zhang Y. F., « Reproducing Kernel Particle Methods », *International Journal for Numerical Methods in fluids*, vol. 21, p. 1081-1106, 1995.

Martínez M. A., Cueto E., Alfaro I., Doblaré M., Chinesta F., « Updated Lagrangian free surface flow simulations with Natural Neighbour Galerkin methods », *International Journal for Numerical Methods in Engineering*, vol. 60(13), p. 2105-2129, 2004.

Piper B., « Properties of local coordinates based on Dirichlet tessellations », *Computer Suppl*, vol. 8, p. 227-239, 1995.

Yvonnet J., Chinesta F., Ryckelynck D., Lorong P., « The Meshless Constrained Natural Element Method (C-NEM) for Treating Thermal Models Involving Moving Interfaces », *International Journal of Thermal Sciences*, vol. 44, p. 559-569, 2005a.

Yvonnet J., Ryckelynck D., Lorong P., Chinesta F., « A New Extension of the Natural Element Method for Non Convex and Discontinuous Problems : the Constrained Natural Element Method », *International Journal for Numerical Methods in Engineering*, vol. 60(8), p. 1451-1474, 2004.

Yvonnet J., Ryckelynck D., Lorong P., Chinesta F., « Simulating thermo-elasto-plasticity in large transformations with adaptive refinement in the natural element method : application to shear banding », *International Journal of Forming Processes*, vol. 8, p. 347-373, 2005b.
Yvonnet J., Villon P., Chinesta F., « Approximations involving bubbles for treating incompres-
sible media », *International Journal for Numerical Methods in Engineering*, n.d.