r3.03.09 Raccord Arlequin 3D – Poutre

Résumé :

Ce document explique la méthode Arlequin développée dans Code_Aster pour relier une modélisation milieu continu 3D et une modélisation poutre.

En 3D, ce raccord se traduit par des relations linéaires reliant les déplacements de l’ensemble de nœuds 3D (3 degrés de liberté par nœud) liés avec l’ensemble des nœuds de poutre (6 degrés de liberté par nœud).

Le raccord Arlequin 3D-poutre

Objectifs et hypothèses

Contrairement aux méthodes adaptatives classiques, le cadre Arlequin permet de superposer des modèles locaux adaptés, plutôt que de corriger les modèles globaux existants. La solution physique est alors reconstituée dans la zone de superposition comme une partition des champs solutions.

Dans ce qui suit, les qualifications « fin » ou « local » et « grossier » ou « global » sont utilisées ici par anticipation des hétérogénéités des modélisations utilisées (poutre 1D et volumique 3D).

La méthode Arlequin est une méthode de couplage avec recouvrement basée sur la combinaison de modèles de finesses et/ou de modélisations différentes. Elle permet le mélange ainsi que le raccord “en volume” des formulations de comportements hétérogènes, et ce sans imposer a priori des contraintes sur les maillages à recoller [1] . Elle est fondée sur deux idées principales :

• Le raccord des sous domaines par l’intermédiaire d’une formulation faible : l’introduction des multiplicateurs de Lagrange dans la zone de collage garantit le couplage des modèles, la continuité des quantités cinématiques, ainsi que le contrôle des écarts des contraintes et des déformations entre les zones couplées;

• La distribution de l’énergie entre domaines et modèles : dans le but de ne pas compter deux fois l’énergie du système global dans la zone de recouvrement, les travaux virtuels associés aux deux modèles sont distribués entre les sous-domaines couplés à travers la zone de collage par le biais de fonctions de pondération \((\gamma ,1-\gamma )\) qui forment une partition de l’unité (la somme des deux fonctions est égale à 1) sur l’ensemble du domaine d’étude. On note que pour deux sous-domaines \({\Omega}_{\mathrm{1D}}\) et \({\Omega}_{\mathrm{3D}}\) , la zone de recouvrementest définie par leur intersection.

Figure 2.1-a fonctions de pondération dans le cadre Arlequin

Cette pondération permet d’éviter de prendre en compte l’énergie plusieurs fois dans la même zone et autorise ainsi une certaine liberté à l’utilisateur pour le choix du modèle prédominant. En effet, elle permet de mettre le poids sur le modèle que l’on souhaite faire exprimer.

Dans le cadre Arlequin, dès que des modélisations hétérogènes sont couplées, il est préconisé de discrétiser les multiplicateurs de Lagrange sur le modèle grossier en espace. Ceci permet d’éviter le phénomène de verrouillage numérique et, par conséquent, d’autoriser plus de liberté à chacun des modèles en permettant au modèle fin d’exprimer sa richesse. Il en résulte que l’espace éléments finis de ces multiplicateurs de Lagrange est construit sur le modèle grossier (modèle de poutre pour ce cas de figure).

Le raccord 3D-Poutre Arlequin, tel que développé dans Code_Aster, doit satisfaire une exigence hautement importante. En effet, les maillages 1D et 3D doivent être « hiérarchiquement compatibles », au sens où tous les éléments 3D sont inclus dans l’espace cylindrique de l’élément 1D en vis-à-vis (pas d’éléments 3D à cheval entre deux éléments poutre). La seule contrainte que l’on impose à l’utilisateur consiste donc à faire un effort lors de l’élaboration du maillage, ce qui n’est pas une mauvaise chose car le modèle final sera mieux maîtrisé.

Figure 2.1-b représentation schématique d’un raccordement d’un modèle mixte 1D-3D de poutre cantilever avec « maillages hiérarchique »

La méthode Arlequin est ici utilisée pour raccorder deux modèles 1D et 3D (opération de jonction). Par conséquent, la zone de couplage correspond à la zone de collage. La gestion des pondérations pour assurer la partition de l’unité dans la zone libre (des conditions Arlequin) est gérée par l’utilisateur.

