Exemple d'utilisation du code probThermoHygroMecanique

Le programme probThermoHygroMecanique est la dernière version du code. On retrouve une copie du code ici mais il se retrouve aussi dans "app/src/Test.Formulation/probThermoHygroMecanique.cc" de l'arborescence MEF++. Il diffère des versions précédentes par le fait que l'on peut maintenant résoudre les problèmes bidimensionnnelle ou tridimensionnelle. Le point de départ de la modification du code original est que je voulais pouvoir utiliser les mêmes fichiers champs, à quelques modifications près, dans le cas ou on voudrait faire les deux résolutions. Puisque je m'étais fixé cette contrainte il semblait logique de ne faire qu'un seul code qui lirait, suivant la dimension du problème, les bonnes données.

Pour illustrer l'utilisation du code en 2d et 3d j'ai choisi un problème hygro-mécanique (en M et U). Il s'agit d'un problème tel que présenté par Ganev: une plaque "carré" de MDF imperméable sur la base et les rives et avec un échange sur la surface du dessus. Sous l'effet de l'humidité ambiante (qui est supérieure à l'humidité initiale du panneau) le panneau se déforme. Évidemment le 2D dans ce cas est discutable mais je veux seulement montrer comment on résout les deux types de problèmes.

Matériaux

On suppose un panneau de MDF "symétrique" formé d'un noyau de densité 650 kg/m³ entre deux couches de surfaces de densité 800 kg/m³.

Modèle

Le modèle usuel, diffusion pour M avec la condition d'échange sur la face supérieure comme "moteur". Pour la mécanique c'est la dilatation dû à M qui engendre la déformation, il n'y a pas d'autre force. Le panneau est "simplement posé" alors on a seulement une condition de contact pour la surface du dessous.

C.L.

Pour M on a une seule condition limite : échange sur le dessus.

Pour U : des symétries car on ne prend que le quart de la plaque et une condition de contact sur la base.

Marche à suivre

Voyons comment on résout en 2d et en 3d

étape 1.

Création du maillage. On met des étiquettes sur différentes parties, "entités", de la plaques qui nous serviront par la suite pour les C.L. et pour la description des matériaux. Pour simplifier les choses j'ai conservé les noms des entités que l'on soit en 2d ou 3d. Les entités mises en couleur sont celles communes aux deux maillages.

On suppose que le panneau se compose de 3 "volumes" (un "biscuit") chacun de ces volumes ayant des propriétés physiques différentes. On les nommera volume_dessus, volume_milieu et volume_dessous.

On identifiera la surface du dessus : surface_dessus pour M, pour U : surface_dessous pour le contact et les surfaces de symétrie SYM_X, SYM_Y

Différences entre le 2d et le 3d :

Naturellement les volumes sont des surfaces en 2d, les surfaces des courbes et les courbes des points.
Une des conditions de symétrie devient inutile.

On choisi de prendre une coupe transversale dans le plan x0z alors la condition de symétrie par rapport à l'axe des x devient inutile.

étape 2.

Création du fichier champs. Voir les commentaires dans le fichier champs donnees_exemple/article4.champs. Le fichier est coupé en trois parties :

Caractérisation du type de problème et choix du solveur. On spécifie la dimension de l'espace, les équations à résoudre (si on résout en thermique et/ou mécanique), si des mouvements de corps rigide sont prévus, hystérèse par rapport à M. Le type de solveur instationnaire pour M (et T possiblement). Le solveur linéaire pour U. Imposition de valeur de temps pour lesquels on doit résoudre le problème instationnaire (pour les discontinuités dans le temps des données), points particuliers dans le domaine pour lesquels on veut sauvegarder une évaluation dans un fichier, etc.
Choix de l'interpolation et des fonctions utilisé pour définir les conditions initiales et limites. Définition de la valeur initiale de M et U, la définition de l'humidité ambiante, etc.
Définition des différents matériaux composant le panneau et assignation des ces matériaux aux entités formant le panneau. D'abord la définition des propriétés des matériaux, puis on assigne celle-ci aux différentes entités avec un "champssurgeom"

C'est surtout dans la partie 3 que l'on distinguera le cas bi- et tri- dimensionnelle. Notons que pour "couvrir" tous les cas on ne commente pas les déclarations de champs inutiles (les champs 2d quand on fait du 3d et vice versa). Cette généralité se fait au dépend de l'efficacité puisque l'on à ainsi une série de pointeur qui sont initialisés sans être utilisés.

Différences entre le 2d et le 3d :

Peu importe la dimension de l'espace, dans le cas orthotrope ou isotrope transverse (comme c'est le cas ici) il faut donner les 3 modules de Young et les 3 coefficients de Poisson (dans le cas des déformations planes). La différence entre le cas 2d et 3d est donc minime pour M et pour U. Il faut faire la permutation appropriée pour tenir compte des changement d'axes (pour le coefficients nécessaires à M et U mais aussi pour éliminer les corps rigides). Pour U il faut choisir l'hypothèse plane voulue (contraintes planes ou déformations planes) et faire appel au champs "coefelasticite2D" pour le cas 2d plutôt que "to4elasticite" qui défini le tenseur des coefficients d'élasticité dans le cas 3d. Dans le fichier champs cela se traduira par un nombre restreint d'entrée spécialisée pour le 2d ou 3d.

