Commit 2011-07-25

Commit des premiers développements ThixoWall permettant de s’affranchir de la description topologique des frontières lors d’un calcul de mise en forme en grandes déformations.

La frontière de la structure est décrite par l’iso-zéro d’une level set (fonction distance signée - cf l’exemple si dessous pour une frontière circulaire). Dans les zones ou la level set est négative, on suppose qu’on est dans la matière, inversement si elle est positive on est dans le vide. La formulation éléments finis doit donc être enrichie de manière à pouvoir considérer les cas ou l’élément est vide, rempli, ou partiellement rempli. Les deux premiers cas sont triviaux : on les enlève de l’analyse ou on utilise une formulation classique sur l’élément. Pour le dernier cas, on va intégrer par partie les forces internes et la matrice de raideur en utilisant des sous-cellules triangulaires d’intégration. Ces sous-cellules n’ont pas de degrés de liberté associés. Ce sont uniquement des supports des quantités du quadrangle sur les section de l’élément ou il existe de la matière.

On introduit deux nouvelles dll :

• mtXFEM : cette dll contient les 3 classes utilisées lors d’un calcul XFEM :
  • L’élément XFEM (XFEMElement) est un élément fini de type quadrangulaire qui possède les informations permettant de décrire la frontière discrétisée ainsi que les sous-cellules triangulaires utilisées pour l’intégration.
  • La méthode d’intégration (XFEMVolumeIntegrationMethod) permet la prise en compte des parties « vide » et « matière » sur un élément ainsi que le calcul des forces internes et de la matrice de raideur en fonction du positionnement de la frontière.
  • La nouvelle « interaction » (XFEMFieldApplicator) permet de générer des éléments XFEM et des sous-cellules d’intégration. Elle contient également les informations nécessaires à l’affichage.
• mtXFEMDrawable : cette dll contient uniquement le drawable spécifique aux éléments XFEM (XFEMFieldApplicatorDrawable). Il est ajouté dans une dll indépendante de mtDrawable et de mtXFEM car :
  • Si on met le drawable des XFEM directement dans le mtXFEM, on ne pourra plus faire de XFEM sans visualisation
  • Si on met le drawable XFEM dans mtDrawable, la lib XFEM sera systématiquement appelée dès qu’on fera de la visualisation

Pour l’élément XFEM, j’ai divisé les sources en 3 :

• _Basics contient les fonctions de mise en donnée de l’élément classiques auquelles doit répondre une classe de type élément.
• _Enrichment contient tout ce qui concerne la level set et les sous-cellules d’intégration.
• _io regroupe les fonctions utilisées pour l’affichage.

Plus précisément, l’élément XFEM est définit un peu comme un quadrangle classique mais on lui associe quelques données supplémentaires :

• Une fonction d’interpolation XFEM (utilisé pour interpoler les données des sous-cellules d’intégration sur l’élément)
• Une variable booléenne pour savoir si c’est un élément enrichi (si il est vide on utilise le setEnable(false) utilisé en rupture, si il est plein il a un comportement équivalent à celui d’un Volume2DElement)
• Des std::vector de pointeurs vers
  • les points d’intersection entre la level set et les frontières de l’élément
  • les sides définissant les sous cellules d’integration
  • les points définis de manière découplée des éléments (pour mettre à jour quand il faut la position des points d’intersections et du centre de gravité de la zone matière)
• La fonction distance signée est stockée dans Field(HX)

… et bien sûr on définit dans l’élément quelques fonctions pour calculer tout ca…

Le XFEMFieldApplicator a été définit à l’origine pour avoir un conteneur de mes sous cellules d’intégrations (getXFEMintegrationCells(), qui permet un affichage via win.add(app.getXFEMintegrationCells().getWireSet())). Ces sous-cellules sont vues comme des TriangleSide, et associés en une géométrie. La fonction getXFEMElementSet() regroupe les éléments non-découpés et substitue les éléments triangulaires associés aux sous-cellules d’intégration aux éléments quadrangulaires découpés. Cet ElementSet est utilisé pour l’affichage. On ne peut pas sur-définir l’ElementSet classique car celui-ci est utilisé au moment de l’assemblage (et donc ne doit pas voir les éléments triangulaires). Le XFEMFieldApplicator contient également les fonctions qui permettent de définir les sides triangulaires et les éléments associés a partir des coordonnées des nœuds et des points d’intersection. Les sous-cellules sont systématiquement crées autour du centre de gravité de l’élément :