Comme il s’agit d’un couplage de type L2, il suffit d’imposer les coefficients de pondération dans la définition des matériaux relatifs aux différentes zones libre et de collage.

Remarque :

Une procédure générale pour gérer les maillages incompatibles serait très lourde et coûteuse à programmer, elle nécessiterait une gestion compliquée des fonctions de pondération dans la zone de collage ainsi qu’une intersection ou un découpage en sous-éléments des maillages 1D et 3D pour gérer l’intégration des termes de couplage.

Notations

On notera :

• \(\Omega\) la zone de couplage

• \({\mathrm{u}}_{\mathrm{1D}}\) et \({\mathrm{u}}_{\mathrm{3D}}\) les champs de déplacements correspondant aux modèles 1D poutre et 3D

• \(\lambda\) le multiplicateur de Lagrange associé au domaine 1D

• \(\text{M}\) la matrice classique de raideur poutre

Expression de la condition de raccord Arlequin

L’équation de couplage L2 en milieu continu s’écrit :

\({\int}_{\Omega}\lambda :({\mathrm{u}}_{\mathrm{1D}}-{\mathrm{u}}_{\mathrm{3D}})\text{d}\Omega ={\int}_{\Omega}\lambda :{\mathrm{u}}_{\mathrm{1D}}\text{d}\Omega -{\int}_{\Omega}\lambda :{\mathrm{u}}_{\mathrm{3D}}\text{d}\Omega\) éq 2.4-1

\(\lambda\) et \({\mathrm{u}}_{\mathrm{1D}}\) ont exactement le même type de degrés de liberté et le même type de fonctions de forme car ils sont définis dans le même espace des champs cinématiquement admissibles. En effet, pour un élément 1D du multiplicateur de Lagrange correspond un élément 1D poutre.

On suppose, pour simplifier, que l’axe de la poutre est selon X.

Figure 2.3-a représentation schématique d’un élément de poutre à 6 degrés de liberté

Alors, le déplacement de type poutre s’écrit classiquement comme suit :

\(\lbrace \begin{array}{cccc}u=u(x)& +z{\theta}_{y}(x)& -y{\theta}_{z}(x)& +0\\ v=0& +0& +v(x)& -(z-{z}_{c}){\theta}_{x}(x)\\ w=0& +w(x)& +0& +(y-{y}_{c}){\theta}_{x}(x)\end{array}\)

Le multiplicateur de Lagrange s’écrira alors selon une cinématique de poutre comme suit :

\(\lbrace \begin{array}{cccc}{\lambda}_{u}={\lambda}_{u}(x)& +z{\lambda}_{{\theta}_{y}}(x)& -y{\lambda}_{{\theta}_{z}}(x)& +0\\ {\lambda}_{v}=0& +0& +{\lambda}_{v}(x)& -(z-{z}_{c}){\lambda}_{{\theta}_{x}}(x)\\ {\lambda}_{w}=0& +{\lambda}_{w}(x)& +0& +(y-{y}_{c}){\lambda}_{{\theta}_{x}}(x)\end{array}\)

Sur un élément 1D noté e1D le premier terme de couplage 1D-1D s’écrit dans le même espace 1D :

\({\int}_{{\text{e}}_{\mathrm{1D}}}\lambda :{\mathrm{u}}_{\mathrm{1D}}\text{d}\Omega ={\left\lbrace \lambda \right\rbrace }^{t}{\text{C}}_{\mathrm{1D}-\mathrm{1D}}\left\lbrace {u}_{\mathrm{1D}}\right\rbrace ={\left\lbrace \lambda \right\rbrace }^{t}\text{M}\left\lbrace {u}_{\mathrm{1D}}\right\rbrace\) éq 2.4-2

où \(\left\lbrace \lambda \right\rbrace\) est le vecteur des 12 degrés de liberté du multiplicateur de Lagrange et \(\left\lbrace {u}_{\mathrm{1D}}\right\rbrace\) est le vecteur des 12 degrés de liberté des éléments poutre :

