Tutoriel Code_Aster

Une fois n’est pas coutume, il pleut en Haute‐Garonne ! Alors j’ai tué le temps jetant un œil à Code_Aster, qui est un code de calcul de structure thermomécanique par la méthode des éléments finis isoparamétriques. Il est développé par EDF sous licence GNU GPL v3. Il permet de faire à peu près tout ce qui est imaginable en mécanique, voir à ce propos la plaquette de présentation (PDF).

Ce code de calcul est intégré à la suite de logiciels libres Salomé‐Méca, qui contient un préprocesseur, Code_Aster, et un post‐processeur/visionneur pour voir les résultats. Aujourd’hui, nous allons utiliser le code en version autonome (stand alone) et nous utiliserons notre éditeur de texte préféré, gmsh, astk, puis de nouveau gmsh pour voir les résultats de nos calculs.

Sommaire

• Installation
• Calcul de poutre
  • Création de la géométrie avec Gmsh
  • Fichier de commandes
  • Astk
  • Résultat
• Aller plus loin

Installation

Cela se passe ici. Deux options :

1. soit on compile code_aster ;
2. soit on n’aime pas compiler et l’on télécharge le binaire de Salomé‐Méca qui contient Code_Aster de façon préinstallé et quelques Gio d’outils dont nous ne nous servirons pas.

La compilation se passe assez bien, et les paquets prérequis (voir les instructions de compilation sur leur site) se trouvent assez facilement.

Calcul de poutre

Nous allons étudier le comportement mécanique d’une poutre encastrée d’un côté et soumise à un effort ponctuel de l’autre côté :

Nous allons le faire de trois façons :

1. En modélisation poutre 1D ;
2. En modélisation plaque 2D ;
3. En modélisation complète 3D.

Création de la géométrie avec Gmsh

Pour fonctionner, Code_Aster a besoin d’un fichier de commandes, et d’un fichier décrivant un maillage : une liste de nœuds et une liste d’éléments reliant ces nœuds. On pourrait lui fournir un fichier texte contenant les coordonnées géométriques de chaque nœud du maillage, mais vu qu’on a la flemme et que cela peut être assez ennuyeux pour des problèmes complexes, on va demander à Gmsh de le faire pour nous.

On crée tout d’abord la géométrie de notre problème à l’aide de points, de lignes, surfaces et volumes, on doit aussi définir des groupes d’intérêts (la poutre entière, la partie encastrée et la partie sur laquelle on applique la force). On peut jouer cinq minutes avec la partie interface graphique de Gmsh, pour lequel on trouvera de nombreux tutoriaux sur le Web, mais on en revient vite à un fichier texte.

Voici donc poutre1d.geo :

//== paramètres ==
//taille maillage autour des nœuds. 2,5 mm entre deux nœuds.
cl__1 = 2.5;

// == géométrie ==
//points
Point(1) = {0,0, 0, cl__1}; // extrémité encastre de ma poutre
Point(2) = {100,0, 0, cl__1}; // extrémité libre de ma poutre, soumise à une force
//lignes
Line(1) = {1, 2}; // on crée la ligne entre Point1 et Point2

//== groupe ==
Physical Point("encastr") = {1};
// on encastrera le nœud correspondant à ce point

Physical Point("force") = {2};
//on appliquera la force sur le nœud correspondant à ce point

Physical Line ("poutre") = {1};
// notre poutre sera constituée de tous les nœuds et éléments correspondant à cette ligne

Une fois ce fichier poutre1d.geo créé, on l’ouvre avec gmsh (terminal : gmsh poutre1d.geo). On clique sur Mesh → 1D, le maillage est fait, mais on ne le voit pas, car seule la géométrie est affichée ! Donc, Tools → Options, onglet Visibily de Geometry, on décoche Points et Lines, et dans Mesh, on coche Nodes et Lines. Cela donne ceci :

Notez qu’avec Tools → Options, dans l’onglet list bowser, on peut visualiser ou masquer (taper la touche Entrée au clavier après avoir cliqué sur le nom du groupe dans la fenêtre) les groupes que nous avons créés et leurs affectations. C’est très pratique. On voit par exemple bien que notre groupe « poutre » est constitué de tous les éléments de la poutre.

