Journal Tutorial Code_Aster

Sommaire

• Code_Aster, qu'est-ce c'est ?
  • Installation
• Calcul de Poutre
  • Création de la géométrie avec Gmsh
  • Fichier de commande
  • Astk
  • Résultat
• Aller plus loin

Bonjour Nal', une fois n'est pas coutume, il pleut en Haute-Garonne ! Alors j'ai tué le temps jetant un œil à Code Aster. Je laisse ici mes notes.

Code_Aster, qu'est-ce c'est ?

Code_Aster est un code de calcul de structure thermo-mécanique par la méthode des éléments finis isoparamétriques.
Il est développé par EDF sous licence GNU/GPLv3.
Il permet de faire à peu près tout ce qui est imaginable en mécanique, voir à ce propos-là plaquette de présentation

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/vieweur pour voir les résultats.

Aujourd'hui, nous allons utiliser le code en version stand_alone, et nous utiliserons notre éditeur de texte préferé, gmsh, astk, puis de nouveau gmsh pour voir les résultats de nos calculs.

Installation

Cela se passe ici.
Deux options:
1. Soit on compile code_aster
1. Soit on n'aime pas compiler et on télécharge le binaire de Salome-Méca qui contient code_aster de façon préinstallé, et quelques Gio d'outils dont 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 3 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 commande, 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'aides 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 5 minutes avec la partie GUI de Gmsh pour lequel on trouvera de nombreux tutoriaux sur le web, mais on en revient vite à un fichier texte.

Voici donc poutre1d.geo:

//== parametres ==
//taille maillage autour des noeuds. 2.5mm entre deux noeuds.
cl__1 = 2.5;

// == geometrie ==
//points
Point(1) = {0,0, 0, cl__1}; // extremite encastre de ma poutre
Point(2) = {100,0, 0, cl__1}; // extremite libre de ma poutre, soumise à une force
//lignes
Line(1) = {1, 2}; // on cree la ligne entre Point1 et Point2

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

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

Physical Line ("poutre") = {1};
// notre poutre sera constitue des tous les noeuds et elements correspondant à cette ligne

Une fois ce fichier poutre1d.geo crée, 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 seul la géometrie est affichée ! Donc Tools > Options, Onglet Visibily de Géometrie, on décoche Points et Lines et dans Mesh, on coche Nodes and Lines. Cela donne ceci:

Notez qu'avec Tools > Options, dans l'onglet list bowser, on peut visualiser/désafficher (toucher enter du clavier une fois cliqué sur le nom du groupe dans la fenêtre) les groupes que nous avons crée 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 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 pop-up qui apparaît et on clique rapidement sur OK.

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

//== parametres: ==
//taille maillage autour des noeuds. 2.5mm entre deux noeuds.
cl__1 = 2.5;
L=100; //longeur poutre de 100mm
R=5; // ratio longueur/largeur
l=L/R;

//== geometrie ==
//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 cree 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 mesh en 3D avec gmsh:

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

//== geometrie ==
//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 extdure la surface 1 de -3mm selon l''axe Z
//en créant 3 éléments dans l''épaiseur avec l''aide de layers


//== 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''import 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

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

Nous voici maintenant avec 3 maillages au format.med Il nous faut maitenant créer notre fichier de commande !

Fichier de commande

#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 grée un groupe de noeud qui contient tous les noeuds de toutes les mailles.
    # Il faut le faire quand le maillage provient de Gmsh, car Gmsh transforme les noeuds en maille, on les retransforme ici en noeud
    CREA_GROUP_NO=_F(TOUT_GROUP_MA='OUI',),
);

#on affecte au groupe de maille 'poutre' crée 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 20mm et d'épaisseur 3mm

elemcar=AFFE_CARA_ELEM(
    MODELE=model,
    INFO=2,
    POUTRE=(
        _F(
            GROUP_MA=('poutre',),
            SECTION='RECTANGLE',
            CARA=('HY','HZ',),
            VALE=(3,20),
        ),
    ),
);
#on interdit toutes rotations et translations aux noeuds du groupe 'encastr' (1 seul noeud 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 500N selon la direction -Z au noeud 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 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 à Logique Unité, c'est l'endroit de la mémoire où je ne sais quoi où on s'attend à trouver le fichier, dans le 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 Datas, Récrire, 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 était crée, 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 renvoit 0, c'est que cela a fonctionné ! S'il renvoit 1, c'est que ça a planté et qu'il faut débuger…

Résultat

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

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

Voici les paramètres sur lequels 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/moment aux contraintes de la poutre.

C'est grossomodo un vecteur de 6 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. 0 Correspond à SN et 5 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 présente par défaut la contrainte maximum selon XX.

Voici d'autre visu 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 calcul 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 par contre 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 journal, n'hésite pas à t'amuser !

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