\(\left\lbrace {u}_{\mathrm{1D}}\right\rbrace =({u}^{1},{v}^{1},{w}^{1},{{\theta}_{x}}^{1},{{\theta}_{y}}^{1},{{\theta}_{z}}^{1},{u}^{2},{v}^{2},{w}^{2},{{\theta}_{x}}^{2},{{\theta}_{y}}^{2},{{\theta}_{z}}^{2})\)

\(\left\lbrace \lambda \right\rbrace =({\lambda}_{u}^{1},{\lambda}_{v}^{1},{\lambda}_{w}^{1},{\lambda}_{{\theta}_{x}}^{1},{\lambda}_{{\theta}_{y}}^{1},{\lambda}_{{\theta}_{z}}^{1},{\lambda}_{u}^{2},{\lambda}_{v}^{2},{\lambda}_{w}^{2},{\lambda}_{{\theta}_{x}}^{2},{\lambda}_{{\theta}_{y}}^{2},{\lambda}_{{\theta}_{z}}^{2})\)

Dans le cas particulier du couplage L2, la matrice \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}\) est égale à la matrice classique de masse d’un élément de poutre.

Toujours sur un élément 1D noté e1D ,le deuxième terme de couplage 1D-3D s’écrit dans le même espace 1D :

\(-{\int}_{{\text{e}}_{\mathrm{1D}}}\lambda :{\mathrm{u}}_{\mathrm{3D}}\text{d}\Omega =-{\left\lbrace \lambda \right\rbrace }^{t}{\text{C}}_{\mathrm{1D}-\mathrm{3D}}\left\lbrace {u}_{\mathrm{3D}}\right\rbrace\) éq 2.4-3

et \(\left\lbrace {u}_{\mathrm{3D}}\right\rbrace\) est le vecteur des (3 x \({N}_{\mathit{no3D}}\) ) degrés de liberté des éléments 3D (\({N}_{\mathit{no3D}}\) étant le nombre de nœuds de l’élément 3D) :

\(\left\lbrace {u}_{\mathrm{3D}}\right\rbrace =({u}^{1},{v}^{1},{w}^{1},\text{...},{u}^{{N}_{\mathit{no3D}}},{v}^{{N}_{\mathit{no3D}}},{w}^{{N}_{\mathit{no3D}}})\)

La démarche de construction de la matrice \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}\) sur le modèle fin 3D est plus délicate que celle du modèle grossier de Poutre. La complexité de ce couplage réside dans la coexistence de deux champs appartenant à des espaces différents. Elles’obtient en considérant tous les éléments 3D contenus dans l’espace cylindrique de l’élément 1D.

\({\text{C}}_{\mathrm{1D}-\mathrm{3D}}=\sum_{k=1}^{{N}_{\mathit{el3D}}}{\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}\)

Il est rappelé que, par construction même du maillage, tous les éléments 3D doivent être inclus dan cet espace. Ces éléments sont déterminés par une procédure d’appariement (graphe de parenté) qui consiste à projeter leurs barycentres dans le support géométrique de l’élément 1D associé (projection d’un point sur un segment).

Pour chaque élément 3D numéro k (k allant de 1 à Nel3D), on calcule la matrice \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}\) de la manière suivante :

