La méthode Area In Contact associe une aire à chaque noeud esclave (élément de contact). Ainsi, on peut tenir compte de la densité du maillage sur la frontière de contact (petite maille à côté de grosse maille). Dès lors, la force de contact agissant sur le noeud esclave, dans le cas d'un contact sans frottement (pour faire simple), est donnée par

$\vec{F}^{(s)} = - t_n A_s \vec{n}^{(m)}$

où

- (s) indique l'esclave.
- (m) indique le maître.
- $A_s$ est l'aire de contact nodale.
- $\vec{n}$ est la normale extérieure unitaire.
- $t_n$ est la pression de contact.

Dans le cadre de la méthode de la pénalité classique, la pression de contact est donnée par

$t_n = C_n g_n $

où

- $C_n$ est le coefficient de pénalité.
- $g_n$ est la distance signée du noeud esclave par rapport à la surface maître.

Lors de la définition d'une interaction de contact, nous devons spécifier l'entité géométrique esclave sur laquelle les éléments de contact sont générés.

Lors de l'appel à la fonction toDofSet() dans la classe ContactInteraction, nous remplissons une table de hachage contenant l'adresse mémoire de tous les noeuds de l'entité géométrique esclave, voire maître dans le cas du contact déformable-déformable avec double passe.

Ainsi, le test suivant est ajouté aux tests existants pour corriger les effets de bord indésirables. Lors l'appel à la fonction computeSlaveArea dans la classe patron ContactElement, pour chaque noeud esclave (i.e l'élément de contact), nous considérons que l'entité topologique du maillage incidente à ce noeud (Curve en 2D et Side en 3D) contribue à l'aire de contact nodale à condition que l'ensemble des noeuds appartenant à cet entité géométrique (2 noeuds en 2D et 4 noeuds en 3D) possède une occurrence dans la table de hachage, c'est à dire que la “Curve” en 2D ou la “Side” en 3D est contenue entièrement dans l'entité géométrique esclave sur laquelle est définie l'interaction de contact.

Finalement, la correction des effets de bord fonctionne quelque soit le type de l'interaction de contact (auto-contact, contact rigide déformable, contact déformable-déformable avec simple ou double passe)

Jusqu'à présent, il faut définir dans l'ElementProperties des éléments de contact, le paramètre AREAINCONTACT :

    prp1 = ElementProperties(Contact[2|3]DElement)   
    prp1.put(AREAINCONTACT,AIC_ONCE)

Si le paramètre AREAINCONTACT vaut AIC_ONCE, l'aire de contact nodale est calculée sur la configuration initiale et n'est plus mise à jour au cours de l'intégration temporelle. Cependant, si le paramètre AREAINCONTACT vaut AIC_ONCEPERSTEP, l'aire de contact nodale est calculée sur la configuration courante obtenue au pas de temps précédent = la dernière configuration équilibrée.

Cependant, il est nécessaire de clarifier les situations complexes suivantes :

• Rupture d'un élément :

Lors de l'appel à la fonction checkFreeNode, nous devons recalculer l'aire de contact nodale sur la configuration initiale, si l'option AREAINCONTACT est égale AIC_ONCE; ou sur la dernière configuration équilibrée, si l'option AREAINCONTACT est égale AIC_ONCEPERSTEP, lors de la rupture d'un élément fini volumique.

• Formalisme ALE :

Lors de l'appel à la fonction endALE, on recalcule l'aire de contact nodale sur la dernière configuration équilibrée uniquement si l'option AREAINCONTACT est égale AIC_ONCEPERSTEP.

Il faut définir dans l'ElementProperties des éléments de contact, le paramètre AREAINCONTACTMETHOD :

    prp1 = ElementProperties(Contact[2|3]DElement)   
    prp1.put(AREAINCONTACTMETHOD, AIC_METHOD_CONSISTENT)