Une dernière fonction « initIteration » est utilisée dans le XFEMFieldApplicator pour la mise a jour des positions des points d’intersection et du centre de gravité, qui ne sont pas mis à jour au cours des itérations mécanique quand on calcule l’équilibre. Elle remplace dans le fieldApplicator et les classes héritées la fonction « contactDetection » qui est appelée à chaque itération mécanique.

La méthode d’intégration permet de calculer les forces internes et la matrice de raideur en utilisant les sous-cellules d’intégration. Il faut donc connaitre les déformations et les contraintes sur les sous-cellules. Dans un premier temps je les calcule par sous-cellule en utilisant la déformation des sous-cellules d’intégration (il faudrait surement les calculer au niveau de l’élément pour bien faire, en ramenant le point de Gauss du triangle dans le quadrangle). Une fois qu’elles sont connues sur chaque sous cellule, on utilise une intégration « par parties ». La matrice gradient sur le quadrangle est calculée en chaque point de Gauss des sous-cellules en utilisant les fonctions de forme du quadrangle. De même que la matrice jacobienne de la transformation est calculée sur les points de Gauss à partir des données sur le quadrangle.

Le drawable XFEM est assez simple: il parcourt l’ElementSet getXFEMElementSet() au lieu du getElementSet() classique. En présence d’un élément enrichi, il trace donc les données connues sur les sous-cellules d’intégration. Si on trace l’élement set classique, les élements qui sont enrichis ne possèdent aucune donnée en propre (contrainte et déformations sont naturellement nulles car calculées sur les sous-cellules).

En entête de fichier, rajouter :

from wrap.mtXFEM import *
try:
    from wrap.mtXFEMDrawables import *  #importation de la vizu si on est en GUI
except:
    pass

Dans ElementProperties, definir un XFEMElement. Les propriétés additionnelles sont :

• BOUNDARY_WIRE : Numéro de la Wire décrivant la frontière. Dans un premier temps on ne peux en définir qu’une seule. Il faudrait voir pour en définir plusieurs et gérer les inclusions multiples. Valeur par défaut 0 (dans ce cas, on utilise systématiquement l’intégration XFEM sur un élement quadrangulaire complet découpé en 3 sous-cellules).
• NB_PG_INTCELL : Nombre de points de Gauss par sous-cellule triangulaire d’intégration. Valeur par défaut 1.
• INVERSE_WIRE : Sens de parcours de la wire (-1 si la wire est définie avec le vide « a gauche »). Valeur par défaut (-1).
• TOL_LEVELSET : Définit la distance minimale (en pourcentage) à partir de laquelle on considère la distance comme étant nulle (si d < %.sqrt(Aire_Quadrangle), d = 0). Valeur par défaut 0.15.

prp1 = ElementProperties(XFEMElement)
prp1.put(BOUNDARY_WIRE, 3)
prp1.put(NB_PG_INTCELL, 4)
prp1.put(INVERSE_WIRE, -1)
prp1.put(TOL_LEVELSET, 0.15)

Enfin, le FieldApplicator doit être remplacé par un XFEMFieldApplicator :

 app = XFEMFieldApplicator(1) 

Si tout se passe bien, ca devrait fonctionner. A noter que ca ne fonctionne pour le moment pas dans toutes les configurations (les élements trop distordus limitent la convergence) et qu’il faut souvent baisser la précision du solveur itératif (pour le moment, ca passe sur les cas test avec 1e-2).

Les hypothèses de calcul sont dans un premier temps :

• Quasi-statique (les matrices de masse ne sont pas implémentées)
• Pas d’axisymétrique
• Matériau élastique (ca doit marcher avec un autre type de matériau théoriquement, mais l’approximation pour le calcul des variables au points de gauss risque d’être bof).
• Methode d’intégration de Cauchy sur les sous-cellules. On peut par contre choisir n’importe laquelle dans le ElementProperties pour les éléments non sous-découpés.
• Pas de conditions limites sur la wire de la frontière (tous les nœuds supportant une condition limite doivent appartenir au maillage quadrangulaire).

