Introduction à FreeFem++, une plateforme Eléments Finis libre

FreeFem++ : un outil libre pour la
résolution d’EDP par la méthode des
Eléments Finis
Arnaud Lejeune
Pôle Calcul Scientifique – Département Mécanique Appliquée
20 Février 2014
Logiciels Libres et/ou OpenSource
https://www.projet-plume.org/relier
Logiciels Libres
http://www.gnu.org/philosophy/fre
e-sw.fr.html
« Logiciel libre » [free software]
désigne des logiciels qui respectent
la liberté des utilisateurs. En gros,
cela veut dire que les utilisateurs ont
la liberté d'exécuter, copier,
distribuer, étudier, modifier et
améliorer ces logiciels.
OpenSource
http://opensource.org/
Un logiciel OpenSource est un
logiciel qui peut être librement
utilisé, modifié et partagé (avant ou
après modifications) par quiconque
Solveurs Eléments finis
Outils de CAO
•
FreeCAD
http://www.freecadweb.org/
•
•
http://www.freefem.org/ff++/
LibreCAD (2D)
FreeFEM++,
http://librecad.org/cms/home.html
•
Mailleurs
•
GMSH
http://geuz.org/gmsh/
•
MESHLAB (3D)
http://meshlab.sourceforge.net/
POST-TRAITEMENT
•
http://fenicsproject.org/
•
GETFEM++,
http://download.gna.org/getfem/html/ho
mepage/
•
ELMER,
http://www.csc.fi/english/pages/elmer
Paraview
http://geuz.org/gmsh/
•
FENICS,
VISIT
•
FEEL++,
http://www.feelpp.org/
https://wci.llnl.gov/codes/visit/home.html
…..
Les licences Gnu GPL, BSD, …. sont des licences OpenSource créées initialement aux USA
Les licences CeCILL sont une adaptation au droit français (http://www.cecill.info/)
FREEFEM++: Tour d’horizon
Caractéristiques
• Environnement intégré FreeFem++-cs
•Langage de script basé sur la syntaxe du C/C++
•Librairie Eléments Finis avec éléments standards et non
standards
• Résolution de systèmes linéaires performante
( SuperLU, MUMPS,UMFPACK,….)
• Raffinement adaptatif et partition de maillage
• Calcul de valeurs propres
• Interface MPI dispo
• Post-traitement inclus en version allégée
• Disponible sous Linux/Mac OSX et Windows
Processus Elément Fini
Un élément fini est défini par
(Raviart) :
• Un domaine (connexe,
compact , d’intérieur non vide)
• Un ensemble fini de points
• Un ensemble de fonctions
d’interpolation
Un processus EF nécessite la connaissance de :
• L’espace de définition de la solution recherchée (test function) u et de la fonction test v (trial function)
• Géométrie de l’élément de référence et le type d’interpolation (ordre, nœuds, continuité)
• La méthode d’intégration numérique
Principe de base d’un solveur éléments finis
Création du modèle
Nœud 3
Trouver u  V tel que
KU  F
 u.d d  
 k11
x



x
k n1
x
x
x
k1n 

x

   x


 x
x

k nn 
x
x
Nœud 1
f .v d
Nœud 2
Choix de l’élément
Pt Intégration 3
x
x
x
x
x 
x 
Pt Intégration 2
Pt Intégration 1
Choix de l’intégration
L’utilisateur se concentre sur :
• la formulation du problème
• le choix des éléments (espace fonctionnel, ordre, type d’interpolation)
•Le choix des méthodes d’approximation et de résolution
FREEFEM++: Tour d’horizon
Sorties
graphiques/co
nsole/dans un
fichier
Fichier .edp
Création/ Import
Géométrie/Maillage
FREEFEM++
Résolution
Traitement
Numérique (lié à
la discrétisation)
A la différence des logiciels commerciaux,
• Pas de mode « Clic-bouton »
• l’utilisateur doit maîtriser son modèle
• l’utilisateur maîtrise / a conscience des techniques de résolution utilisées
Analyse
Formulation Faible
FREEFEM++: Tour d’horizon
Construction/ Définition d’un maillage
Forme prédéfinie:
mesh nomMaillage = square ( n ,m , [x0+(x1-x0)*x, y0+(y1-y0)*y])
: génération d’un
maillage régulier sur un carré [x0,x1]x[y0,y1]
A partir de frontières définies par des courbes paramétrées:
border nomBord1 (t=debut, fin){x=f(t), y=g(t), label=numLabel }
mesh nomMaillage = buildmesh(nomBord1(valeur1)+nomBord2(valeur2)+…)
1: border a(t=0,2*pi){ x=cos(t); y=sin(t);label=1;}
2: border b(t=0,2*pi){ x=0.3+0.3*cos(t); y=0.3*sin(t);label=2;}
3: plot(a(50)+b(+30)) ; // to see a plot of the border mesh
4: mesh Thwithouthole= buildmesh(a(50)+b(+30));
5: mesh Thwithhole = buildmesh(a(50)+b(-30));
6: plot(Thwithouthole,wait=1,ps="Thwithouthole.eps");
7: plot(Thwithhole,wait=1,ps="Thwithhole.eps");
Note : Exemple repris du manuel utilisateur FreeFEM++
FREEFEM++: Tour d’horizon
Définition de l’espace d’interpolation
fespace nomEspace (nomMaillage, type Elément Finis)
Familles
FreeFEM ne connaît pas :
•Les éléments coque/plaque
Crouzeix-Raviart
Lagrange Discontinu
Lagrange
Nédélec
Définition du problème variationnel
Raviart-Thomas
problem nomProblem (u,v, [solver])
a(u,v)-l(v) + (cdl)
Note : On peut également utiliser la commande varf pour définir séparément les formes bilinéaires et
linéaires
FREEFEM++: Exemple Simpliste
Exemple Basique :
u  f sur le domaine Ω
u  u d sur le bord Γ d
u

 g sur le bord Γ n
n
Formulation faible:
Trouver u V tel que
 u.v d  


f .v d   N g.v d v V0

Elément Fini utilisé:
Lagrange - 3 nœuds – ordre 1
// Creation du maillage
mesh nomMail= square(10,20,[1+2*x,-1+3*y]);
// Creation de l’espace d’interpolation
fespace Vh(nomMail, P1);
// Definition de l’inconnue et de la variable virtuelle
Vh u,v;
// Definition de la formulation variationnelle
problem problem1(u,v) = int2d(nomMail)(dx(u)*dx(v)+dy(u)*dy(v))
– int2d(nomMail) (1*v) + on(2, u=5);
// Résolution
problem1;
// Sortie graphique
plot(u,fill=1,wait=1,value=true,ps="simpliste.eps");
ud  5
g x   0
f  x, y   1
Elasticité linéaire
Trouver u  V tel que v  V0
div  f sur Ω
u  u sur Γ
d
   u . v  d   f .v d
  g.v d

d

 n0
  .n  g sur le bord Γ n0
  .n  0 sur le bord Γ n1
ux  0
d
20
g
g x   g ( 20  y )
// Création du maillage
mesh barrage= buildmesh(horizontal(10)+diagonal(15)+vertical(-20));
plot(barrage,wait=1,ps="barrage.eps");
gy  0
f  x, y   0
E  40e9
  0.2
// Création du domaine
int Dirichlet=1,Neumann=2;
border horizontal(t=0,1){ x=10*t; y=0;label=Dirichlet;}
border vertical(t=0,1){ x=0; =20*t; label=Neumann;}
border diagonal(t=0,1){ x=10*t; y=-20*t+20;}
10
// Création de l’espace variationnel
fespace Vh(barrage,P1);
Vh ux,uy,vx,vy;
// Définition des paramètres
real rho=1000.0;
real gravite=9.81;
real longueur=10.0;
func gx= rho*gravite*(20-y);
real f=0.;
real E = 40e9, nu = 0.2, mu= E/(2*(1+nu));
real lambda = E*nu/((1+nu)*(1-2*nu));
Elasticité linéaire
Trouver u  V tel que v  V0


// Définition de la fonction de calcul des déformations
macro epsilon(ux,uy) [dx(ux),dy(uy),(dy(ux)+dx(uy))]
// Définition de la fonction divergence
macro div(ux,uy) ( dx(ux)+dy(uy) )
// Definition et resolution de la formulation variationnelle
solve problemeBarrage([ux,uy],[vx,vy])=
int2d(barrage)(lambda*div(ux,uy)*div(vx,vy) +2.*mu*( epsilon(ux,uy)‘ *
epsilon(vx,vy) ) )
- int2d(barrage)(f*vy) - int1d(barrage,Neumann)(gx*vx)
+ on(Dirichlet,ux=0,uy=0);
// sauvegarde des deplacement en x et y
plot(ux,wait=1,value=1,fill=1,ps="barrageUx.eps");
plot(uy,wait=1,value=1,fill=1,ps="barrageUy.eps");
// sauvegarde du maillage déformé
real coef =10000 ;
mesh barrage1 = movemesh(barrage, [x+ux*coef, y+uy*coef]);
plot(barrage1,wait=1,value=1, fill=1, ps=« deform.eps");
real dxmax = ux[].max;
real dymax = uy[].max;
// Sortie console des deplacements max
cout << " - ux max = "<< dxmax<< " uy=" << dymax << endl;
  u . v  d 


f .v d  
 n0
g .v d
Conclusions et Regrets
Bilan
• Résolution d’edp par MEF:
• linéaire, non linéaire,
• statique, quasi-statique,
d’évolution,
• multi-domaines,
• maillage adaptatif,
•…
• Solveurs actuels et parallélisation
disponibles
Point fort
• descendance de Pironneau 1987
• http://www.um.es/freefem/ff++/pmwiki.php
•
• Large Communauté
d’utilisateurs/développeurs
• Evolution pernanente : freefem++-cs
Applications :
• Fluides incompressibles
• Navier-Stokes incompressible
• Elasticité
• Interaction fluide-structure
• …..
REGRETS
• Absence d’éléments quadrangle, hexaèdre
• Manque éléments coque et plaque (ajout
non trivial)
• Post-traitement parfois délicat

Download Report