Si le paramètre AREAINCONTACTMETHOD vaut AIC_METHOD_GEOMETRIC, le calcul de l'aire de contact nodale est basée sur une approche purement géométrique. Cependant, si le paramètre AREAINCONTACTMETHOD vaut AIC_METHOD_CONSISTENT, le calcul de l'aire de contact nodale est basée sur une approche énergétiquement consistante, c'est à dire sur base du principe des travaux virtuels (ici le travail virtuel des forces de contact).

Cette dernière méthode permet de retrouver une pression de contact nodale uniforme lorsque l'on écrase un disque creux/plein, un bloc en état plan déformation ou en 3D au moyen d'un plan rigide. Attention que la force de contact n'est pas uniforme sur toute la surface de contact !!!

Sur les figures, ci-dessus, vous pouvez voir la différence entre la méthode géométrique ou consistante dans le cas de l'écrasement d'un bloc en 3D au moyen d'un plan rigide :

CONFIGURATION INITIALE AIC_METHOD_CONSISTENT AIC_METHOD_GEOMETRIC

Les extracteurs de grandeur d'une interaction de contact possèdent un filtre qui permet d'archiver les valeurs du nom de l'extracteur en un sous-ensemble de noeuds esclave (éléments de contact) de l'interaction de contact.

valueExtractor = ContactStatusValueExtractor(contInt, sOp=None, maxV=-1)
valueExtractor.setGeoFilter(gObject)
gObject = Référence vers un objet géométrique maillé, qui est inclu dans l'objet géométrique sur lequel l'interaction de contact a été définie.  

Par défaut, nous considérons tous les éléments de contact de l'interaction de contact (Dans le cas du contact déformable-déformable en double passe, nous avons les éléments de contact du maître et de l'esclave.).

Les opérations se font dans l'ordre suivante :

- Opération de filtrage
- Opération de tri
- Opération de “resize” c'est à dire le nombre de valeur filtrée puis triée, retenue.

Il est possible d'extraire l'aire de contact nodale :

areaInContactValues=AreaInContactValueExtractor(contInt, sOp=None, maxV=-1)

De là, il est possible de déduire l'aire de contact (potentiellement en contact) :

valuesmanager.add(no, areaInContactValues, SumOperator(), 'AreaInContact')

Il est possible d'extraire de l'aire de contact nodale des noeuds eslave en contact :

activeAreaInContactValues=AreaInContactValueExtractor(contInt, sOp=None, maxV=-1)
activeAreaInContactValues.setOnlyInContactStatus()

De là, il est possible de déduire l'aire de contact active (réellement en contact) :

valuesmanager.add(no, activeAreaInContactValues, SumOperator(), 'ActiveAreaInContact')

Dans le fichier source “mtGeoBilinearCoons.cpp”, il y a une erreur dans le calcul de l'aire à la ligne 71 et 72 (Fonction “area”):

    Vect3 sum2 ( x3+x4-x2-x1 );
          sum2/=4.0;

J'ai changé sum1/=4.0 par sum2/=4.0.

Dans le cas test apps_imp_revolutionSurfaceTRSym, il y a une erreur sur l'affichage d'une partie des surfaces de révolution. Sous les recommandations de Luc, j'ai ajouté le gardien suivant :

    if(!wire.getSens(ns))
        ksi = 1.0-ksi;

dans le fichier source SurfaceCloud.cpp à la ligne 163 (Fonction “getKsi”).

Les effets de bord indésirables en 2D et en 3D sont corrigés désormais dans la version de Metafor pour le calcul de l'aire de contact nodale. Cependant, en discutant avec luc et romain, j'ai adapté les cas tests suivant pour étendre l'entité géométrique esclave de l'interaction de contact rigide-déformable afin de reproduire au mieux les résultats dans la batterie. Ces cas tests concernent l'interaction de contact entre une aube déformable et le carter rigide. La zone potentiellement en contact de l'aube est désormais son sommet et une couche de noeuds sur sa surface latérale.