Pour sauvegarder notre maillage, on fait File → Export et l’on choisit le format de maillage appelé .MED, on obtient donc un beau mesh1d.med. Surtout, on veille à ce que tout soit décoché dans la fenêtre surgissante qui apparaît et l’on clique rapidement sur OK.

De même, voici poutre2d.geo, qu'on maille en 2D avec gmsh:

//== paramètres: ==
//taille maillage autour des nœuds. 2,5 mm entre deux nœuds.
cl__1 = 2.5;
L=100; //longueur poutre de 100 mm
R=5; // ratio longueur/largeur
l=L/R;

//== géométrie ==
//points
Point(1) = {0, 0, 0, cl__1};
Point(2) = {L, 0, 0, cl__1};
Point(3) = {L, l, 0, cl__1};
Point(4) = {0, l, 0, cl__1};
Point(5) = {L, l/2, 0, cl__1};
//lignes
Line(1) = {1, 2};
Line(2) = {2, 5};
Line(3) = {5, 3};
Line(4) = {3, 4};
Line(5) = {4, 1};

//surface
Line Loop(1) = {1, 2, 3, 4, 5}; //on crée une boucle de lignes
Plane Surface(1) = {1}; // on crée une surface à partir de la boucle

//== groupe ==
Physical Line("encastr") = {5}; // on encastrera cette ligne
Physical Point("force") = {5}; // lieu application force
Physical Surface("poutre") = {1}; // notre poutre 2D

Et poutre3d.geo qu’on maille en 3D avec gmsh :

//== paramètres: ==
//taille maillage autour des nœuds.
cl__1 = 5;
L=100; // longueur poutre
R=5; // ratio longueur/largeur
l=L/5;

//== géométrie ==
//points
Point(1) = {0,0, 0, cl__1};
Point(2) = {L,0, 0, cl__1};
Point(3) = {L,l, 0, cl__1};
Point(4) = {0,l, 0, cl__1};
Point(5) = {L, l/2, 0, cl__1};
//lignes
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
//surface
Line Loop(1) = {1, 2, 3, 4};
Plane Surface(1) = {1};
Point{5} In Surface{1}; // pour que le point 5 soit contenu dans la surface

//volume
Extrude {0,0,-3}{Surface{1};Layers{3}; Recombine;}
//on extrude la surface 1 de -3 mm selon l’axe Z
//en créant 3 éléments dans l’épaisseur avec l’aide de calques


//== groupe ==
//on sait que c’est la surface 25 parce qu’on le visualise sous gmsh en affichant « surface label ».
//il peut y avoir une erreur lors de l’importation si le numéro de la surface créée par l’extrusion n’est pas 25.
//      C’est pas grave, on regarde à quoi correspond la surface à encastrer, on trouve son label, et mon modifie les lignes ci-dessous.
Physical Surface("encastr") = {25}; // on encastrera cette surface
Physical Point("force") = {5}; // lieu application force
Physical Volume("poutre") = {1}; // notre poutre 3D

//== maillage ==
Transfinite Line{1,3}=8*R+1; // 8*R élém dans la longueur = 41 nœuds entre lignes 1 et 3
Transfinite Line{4,2}=8+1; // 8 élém dans la largeur =  9 nœuds entre lignes 4 et 2
Transfinite Surface "*"; // on veut un maillage propre
Recombine Surface "*"; // on veut un maillage quadra

Nous voici maintenant avec trois maillages au format .med. Il nous faut maintenant créer notre fichier de commandes !

Fichier de commandes

#U1.03.02
DEBUT();

#on charge le fichier de maillage .MED, unité logique 20
mesh=LIRE_MAILLAGE(
    INFO=2,
    INFO_MED=2,
    UNITE=20,
    FORMAT='MED',
);

#on a importé le maillage et ses groupes, on crée d'autres groupes:

mesh=DEFI_GROUP(
    reuse =mesh,
    INFO=2,
    MAILLAGE=mesh,
    #on crée un groupe nommé TOUT qui contient toutes les mailles du maillage.
    #on ne va pas s'en servir, mais ça peut être utile
    CREA_GROUP_MA=_F(NOM='TOUT',TOUT='OUI',),
    #on crée un groupe de nœuds qui contient tous les nœuds de toutes les mailles.
    # Il faut le faire quand le maillage provient de Gmsh, car Gmsh transforme les nœuds en maille, on les retransforme ici en nœuds
    CREA_GROUP_NO=_F(TOUT_GROUP_MA='OUI',),
);

#on affecte au groupe de mailles 'poutre' créé avec gmsh,
#   des éléments finis de types Poutre, ici POU_D_T
model=AFFE_MODELE(
    MAILLAGE=mesh,
    AFFE=(
        _F(
          GROUP_MA=('poutre',),
          PHENOMENE='MECANIQUE',
          MODELISATION='POU_D_T',
        ),
    ),
);

#on définit un matériaux, ici de l''acier:
#  Module d'Young' E = 210000 N/mm2
#  Coefficient de Poisson, nu = 0.3
#  masse volumique = 8e-9 tonne/mm3
steel=DEFI_MATERIAU(ELAS=_F(E=210000.,NU=0.3,RHO=8e-9),);

#U4.43.03
#on assigne notre matériaux à nos mailles du groupe 'poutre'
material=AFFE_MATERIAU(
    MAILLAGE=mesh,
    AFFE=_F(GROUP_MA=('poutre',), MATER=steel,),
);


#U4.42.01
#On assigne à nos éléments poutre POU_D_T une section rectangulaire de largeur 20 mm et d’épaisseur 3 mm

elemcar=AFFE_CARA_ELEM(
    MODELE=model,
    INFO=2,
    POUTRE=(
        _F(
            GROUP_MA=('poutre',),
            SECTION='RECTANGLE',
            CARA=('HY','HZ',),
            VALE=(3,20),
        ),
    ),
);
#on interdit toute rotation et translation aux nœuds du groupe 'encastr' (1 seul nœud ici).
#   cela simule l'encastrement
encast=AFFE_CHAR_MECA(
    MODELE=model,
    DDL_IMPO=(
        _F(
            GROUP_NO=('encastr',),
            DX=0,DY=0,DZ=0,DRX=0,DRY=0,DRZ=0,
        ),
    ),
    INFO=1,
);

# on applique 500 N selon la direction -Z au nœud de notre groupe 'force'
force_f=AFFE_CHAR_MECA(
    MODELE=model,
    FORCE_NODALE=_F(
        GROUP_NO=('force',),
        FZ=-500,
    ),
    INFO=2,
);

#U4.51.01
#on compile les précédents concepts pour le calcul
stat=MECA_STATIQUE(
    MODELE=model,
    CHAM_MATER=material,
    CARA_ELEM=elemcar,
    EXCIT=(
        _F(CHARGE=encast,),
        _F(CHARGE=force_f,),

    ),
);

# Par défaut, sont calculés uniquement les déplacements et les réactions nodales aux points de Gauss des éléments, je crois.
# du coup on enrichit le concept "stat" pour lui demander d'autres choses.
# SIEF_ELNO: ici, efforts théorie des poutres au niveau des nœuds des éléments
# SIPO_ELNO: ici, contraintes dans la poutre, au niveau des nœuds des éléments
# SIPM_ELNO: ici, contrainte max dans la poutre
# REAC_NODA: forces/moments aux nœuds limites
stat=CALC_CHAMP(
    reuse =stat,
    RESULTAT=stat,  
    CONTRAINTE=(
        'SIEF_ELNO','SIPO_ELNO','SIPM_ELNO',
    ),
    FORCE=('REAC_NODA',),
);

#on imprime ça dans un fichier de sortie .med, unité logique 80.
#on n'imprime que les déplacements et les contraintes
# (on n'affiche pas tout ce qu'on a calculé, genre SIPM_ELNO ou REAC_NODA pourquoi pas !)
IMPR_RESU(
    FORMAT='MED', 
    UNITE=80,
    RESU=_F(
        RESULTAT=stat,
        NOM_CHAM=(
            'DEPL',
            'SIPO_ELNO',
            'SIPM_ELNO',
        ),
    ),
);

FIN();