• Changement du CMakeList de MtElement pour que VolumeIntegrationMethods apparaisse dans les méthodes et pas dans les éléments au niveau du projet.
• Correction d’un bug qui faisait que quelque soit le nombre de points de Gauss demandés pour des éléments triangulaires, on en avait qu’un seul.
• Correction des points de Gauss des triangles dans optimalGP (le poids était faux dans le cas de 6 points de Gauss) et ajout des options 7 et 12 GP.
• Changement du allocate des méthodes d’intégration volumique. On passe désormais un « const NbIntegrationPoints » plutôt qu’un « NbIntegrationPoints ». Ca évite la copie et il y avait un souci précédemment.
• contactDetection() utilisé pour la mise à jour du contact est renommé en initIteration() de manière à ce que j’utilise moi même la fonction avec un nom cohérent (pour permettre de retrouver facilement contactDetection pour les habitués, j’ai rajouté un commentaire sur chaque ligne remplacée… Je suis pas sur que j’aurais du le remplacer partout, mais ca permet d’avoir un truc homogène)
• J'ai modifié les elements XFEM de Lara (commenté les lignes qui semblaient inutiles(au cas ou elles l'auraient pas été, j'ai préféré les laisser ;), définit proprement les private/public (proprement pour le 2D, rapidement pour le 3D))

• la fonction getInteraction de interactionDrawable et la fonction generateElements de interaction.h deviennent virtuelles pour que je puisse les sur-définir.
• Cloud.h, ElementCloud.h sont maintenant exportés sous windows à l’appel de la dll. De même pour VolumeElement
• Ajout de la fonction « getGKstate() » dans les méthodes d’intégration volumique (j’en avais besoin pour mes subcells et les GKstates sont privés)
• Définition de deux fonctions forceDefine et forcePush_back dans l’elementSet, car j’avais des problèmes avec le « if(!ownsElems) » à la creation de mon XFEMElementSet.
• Définition d’une fonction projectionOf dans mtGeoCurve qui calcule le point le plus proche sur une curve d’un item de coordonnées données (renvoie le projeté si il est sur la curve, et une des extrémités sinon)

Les résultats de compilation de Metafor sur les stations étaient illisibles a cause des nombreux warning (plus de 10k lignes chez moi…). J'ai fait quelques modifs pour les limiter :

• Remplacement des conversions implicites par des conversion explicites
• Ajout de retour chariots en fin de fichier
• La class EventAware était ignorée dans Element, ElementSet, Interaction, et InteractionSet. J'ai ajouté la déclaration “class EventAware” en entête de ces fichiers pour forcer la reconnaissance (j'ai un doute savoir si j'ai bien fait)

A présent les messages d'erreur sont limités aux lignes suivantes que je n'ai pas réussi à résoudre ou que je n'ai pas osé toucher (entre autre SloanConnexionBuilder n'est pas évident avec des Goto dans tous les sens, j'ai eu peur d'introduire une crasse) :

/home/biotteau/b/oo_meta/geniso/_src/geniso.i:57: Warning(509): Overloaded method gIso::Vec3::__mul__(gIso::Vec3 &) is shadowed by gIso::Vec3::operator *(gIso::Vec3 const &) const at /home/biotteau/b/oo_meta/geniso/src/gisoVec3.h:48.
/home/biotteau/b/oo_meta/mtFEMBase/DofSet.h:41: Warning(503): Can't wrap 'operator Dof*' unless renamed to a valid identifier.
/home/biotteau/b/oo_meta/mtFEMBase/DofSet.h:48: Warning(503): Can't wrap 'operator ==' unless renamed to a valid identifier.
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:1442: attention : converting to ‘int’ from ‘float’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:1451: attention : converting to ‘int’ from ‘float’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:1454: attention : converting to ‘int’ from ‘float’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:656: attention : converting to ‘int’ from ‘double’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:666: attention : converting to ‘int’ from ‘double’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:680: attention : converting to ‘int’ from ‘double’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:692: attention : converting to ‘int’ from ‘double’
/home/biotteau/b/oo_meta/mtFEMBase/SloanConnexionBuilder.cpp:725: attention : converting to ‘int’ from ‘double’
/home/biotteau/b/oo_meta/mtFEMBase/StrID.h:20: Warning(401): Maybe you forgot to instantiate 'IDBase< DummyStrID >' using %template.
/home/biotteau/b/oo_meta/mtFEMBase/StrID.h:20: Warning(401): Nothing known about base class 'IDBase< DummyStrID >'. Ignored.
/home/biotteau/b/oo_meta/mtMath/SmallList.inl:62: attention : cannot receive objects of non-POD type ‘class Field’ through ‘...’; call will abort at runtime
/home/biotteau/b/oo_meta/mtPython/PythonMultiParameterFunction.cpp:102: attention : cannot pass objects of non-POD type ‘class Field’ through ‘...’; call will abort at runtime

Les résultats de compilations sont désormais plus propre, avec plus de 8000 lignes de warning en moins…