Pour ces cas tests, j'ai modifié le fichier blisk.py dans newac.tools.snecma. Dans ces cas tests, il y a une ou trois aubes déformable(s) en contact avec le carter rigide. Ces interactions de contact sont définies à l'aide des groupes de noeuds 2, 3 et 4. Les noeuds du groupe sur lequel est défini l'interaction de contact rigide déformable sont sélectionnées en fonction de leur position relative par rapport à une boîte

    xmin = tpointset(1113).getPos0().get1()
    ymin = tpointset(840).getPos0().get2()
    zmin = tpointset(1113).getPos0().get3()
    xmax = tpointset(840).getPos0().get1()
    ymax = tpointset(1113).getPos0().get2()
    zmax = tpointset(840).getPos0().get3()
    groupset(2).addMeshPointsInBox(xmin,xmin+(xmax-xmin)/2.0,ymin+(ymax-ymin)/2.0,ymax+1.0,zmin,zmin+2.0*(zmax-zmin)/3.0)
    groupset(2).addMeshPointsInBox(xmin+(xmax-xmin)/2.0,xmax,ymin,ymax+1.0,zmin+(zmax-zmin)/2.0,zmax)

(lignes de code entre les lignes 267 à 274 du fichier blisk.py) construite sur base la position de deux noeuds du maillage (ici les noeuds 840 et 1113). L'un se trouve sur le bord d'attaque de l'aube et l'autre se trouve sur le bord de fuite de l'aube.

Pour chaque groupe de noeud, j'ai dû aller chercher le numéro du noeud du maillage situé sur le bord d'attaque de l'aube et sur le bord de fuite de l'aube de la couche de noeuds localisée juste en dessous de la couche de noeuds située au sommet de l'aube, à l'aide de la fenêtre de visualisation. Finalement, j'ai complété le groupe de noeuds au dupliquant les commandes utilisées entre la ligne 267 à 274 et en adaptant juste les numéros des noeuds en question.

Finalement, l'interaction de contact collant déformable/déformable entre le blisk et les aubes a été adaptée : le coefficient de pénalité a été augmenté pour compenser le fait que l'aire de contact nodale a diminué en moyenne, lorsque nous ne considérons pas les effets de bords (Cf. commentaires à la ligne 34 du fichier blisk.py)

Pour ces cas tests, j'ai modifié le fichier aube3D18ER.py dans newac.tools.techspace. Les noeuds du groupe sur lequel est définie l'interaction rigide déformable de contact sont sélectionnés en fonction de leur position relative par rapport à un plan

 groupset(1002).addMeshPointsFromObject(groupset(1000),PlaneSelector(sur400,0)

(ligne de code à la ligne 314 du fichier aube3D18ER.py).

Entre les lignes 244 à 268 du fichier aube3D18ER.py, j'ai ajouté la variable dispPlaneContactOffSet pour chaque maillage bacon importé, afin de décaler le plan, i.e. diminuer la variable dispPlaneContact de dispPlaneContactOffSet, de manière à sélectionner une seule couche supplémentaire de noeuds au plus des noeuds au sommet de l'aube.

Pour ces cas tests, j'ai modifié à l'aide de Samcef for machining, les groupes de noeuds correspondant au sommet de l'aube pour chacun des maillages bacon importés (c'est à dire tous les fichiers _sommet.sel dans techspace.meshes), et j'ai modifié aussi dans le fichier meshGDTech.dat le groupe de noeuds correspondant au sommet de l'aube.

Pour ce cas test, avec l'aide de Yves, j'ai reconstruit la géométrie après l'extrusion du maillage pour que l'interaction de contact se fasse entre l'outil rigide et uniquement la face supérieure du maillage.

Ajoutez mes cas tests dans la batterie

[r]:
[a]:oo_meta\mtElements\DdContactInteraction.inl
[a]:oo_meta\mtElements\ContactInteraction.inl
[a]:oo_meta\mtElements\boundaries\AreaInContactValueExtractor.h
[a]:oo_meta\mtElements\boundaries\AreaInContactValueExtractor.cpp

[r]:
[a]:

— Gaëtan WAUTELET 2013/01/24