Notez que les #U4.51.01 ou autres renvoient à la documentation.
On enregistre ce fichier texte en 1d.comm, par exemple, et nous allons lancer le calcul à l’aider d’astk.

Astk

Astk est l’outil permettant de mener à bien un calcul, on le lance via /opt/code_aster/bin/astk (si vous avez installé Code_Aster dans /opt).

On cherche à obtenir une fenêtre qui a cette allure :

Ensuite :

• File → New ;
• on choisit notre dossier de travail (path) ;
• dans la colonne d’icônes au milieu à droite, on clique sur l’icône en forme de dossier bleu, pour aller chercher son mesh1d.med et son 1d.comm ;
• on clique sur l’icône du dessus pour ajouter deux lignes, puis dans type pour la ligne, on choisit mess et rmed, dans name on les appelle ./log1d.mess et ./resu1d.rmed ;

• File → Save_As → 1d.astk.

• File -> New

• on choisit notre path / dossier de travail

• dans la colonne d'icônes au milieu à droite, on clique sur l’icône en forme de dossier bleu, pour aller chercher son mesh1d.med et son 1d.comm

• on clique sur l’icône du dessus pour ajouter deux lignes, puis dans type pour la ligne, on choisit mess et rmed, dans name on les appels ./log1d.mess et ./resu1d.rmed

• File -> Save_As -> 1d.astk

La colonne LU correspond à Logical Unit ou Unité Logique, c’est l’endroit de la mémoire où je ne sais quoi où l’on s’attend à trouver le fichier ; dans fichier.comm, on a précisé que l’unité logique était 20 pour le maillage .med et 80 pour le résultat .med. Les colonnes DRC veulent dire Données, Résultats, Compressé.

Une fois que cela est fait on clique sur Run ! Le calcul est lancé. Il se termine, on va voir le log1d.mess qui a été créé, il contient toutes les infos relatives au calcul. L’information la plus importante étant la dernière ligne.

Chez moi j’ai : EXECUTION_CODE_ASTER_EXIT_13223-localhost=0. Si le code renvoie 0, c’est que cela a fonctionné ! S’il renvoie 1, c’est que ça a planté et qu’il faut déboguer…

Résultat

Normalement tout a fonctionné, nous avons un beau resu1d.rmed que nous ouvrons avec gmsh (terminal gmsh resu1d.rmed).

On peut donc voir les déplacements et la contrainte, tout ce dont a besoin un mécanicien pour dimensionner des systèmes mécaniques !

Voici les paramètres sur lesquels agir pour afficher le déplacement multiplié par 10. Il faut afficher des Vectors et non pas l’Original Field. Comme ci‐dessous :

Pour les contraintes, SIPO_ELNO contient la contribution de chaque force ou moment aux contraintes de la poutre.

C’est grosso modo un vecteur de six composantes que voici :

Pour les afficher une par une, on se place dans Options → Visibility et, en bas, la première case à droite de la liste déroulante Original Field/Scalar Force/Vector/Tensor. Zéro correspond à SN et cinq correspond à SNT, par rapport au tableau ci‐dessus (je ne sais pas trop ce que présente SIPO_ELNO par défaut).

SIPM_ELNO, quant à lui, représente par défaut la contrainte maximum selon XX.

Voici d’autres visualisations avec les modèles 2D et 3D :

Aller plus loin

Code_Aster est très vaste, il contient près de 400 types d’éléments finis ! Pour aller plus loin, n’hésitez pas à lire la doc, qui contient aussi des exemples de calculs qui sont livrés avec le code.

Je vous conseille aussi notamment l’excellent livre sous licence libre de Jean‐Pierre Aubry, qui est un passage obligatoire pour prendre en main le code ! Le code date, en revanche, de la version 11 de Code_Aster ; mais une nouvelle version est en cours d’écriture !
On y fait notamment des analyses non linéaires avec du contact entre pièces et du frottement.
Aster Study vient aussi de faire son apparition.

Voilou, cher lecteur. N’hésite pas à t’amuser !

Je poste en commentaire les fichiers .comm de calcul en 2D et 3D.

Aller plus loin

• Journal à l’origine de la dépêche (210 clics)
• Code Aster, sur Wikipédia (200 clics)
• Code Aster, site officiel (279 clics)