\({\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}=\sum_{g=1}^{\mathit{Nbgauss}}{\mathrm{W}}^{g}{\mathrm{B}}_{\mathrm{1D}}^{gt}({\xi}_{\mathrm{1D}}){\mathrm{G}}^{g}{\mathrm{B}}_{\mathrm{3D}}^{g}({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\det({\mathrm{J}}^{\mathrm{g}})\)

avec G la matrice de passage Élément de référence – Élément physique, contenant des zéros et des éléments de l’inverse du jacobien calculés au point de gauss g.

La matrice \({\mathrm{B}}_{\mathrm{1D}}^{g}\) est de dimension 3 x (3 x nombre de nœuds de l’élément 1D). Ses termes sont calculés aux points de Gauss des éléments 1D.

\({\mathrm{B}}_{\mathrm{1D}}^{g}=\left[\begin{array}{cccccccccccc}{N}_{1}^{p}& y{\xi}_{5}& z{\xi}_{5}& 0& z{\xi}_{6}& -y{\xi}_{6}& {N}_{2}^{p}& y{\xi}_{7}& z{\xi}_{7}& 0& z{\xi}_{8}& -y{\xi}_{8}\\ 0& {\xi}_{1}& 0& -z{N}_{1}^{p}& 0& -{\xi}_{2}& 0& {\xi}_{3}& 0& -z{N}_{2}^{p}& 0& -{\xi}_{4}\\ 0& 0& {\xi}_{1}& y{N}_{1}^{p}& {\xi}_{2}& 0& 0& 0& {\xi}_{3}& y{N}_{2}^{p}& {\xi}_{4}& 0\end{array}\right]\)

où x, y et z sont les variables locales de l’élément 1D, \({N}_{1}^{p}\) et \({N}_{2}^{p}\) les fonctions de forme pour la traction-compression-torsion, \({\xi}_{\mathrm{1...8}}\) sont les fonctions de forme de flexion de la poutre (cf. document [r3.08.01] de Code-Aster).

La matrice \({\mathrm{B}}_{\mathrm{3D}}^{g}\) est de dimension 3 x (3 x nombre de nœuds de l’élément 3D). Par exemple, la matrice \({\mathrm{B}}_{\mathrm{3D}}^{g}\) associée à l’élément volumique linéaire tétraèdre (TETRA4)s’écrit :

\({\mathrm{B}}_{\mathrm{3D}}^{g}=\left[\begin{array}{cccccccccccc}{N}_{1}& 0& 0& {N}_{2}& 0& 0& {N}_{3}& 0& 0& {N}_{4}& 0& 0\\ 0& {N}_{1}& 0& 0& {N}_{2}& 0& 0& {N}_{3}& 0& 0& {N}_{4}& 0\\ 0& 0& {N}_{1}& 0& 0& {N}_{2}& 0& 0& {N}_{3}& 0& 0& {N}_{4}\end{array}\right]\)

Les variables x, y, z ainsi que les fonctions de forme et leurs dérivées sont calculées au point de Gauss g de l’élément 3D. Ainsi, par exemple, on a :

\(\begin{array}{c}x=\langle N\rangle \left\lbrace {x}_{n}\right\rbrace ,{u}_{x}=\langle N\rangle \left\lbrace {u}_{n}\right\rbrace \\ y=\langle N\rangle \left\lbrace {y}_{n}\right\rbrace ,{u}_{y}=\langle N\rangle \left\lbrace {v}_{n}\right\rbrace \\ z=\langle N\rangle \left\lbrace {z}_{n}\right\rbrace ,{u}_{z}=\langle N\rangle \left\lbrace {w}_{n}\right\rbrace \end{array}\)

où :

• xn, ynet znsont les coordonnées des nœuds de l’élément 3D,

• un, vnet wn sont les degrés de liberté associés aux nœuds des élément 3D,

• :math:`langle Nrangle ` est le vecteur contenant toutes les fonctions de forme de l’élément 3D suivant l’ordre des nœuds 3D (cf. document r3.01.01 de Code-Aster).

Pour chaque point de Gauss de coordonnées \(({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\) dans l’élément de référence 3D, nous introduisons sa représentation physique sous forme des coordonnées interpolées \((\stackrel{ˉ}{X},\stackrel{ˉ}{Y},\stackrel{ˉ}{Z})\) suivantes :

\(\stackrel{ˉ}{X}=\sum_{i=1}^{{N}_{\mathit{el3D}}}{N}^{i}{x}^{i}({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\)

\(\stackrel{ˉ}{Y}=\sum_{i=1}^{{N}_{\mathit{el3D}}}{N}^{i}{y}^{i}({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\)

\(\stackrel{ˉ}{Z}=\sum_{i=1}^{{N}_{\mathit{el3D}}}{N}^{i}{z}^{i}({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\)

Figure 2.3-b opération de localisation des points de Gauss d’un élément 3D dans l’unique élément de poutre qui lui est apparié

Pour simplifier l’illustration, nous supposons que l’axe de la poutre est selon X. Pour chaque point de Gauss de l’élément 3D, nous introduisons sa projection sur l’unique élément réel 1D qui lui est apparié comme suit :

\({\stackrel{ˉ}{\xi}}_{\mathrm{1D}}=\frac{\stackrel{ˉ}{X}-{X}_{inf}}{\text{Le}}\)

Les coordonnées de la projection dans les deux autres directions de l’épaisseur de la poutre sont données par \(y=\stackrel{ˉ}{Y}\) et \(z=\stackrel{ˉ}{Z}\) .

Cette technique permet de discrétiser le multiplicateur de Lagrange sur les éléments de poutre et de calculer les intégrales dans les expressions des matrices de couplage \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}\) . Ainsi, pour chaque élément 3D on calcule la matrice \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}\) comme suit :

\({\text{C}}_{\mathrm{1D}-\mathrm{3D}}^{k}=\sum_{g=1}^{\mathit{Nbgauss}}{\mathrm{W}}^{g}{\mathrm{B}}_{\mathrm{1D}}^{gt}({\stackrel{ˉ}{\xi}}_{\mathrm{1D}}){\mathrm{G}}^{g}{\mathrm{B}}_{\mathrm{3D}}^{g}({\xi}_{\mathrm{3D}},{\eta}_{\mathrm{3D}},{\zeta}_{\mathrm{3D}})\det({\mathrm{J}}^{\mathrm{g}})\)

où les fonctions de forme de type poutre sont interpolées aux coordonnées paramétriques des points projetés obtenus.

Implantation de la méthode de raccord

Pour chaque raccord Arlequin, l’utilisateur doit définir les opérandes suivantes sous l’option 3D_POU_ARLEQUIN du mot clé LIAISON_ELEM:

Le groupe de mailles 3D de la zone de collage (mot clé GROUP_MA_1)

Le groupe de mailles 1D de la zone de collage (mot clé GROUP_MA_2)

Le concept CHAM_MATER définissant les matériaux (servant à assurer la partition de l’unité)

Le concept CARA_ELEM définissant les caractéristiques élémentaires servant au calcul des matrices de couplage \({\text{C}}_{\mathrm{1D}-\mathrm{1D}}\) et \({\text{C}}_{\mathrm{1D}-\mathrm{3D}}\)

Pour chaque élément poutre, Code_aster calcule les relations linéaires issues de l’écriture des matrices de couplage assemblées [éq 2.4-1]. Ces relations cinématiques relient :

• les 12 degrés de liberté des deux nœuds de l’élément poutre,

• avec les degrés de liberté des nœuds de tous les éléments appariés à l’élément poutre.

Enfin, les relations linéaires résultantes seront dualisées, comme toutes les relations linéaires issues par exemple du mot clé LIAISON_DDL de AFFE_CHAR_MECA.

[1]

En pratique, on demandera à l’utilisateur de Code_Aster de fournir des maillages hiérarchiques.

Quelles utilisations peut-on faire de cette modélisation ?

Actuellement, il est possible d’envisager trois types d’utilisation de la méthode Arlequin :

• Le zoom : Ce type d’application vise la détection et l’étude d’un phénomène local (contact, fissure,…) au niveau d’une zone précise appartenant à un domaine global. Un patch fin est alors superposé au modèle grossier existant, dans le but de raffiner la solution aux voisinages des zones d’intérêt. Dans ce cas de figure, la zone de collage est égale à un des sous-domaines.

• La jonction : Le but principal visé par cette application consiste à relier deux modèles qui peuvent être de natures et/ou de finesses différentes. Dans cette configuration, la zone de collage est strictement égale à la zone de recouvrement des deux domaines.

• La substitution : Cette opération consiste à substituer localement un modèle existant par un autre plus fin. Cette technique s’avère très utile pour introduire avec une grande flexibilité, un défaut dans un modèle sain tel que la fissuration. Contrairement aux deux cas précédents, en appliquant la substitution, la zone de collage constitue un sous-domaine de la zone de recouvrement des deux modèles.

Figure 3-a différentes utilisations du raccord 1D-3D dans le cadre Arlequin

Bibliographie

1. 1. Ghanem, «Contribution à la modélisation avancée des machines tournantes en dynamique transitoire dans le cadre Arlequin», thèse de l’INSA de Lyon.

Description des versions du document

Version Aster	Auteur(s) Organisme(s)	Description des modifications
12	M. TORKHANI (EDF/R&D/AMA)	Texte initial