Données distinctes pour le 3d :

scalaire dimension 3

v3dquad "U3D" [0, 0, 0]

empilementto2 "K_M3D" [zero,zero,D,zero,zero,zero]
empilementto2 "Beta3D" [beta1,beta1,beta3,zero,zero,zero]
to4elasticite "Eijkl3D" [E1,E1,E3,Nu12,Nu13,Nu13,G1,G3,G3]

alias K_M=K_M3D
alias BetaM=Beta3D
scalaire DirectionRigide 2

Données distinctes pour le 2d:

scalaire dimension 2

v2dquad "U2D"    [0, 0]

empilementto2 "Beta2D" [beta1,beta3,beta1,zero,zero,zero]       <------- ATTENTION AUX AXES (PERMUTATION Y-Z)
empilementto2 "K_M2D" [zero,D,zero,zero,zero,zero]                 <------- ATTENTION AUX AXES (PERMUTATION Y-Z)
coefelasticite2D "Eijkl2D" [DP,E1,E3,E1,Nu13,Nu13,Nu12,G3]     <------- ATTENTION AUX AXES (PERMUTATION Y-Z)

alias K_M=K_M2D
alias BetaM=Beta2D
scalaire DirectionRigide 1

étape 3.

Création des fichiers de conditions limites et du fichier champs pour le contact. Pour notre exemple ces fichiers sont simples :

Pour M :
neumannnonlin scalaire "surface_dessus" "M" fonction(M,h,M_infini)=[h*(M-M_infini)]

Pour U en 2d:
dirichlet scalaire "SYM_Y" "U2DX" "zero"

Pour U en 3d :
dirichlet scalaire "SYM_Y" "U3DX" "zero"
dirichlet scalaire "SYM_X" "U3DY" "zero"

Notre problème de mécanique est mal posé puisque nous permettons des mouvements de corps rigide. Néanmoins ces déplacements sont "connus" i.e. on connaît la direction des déplacements et on peut les éliminer lors de la résolution du problème de contact.

Pour le contact on doit définir :
1) Une valeur de défaut (GapMax) pour "l’écart" entre le panneau et le support. Valeur sans importance puisque l'on travaille uniquement sur "surface_dessous"
2) La fonction écart (GapPeau) défini sur toute la peau
3) Le champ décrivant la normale extérieure au support
4) Le type du champ du multiplicateur Lamdba utilisé pour résoudre le problème d'optimisation sous-jacent.

Pour le contact 2-d

scalaire GapMax 10
scalsurgeom GapPeau defaut=GapMax surface_dessous=U2DY
vectoriel2d normale_ext2D [0,1]
scallin LambdaCtc 0

Pour le contact 3-d

scalaire GapMax 10
scalsurgeom GapPeau defaut=GapMax surface_dessous=U3DZ
vectoriel2d normale_ext2D [0,0,1]
scallin LambdaCtc 0

étape 4.

Dans le cas présent on a une version 2d et une version 3d du maillages article42d et article43d (les données sont dans donnees_exemple). Puisque l'on a un seul fichier champ, donnees_exemple/article4.champs, permettant de résoudre les deux problème on passe de l'un a l'autre en utilisant des liens symbolique : ont fait un lien entre les fichiers de la dimension voulue et des fichiers de préfixe article4.

Par exemple pour le 2d :
ln -s article42d.enti article4.enti
ln -s article42d.geom article4.geom
ln -s article42d.mail article4.mail
ln -s article42d_hyd.CL article4_hyd.CL
ln -s article42d_mec.CL article4_mec.CL
ln -s article42d-peau.champs article4-peau.champs

On peut maintenant exécuter le programme : probThermoHygroMecanique article4 surface_dessous

On aura ainsi la solution du problème 2d. Le second paramètre de l'appel servant à préciser l'entité du maillage de peau sur laquelle on a du contact. Cela permet d'optimiser le problème en réduisant la taille du problème de contact à résoudre.

À la fin de la résolution on aura :

Deux (trois si on ajoute la température) séries de fichiers contenant M et U à chaque pas de temps dans le format MEF++. On pourrait visualiser chacun des fichiers via iMEF++ ou, en les convertissant, dans un autre logiciel (VU par exemple).
Deux fichiers (trois si on ajoute la température au problème), article4.M et article4.U, contenant les évaluations particulières demandées dans le fichier champ. Ces fichiers auront un format qui est facile a importer dans des logiciels tel que MAPLE, MATHEMATICA ou le très simple gnuplot. Le format est (voir article4.M par exemple):

Les solutions se retrouvent dans donnees_exemple/soln_2d et donnees_exemple/soln_3d

Quelques illustrations produites avec VU:

Deplacements verticaux temps final modele 2D

Déplacements verticaux au temps final modèle 3D

Et pour finir un traitement possible des évaluations ponctuelles :