La Méthode Des Éléments Finis, J. GARRIGUES

La methode des elements nis
Jean GARRIGUES [email protected] Janvier 2002
Ecole Superieure de Mecanique de Marseille
J. Garrigues
Avertissement
Ce texte a pour but de familiariser le lecteur avec la methode des elements nis. On n'y trouvera donc que le strict minimum de fondements theoriques. En particulier, on ne trouvera aucune demonstration de theoremes d'analyse fonctionnelle ni de theoreme de convergence. Le lecteur insatisfait est invite a se reporter aux cours d'analyse numerique. Par contre, on trouvera ici des informations pratiques pour maitriser l'outil et pour comprendre les options proposees par les logiciels et les problemes numeriques qui peuvent se presenter. Dans un souci de generalite, on ne fera allusion a aucun logiciel particulier. La terminologie employee s'e orcera donc d'^ etre generale. Il se peut qu'une certaine ((traduction)) soit necessaire pour adapter les concepts aux idiomes particuliers a chaque logiciel. Les quantites vectorielles et tensorielles sont notees en caracteres gras. Leurs composantes sont en caracteres ordinaires.
Convention typographique
J. Garrigues
Chapitre 1
Introduction
La plupart des problemes de physique peuvent se formuler ainsi : Trouver un champ (scalaire, vectoriel, ou tensoriel) u (M t) satisfaisant a un ensemble d'equations aux derivees partielles et d'equations ordinaires en tout point M d'un domaine et a tout instant t, et respectant des conditions aux limites (eventuellement fonction du temps) sur la frontiere @ du domaine. Les conditions aux limites sont des relations ou des valeurs imposees a u et/ou a ses derivees sur la frontiere. Si le temps t n'appara^t pas comme variable, on dit que le probleme est stationnaire, sinon c'est un probleme d'evolution.
Exemple 1: Le probleme thermique stationnaire

Le champ inconnu est un champ scalaire : la temperature T (M ). L'equation di erentielle est : T (M ) = f (M ) 8M 2 ou f (M ) est une fonction donnee. Les conditions aux limites sont par exemple: { Un champ de temperature stationnaire impose sur la partie de frontiere @ 1 : T (N ) = Tf (N ) 8N 2 @ 1 { un champ de ux thermique stationnaire impose sur la partie de frontiere @ 2 : grad T (N ) n(N ) = f (N ) 8N 2 @ 2 ou n(N ) est la normale exterieure en N a la frontiere. En principe, pour que la solution soit entierement determinee, il faut que @ 1 @ 2 = @ et @ 1 \ @ 2 = , c'est a dire qu'il faut imposer des conditions aux limites sur toute la frontiere. Cependant, la plupart des logiciels ont des conditions implicites appliquee aux parties de frontieres sur lesquelles on ne dit rien 1, et qu'il faut conna^tre en lisant les manuels.
1. Par exemple, dans un logiciel d'elements nis specialise dans la mecanique, ne pas mettre de condition aux limites sur une partie de frontiere, est generalement implicitement interprete comme un partie sur laquelle la contrainte exercee par le milieu exterieur est nulle, c'est a dire (N ) n(N ) = 0. En thermique, la condition similaire correspndrait a une paroi a ux thermique nul, c'est a dire adiabatique, ce qui n'est peut-^ etre pas ce que l'utilisateur non averti sous-entend quand il ne dit rien sur cette frontiere...
J. Garrigues
CHAPITRE 1. INTRODUCTION
Exemple 2: Le probleme elastique lineaire d'evolution

Pour simpli er on ne l'exprime ici que dans le cas des petites perturbations (petits deplacements et petites deformations). Le champ inconnu est un champ vectoriel : le champ de deplacement (M t). Les equations sont : 2 (M t) (Principe fondamental de la mecanique) div + f = @ @t 2 avec = 2 + Tr G (loi de comportement petites deformations et petits deplacements) 1 ;grad (M t) + gradt (M t) (petites deformations) et = 2 Les conditions aux limites sont : { Un champ de deplacements (eventuellement nul), fonction du temps, impose sur la partie de frontiere @ 1 : (N t) = f (N t) 8N 2 @
1
{ Un champ de forces surfaciques (eventuellement nul), fonction du temps, imposes sur la partie de frontiere @ 2 : ou n(N ) est la normale exterieure en N a la frontiere. Compte tenu de la loi de comportement et de la de nition de , la derniere condition aux limites se ramene a des conditions sur les derivees normales de (M t) sur la frontiere. Si le modele physique choisi pour modeliser le probleme reel nous fournit les equations et les conditions aux limites, les mathematiques, par contre, sont souvent impuissantes pour en donner une solution analytique. Les mathematiques prouvent parfois l'existence et l'unicite de la solution, et il est prudent de s'en assurer avant de se lancer dans la recherche d'une solution numerique 2. L'analyse numerique fournit des solutions numeriques. Elles sont plus pauvres qu'une solution analytique. Par exemple, dans le probleme elastique, il faut donner une valeur numerique aux caracteristiques du materiau et , ainsi qu'une de nition numerique des conditions aux limites f (N t) et qf (N t). La solution calculee n'est donc valable que pour ce materiau et pour ces conditions aux limites, alors qu'une solution analytique conserve tous les parametres et permet d'etudier leur in uence sur la solution. Pour avoir une idee de l'in uence de ces parametres, l'ingenieur en est reduit a recommencer le calcul avec d'autres valeurs et estimer leur in uence par comparaison entre les di erents resultats. Pour estimer une solution avec des valeurs de parametres non calcules, on peut tenter d'interpoler (avec une methode d'interpolation arbitraire !) entre les solutions calculees, mais il est generalement imprudent d'extrapoler. Une serie de calculs numeriques ne fournit pas de connaissance au sens strict du terme. En e et, entre les valeurs calculees, les solutions peuvent avoir des variations brusques insoupconnees. Quand on interpole entre deux solutions, on postule une certaine regularite. L'analyse numerique fournit plusieurs methodes de resolution des systemes d'equations di erentielle avec conditions aux limites, parmi lesquelles on peut citer la methode de di erences nies, la methode des volumes nis et en n la methode des elements nis, qui fait l'objet de ce cours.
2. Si la solution n'est pas unique, et que l'algorithme converge vers une solution, il n'est pas possible de se rendre compte de la non unicite. Si on utilise un autre algorithme, il peut tres bien converger vers une autre.
(N t) n (N ) = qf (N t) 8N 2 @
J. Garrigues
Chapitre 2
Grandes lignes de la methode des elements nis

2.1 Expose de la demarche
La methode consiste a rechercher une solution approchee de la solution exacte sous la forme d'un e (M t) de ni par morceaux sur des sous domaines de . Les n sous domaines i doivent champ F ^ etre tels que
n i=1 i=
et e i \ e j = 8i 6= j
ou e i designe l'interieur de i. Autrement dit, les i sont une partition de . ei (M t), de nis sur chaque sous domaines sont des champs choisis parmi une famille Les champs f arbitraire de champs (generalement polyn^ omiaux). La famille de champs locaux est appelee espace des fonctions d'interpolation de l'element. e (M t), obtenus par juxtaposition des champs locaux est appelee La famille de champs globaux F espace des fonctions d'interpolation du domaine . Le champ dans chaque sous domaine i est determine par un nombre ni de valeurs du champ (ou de valeurs de ses derivees) en des points choisis arbitrairement dans le sous domaine, et appeles n uds 1 . Le champ local est une interpolation entre les valeurs aux n uds. Le sous domaine muni de son interpolation est appele element. Chercher une solution par elements nis consiste donc a determiner quel champ local on attribue a e (M t) obtenu par juxtaposition de ces champs chaque sous domaine 2 pour que le champ global F locaux soit proche de la solution du probleme. Parmi les contraintes qu'on impose a la solution approchee cherchee, il y a souvent au moins une continuite simple (C0 ) a la frontiere entre les sous domaines. La gure 2.1 montre une solution approchee discontinue d'un champ scalaire sur un domaine de dimension 1. La famille de champs locaux est la famille des champs constants par morceaux. La gure 2.2 montre une solution approchee continue C0 d'un champ scalaire sur un domaine de dimension 1. La famille de champs locaux est la famille des champs polyn^ omiaux de degre 1.
1. On verra plus loin qu'il vaut mieux choisir au moins une partie des n uds sur la frontiere 2. C'est a dire, quelles valeurs il faut donner aux n uds.
J. Garrigues
CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS
Fig. 2.1 {
Solution approchee discontinue
Fig. 2.2 {
Solution approchee continue C0
Fig. 2.3 {
Solution approchee continue C1
J. Garrigues
2.2. UN EXEMPLE DETAILLE

La gure 2.3 montre une solution approchee continue C1 d'un champ scalaire sur un domaine de dimension 1. La famille de champs locaux est la famille des champs polyn^ omiaux de degre 3. La qualite de la solution approchee depend de la division en sous domaines (nombre et dimensions des sous domaines), du choix de la famille de champs locaux dans chaque sous domaine, et des conditions de continuite qu'on impose aux frontieres des sous domaines (C0 C1 :::). Une fois ces choix faits, il reste a rechercher, une combinaison de champs locaux qui satisfait approximativement les equations. Pour resoudre un probleme par la methode des elements nis, on procede donc par etapes successives : 1. On se pose un probleme physique sous la forme d'une equation di erentielle ou aux derives partielles a satisfaire en tout point d'un domaine , avec des conditions aux limites sur le bord @ . 2. On construit une formulation integrale du systeme di erentiel a resoudre et de ses conditions aux limites: C'est la formulation variationnelle du probleme. 3. On divise en sous domaines : C'est le maillage. Les sous domaines sont appeles mailles. 4. On choisit la famille de champ locaux, c'est a dire a la fois la position des n uds dans les sous domaines et les polyn^ omes (ou autres fonctions) qui de nissent le champ local en fonction des valeurs aux n uds (et eventuellement des derivees). La maille completee par ces informations est alors appelee element. 5. On ramene le probleme a un probleme discret : C'est la discretisation. En e et, toute solution approchee est completement determinee par les valeurs aux n uds des elements. Il su t donc de trouver les valeurs a attribuer aux n uds pour decrire une solution approchee. Le probleme fondamental de la methode des elements nis peut se resumer en deux questions : (a) Comment choisir le probleme discret dont la solution est ((proche)) de la solution exacte? (b) Quelle signi cation donner au mot ((proche)) ? 6. On resout le probleme discret: C'est la resolution 7. On peut alors construire la solution approchee a partir des valeurs trouvees aux n uds et en deduire d'autres grandeurs 3 : C'est le post-traitement. 8. On visualise et on exploite la solution pour juger de sa qualite numerique et juger si elle satisfait les criteres du cahier des charges : C'est l'exploitation des resultats. Les etapes 1,2,3,4 et 5 sont souvent rassemblees sous le nom de pretraitement. Le travail de ces di erentes etapes est assiste par les logiciels. Il reste que pour ma^triser leur utilisation, il est indispensable de comprendre les fondements de la methode, notamment les phases 3 et 4, ne serait-ce que pour comprendre et choisir intelligemment parmi les options qu'ils proposent.
2.2 Un exemple detaille

On se propose de rechercher une solution approchee du probleme suivant: Trouver f (x) dans le domaine = 0 1] satisfaisant l'equation di erentielle d2f (x) + 2 df (x) ; 3 x = 0 dx2 dx
3. par exemple, en elasticite,quand on a obtenu une solution approchee du champ des deplacements, on en deduit le champ des tenseurs de deformation, puis le champ des tenseurs de contraintes
J. Garrigues

avec les conditions aux frontieres de : (x) f (0) = 0 et dfdx =0 x=1 Dans ce probleme, est un domaine de dimension 1 et le temps n'appara^t pas. Pour decrire les points du domaine, on n'a besoin que d'une seule variable x. L'equation a resoudre est donc une equation di erentielle ordinaire. La frontiere @ se reduit a deux points 4. La solution exacte de ce probleme est 3 x2 ; 3 e2 + 3 e2;2 x x + f (x) = ; 3 4 4 8 8 Son graphe est donne gure 2.6 page 15.
2.2.1 Choix du maillage

On divise arbitrairement en trois mailles de m^ eme taille (voir gure 2.4):
1
= 01 3
2 3 1/3
2 = 1 3 3
4 5 2/3
= 2 3 1
6 7 1
1 0
Fig. 2.4 {
Maillage du probleme
2.2.2 Choix des n uds et des champs locaux

On decide de prendre des elements a 3 n uds , et pour la famille de champs locaux des polyn^ omes de degre 2. Les n uds sont choisis aux extremites et au milieu de chaque maille. On peut alors determiner chaque champ local en fonction des valeurs aux 3 n uds. Remarquer que le fait d'avoir utilise des n uds aux extremites de chaque element presente deux avantages: 1. Le nombre de n uds 5 est reduit, car il y a des n uds communs a deux elements. 2. On assure ainsi une continuite C0 de la solution approchee : les champs locaux de deux elements voisins auront la m^ eme valeur a leur n ud commun. Remarquer encore qu'il n'est pas necessaire de prendre les elements identiques : on aurait pu prendre certains a deux n uds et d'autres a trois n uds.
4. Pour des domaines de dimension superieure a 1, on a des equations aux derivees partielles. Si le domaine est de dimension 2, @ est un (ou plusieurs) contour(s), s'il est de dimension 3, @ est une (ou plusieurs) surface(s). Les conditions aux limites sont alors des champs imposes sur @ . 5. et donc le nombre d'inconnues
10
J. Garrigues

Chaque champ local peut donc s'exprimer en fonction des valeurs aux n uds. En e et, il n'existe qu'un seul polyn^ ome du second degre qui satisfasse les conditions suivantes 6 : a1 x2 1 + b1 x1 + c1 = f1 a1 x2 2 + b1 x2 + c1 = f2 a1 x2 3 + b1 x3 + c1 = f3 ou x1, x2 et x3 sont les coordonnees des 3 n uds de l'element, et f1 , f2 et f3 sont des valeurs aux 3 n uds de l'element. On obtient : 2 f3 + x1 f3 ; f1 x3 ; f2 x1 + f2 x3 a1 = x2 f1 ; (x (2.1) x2 ; x3) (x1 ; x3) (x1 ; x2 ) 2 2 f + x 2f ; f x 2 ; f x 2 + f x 2 2 3 2 1 3 3 2 1 2 b1 = ; x1 f3 ; x1 (2.2) (x2 ; x3) (x1 ; x3 ) (x1 ; x2) 2 2 2 2 2 2 3 ; f3 x2 x1 + x3 f2 x1 + f1 x2 x3 ; x3 x2 f1 (2.3) c1 = x1 x2 f3 ; x1 f2 x (x2 ; x3) (x1 ; x3) (x1 ; x2 ) Le polyn^ ome avec les coe cients a1 , b1 et c1 ci-dessus est appele fonction d'interpolation de l'element 1. Les polyn^ omes d'interpolation des 2 autres sont de la m^ eme forme: Pour le second, on remplace les indices 1,2,3 par les indices 3,4,5 et pour le troisieme, on remplace les indices 1,2,3 par les indices 5,6,7. Dans cet exemple, on a 7 n uds. Si on donne une valeur arbitraire a chaque n ud, les fonctions d'interpolation sont determinees et de nissent par morceaux une fonction f~ continue C0 sur le ~ , de nies en 3 morceaux, est de dimension 7. domaine . L'espace des fonctions F Les trois fonctions d'interpolation sur les 3 elements sont : (2.4) f~1 = (18 f3 ; 36 f2 + 18 f1) x2 + (;3 f3 + 12 f2 ; 9 f1 ) x + f1 f~2 = (18 f5 ; 36 f4 + 18 f3) x2 + (;15 f5 + 36 f4 ; 21 f3) x + 3 f5 ; 8 f4 + 6 f3 (2.5) 2 ~ f3 = (18 f7 ; 36 f6 + 18 f5) x + (;27 f7 + 60 f6 ; 33 f5) x + 10 f7 ; 24 f6 + 15 f5 (2.6)
2.2.3 Formulation variationnelle du probleme

On montre en analyse fonctionnelle que
(x) h(x) dx = 0 8 (x) () h(x) = 0 8x 2
ou (x) et h(x) et doivent satisfaire a certaines conditions de regularite qu'on n'explicitera pas ici. Autrement dit, resoudre l'equation di erentielle est equivalent a chercher f (x) tel que : Z 2f (x) df (x) (2.7) (x) d dx 2 + 2 dx ; 3 x dx = 0 8 (x) (x) =0 (2.8) avec f (0) = 0 et dfdx
x=1 6. Cette condition est essentielle: Il faut que le champ local soit completement et uniquement determine par les valeurs imposees aux n uds.
J. Garrigues
11

Cette formulation est une formulation variationnelle (on dit aussi formulation integrale ou encore formulation faible ) du probleme. On obtient une autre formulation variationnelle en utilisant l'integration par parties
Z1
0
Z1
0
f (x) dx = ; (x) d dx 2
2
(x) dx = ; (x) dfdx
Z 1 d (x) df (x) df (x) dx + ( x ) dx dx dx Z01 d (x) 1

0
1 0
dx f (x) dx + (x) f (x)]0
(x) Compte tenu des conditions aux limites, f (0) = 0 et dfdx = 0 , une autre formulation x=1 variationnelle du probleme est donc : Trouver f (x) tel que
Z 1 d (x) df (x) (x) ; dx + 0 ; (0) dfdx dx dx 0 Z 1 d (x)

;2
0
x=0
dx f (x) dx + 2 (1) f (1) ; 0
;3
Z1
0
(x) x dx = 0 8 (x) (2.9)
Il est important de noter que dans les termes de bord de cette derniere formulation, les conditions aux limites sont prises en compte. En comparant les deux formulations variationnelles (2.7) et (2.9), on constate que, dans la seconde forme, la fonction inconnue f (x) intervient avec des derivees d'ordre inferieur, et que par contre les derivees des fonctions (x) apparaissent. L'equivalence entre ces deux formulations variationnelles et le probleme initial est soumise a des conditions sur f (x), (x) et precisees par l'analyse fonctionnelle et qu'on supposera satisfaites. Remarque : La demarche suivie ici pour construire une formulation variationnelle du probleme est purement mathematique. Il arrive souvent que la formulation variationnelle ait une interpretation physique. On peut donc aussi les etablir par des raisonnements physiques 7 .
2.2.4 Discretisation du probleme

Le principe de l'approximation par elements nis est le suivant: On choisit une formulation variationnelle (par exemple (2.9)) et on cherche a satisfaire l'equation (2.9) avec les f~ determines par les valeurs aux n uds de nis precedemment. Il va de soit que ce probleme, n'a generalement pas de solution si on garde la condition ((8 (x))) (sinon f~ serait une solution exacte !). On ne veri era donc cette equation que pour certains i (x) seulement. A chaque i (x) correspondra une equation scalaire fonction des valeurs aux n uds. Il su t de les choisir independants et en nombre egal au nombre d'inconnues. Dans notre exemple, il en faut 7 pour determiner les valeurs aux noeuds f1 ::: f7 8. Les fonctions (x) sont appelees fonctions test ou encore fonctions de ponderation. On est ainsi ramene a un systeme algebrique de 7 equations a 7 inconnues.
7. C'est notamment le cas en mecanique: Le theoreme des puissances virtuelles fournit directement une formulation variationnelle, les champs de vitesses virtuelles n'etant rien d'autre que les fonctions (x). 8. On peut aussi diminuer le nombre d'inconnues en imposant aux f~ de respecter exactement les conditions aux limites. Dans l'exemple, il serait ramene a 5. Si on ne le fait pas, les conditions aux limites ne seront qu'approximativement satisfaites.
12
J. Garrigues

Il y a une in nite de manieres de choisir les 7 i (x). Elles engendrent di erentes variantes de la methode des elements nis. Chacune aboutit a une solution approchee di erente. La seule condition imperative est que le choix des i doit conduire a un systeme algebrique a solution unique. Des theoremes d'analyse numerique suggerent certains choix, pour lesquels la solution approchee f~ est garantie unique et convergente vers la solution exacte 9 quand le maillage se ra ne. Dans notre exemple, si on choisit la formulation (2.9), les fonctions i (x) doivent avoir une derivee non nulle partout, pour que la formulation ne soit pas triviale 10. Par contre, dans cette formulation, ~(x), n'appara^t la derivee seconde de la fonction f (x), qui sera remplacee par l'approximation f plus. On pouvait donc choisir une famille d'approximation f~(x) plus simple11. On choisit par exemple les 7 fonctions test
1 1 2
0
, 1 ,..., 6 de nies par la gure 2.5.

7 1
1/6
2/6
5/6
Fig. 2.5 { Les
sept fonctions i (x)
Leurs derivees sont constantes (+6 ou -6) ou nulles et on a i (xj ) = ij . Pour chaque fonction test i on ecrit une equation scalaire:
" ~ # Z 1 d i(x) df~(x) (x) dx ; i(0) df ; dx dx dx 0 x=0 Z 1 d i (x )

;2
0
~ ~ dx f (x) dx + 2 i (1) f (1)
;3
Z1
0
i (x) x dx = 0
pour i = 1::7
ou les f~(x) sont des fonctions de x et des valeurs aux n uds f1 f2 ::: f7. Les integrales de ces 7 equations sont faciles a calculer : elles se reduisent a des integrales sur le segment ou les derivees de i sont non nulles. La construction de ce systeme d'equations s'appelle l'assemblage. Pour un choix di erent de fonctions i (x), les integrales sur chaque element sont moins simples a calculer 12.
9. 0n ne sait le demontrer que pour certains types de problemes, en particulier les problemes ou l'equation di erentielle est lineaire. Il existe aussi d'autres resultats, mais pour beaucoup de problemes, la methode des elements nis reste une methode heuristique. 10. Alors qu'avec la premiere formulation integrale, on aurait pu prendre des fonctions constantes par morceaux ou m^ eme des distributions de Dirac (methode dite de collocation ). 11. par exemple: des fonctions d'interpolation lineaires par morceaux, avec des elements a deux n uds. La dimene serait 4. On aurait donc moins de ( sion de l'espace des F (liberte) ) pour ajuster la solution approchee,l'approximation trouvee serait plus grossiere. On donne en n de chapitre le graphe du resultat 12. On verra plus loin qu'on peut les calculer de maniere exacte ou approchee.
J. Garrigues
13

Par exemple, la premiere equation (avec 1) est :
" ~ # Z d 1(x) df~(x) df (x) ; dx dx dx ; 1 (0) dx x=0 0 Z d 1(x)

1 6
;2
1 6
~ ~ dx f (x) dx + 2 1(1) f (1)
;3
soit encore :
Z
0
1 6
(x) x dx = 0
Z 16 Z1 ~(x) ~(x) df ~(x) dx ; 3 6 (1 ; 6 x) x dx = 0 ; (;6) df dx ; ; 2 ( ; 6) f dx dx x=0 0 0 0 ou f~(x) = a x2 + b x + c dans le segment 0 1 eme 6 , avec a, b et c de nis en (2.1). On procede de m^ pour les 6 autres equations 13. On obtient ainsi un systeme de 7 equations a 7 inconnues (les valeurs aux n uds f1 f2 ::: f7) : 23 f ; 14 f + 17 f ; 1 = 0 6 1 3 2 6 3 72 1 =0 7 f3 + 5 f1 ; 12 f2 ; 12 22 f + 14 f ; 12 f ; 1=6 f + 1 f ; 1 = 0 3 5 3 4 3 2 6 1 6 7 f5 ; 12 f4 + 5 f3 ; 1 4 =0 22 f + 1 f ; 1 f ; 1 = 0 f ;12 f5 + 14 4+ 3 3 6 6 3 6 7 3 5 =0 5 f5 + 7 f7 ; 12 f6 ; 12 14 f ; 29 f + 1 f ; 17 = 0 3 6 6 7 6 5 72
1 6
"
2.2.5 Resolution
La solution de ce systeme donne les valeurs aux n uds de l'approximation f~(x) cherchee : 4107 f = ; 7063 f = ; 236893 f = ; 25563 f = ; 350779 f = ; 374605 f = ; 127275 f1 = 54872 2 8664 3 164616 4 13718 5 164616 6 164616 7 54872 Les choix faits precedemment (les fonctions tests et les fonctions d'interpolation) conduisent a un systeme d'equations non symetrique. On verra plus loin comment obtenir des systemes symetriques, pour lesquels on dispose d'algorithmes numeriques plus e caces. Les fi etant determinees, on conna^t alors les interpolations dans les trois sous domaines, et donc aussi la solution approchee en juxtaposant ces interpolations.
2.2.6 Examen des resultats

La gure 2.6 montre une comparaison entre la solution exacte et la solution approchee.
~(x) dans les intervalles d'integration! 13. en prenant soin de prendre la bonne de nition de f
14
J. Garrigues

0 -0.5 -1 -1.5 -2 0.2 0.4
0.6
0.8
Fig. 2.6 {
Comparaison des solutions exacte et approchee
On peut noter que puisqu'on n'a pas impose aux fonctions f~(x) de respecter exactement les conditions aux limites, celles-ci ne le sont qu'approximativement: "~ # d(x) f~(0) ' 0:0748 dx x=1 ' 0:0394 Bien que cela n'apparaisse pas clairement sur la gure 2.6, les derivees sont bien discontinues : df~(x) ' ;2:944 lim df~(x) ' ;3:017 lim 1 1 dx x! 3 3 x dx A titre de comparaison on donne gure 2.7 le graphe de la solution approchee avec le m^ eme ~ (x) n'est maillage mais avec des elements a interpolation lineaire a deux n uds. L'espace des F que de dimension 4 et la solution est beaucoup plus grossiere. On donne aussi gure 2.8 le graphe
0.5 0 -0.5 -1 -1.5 -2 0.2 0.4
0.6
0.8
Fig. 2.7 {
Approximation par 3 elements lineaires
avec 3 elements a interpolation lineaire a deux n uds ou on a astreint les interpolations des ~ (x) n'est plus que elements 1 et 3 a respecter exactement les conditions aux limites. L'espace des F de dimension 2. On donne en n gure 2.9 le graphe avec 3 elements a interpolation du second degre a 3 n uds ou on a astreint les interpolations des elements 1 et 3 a respecter exactement les conditions aux ~ (x) est de dimension 5. limites. L'espace des F
J. Garrigues
15
0 -0.5 -1 -1.5 -2
0.2
0.4
0.6
0.8
Fig. 2.8 {
Approximation par 3 elements lineaires,avec C.L. imposees exactes
0 -0.5 -1 -1.5 -2 -2.5
0.2
0.4
0.6
0.8
Fig. 2.9 { Approximation
par 3 elements d'ordre 2,avec C.L. imposees exactes
16
J. Garrigues
Chapitre 3
Le maillage
L'operation de maillage consiste a diviser le domaine en sous domaines appeles mailles. Il faut donc : { De nir le domaine { Le diviser en mailles i telles que l'ensemble des i soit une partition du domaine . Ces operations sont assistees par les logiciels, mais sont rarement completement automatiques.
3.1 De nition du domaine

La de nition du domaine peut ^ etre faite dans le logiciel d'elements nis ou peut ^ etre importee d'un autre logiciel mieux adapte pour ce travail (logiciel de CAO).
3.1.1 De nition du domaine dans le logiciel d'elements nis

Chaque logiciel a ses propres outils de de nition et de visualisation d'entites geometriques. On retrouve neanmoins des concepts communs 1 : notion de repere : (parfois appeles referentiels ). Il existe toujours un repere appele repere global, qui est le repere utilise par defaut quand on n'en de nit pas d'autre. On peut de nir et utiliser d'autres reperes si on les trouve plus commodes pour de nir les objets geometriques qui suivent. notion de systeme de coordonnees : pour designer un point par trois nombres dans un repere, on peut choisir si ces trois nombres sont des coordonnees cartesiennes (la plupart du temps c'est l'option par defaut), cylindriques, spheriques ou autre. notion de point : il peut ^ etre designe par ses coordonnees dans un repere ou par une de nition geometrique 2. notion de ligne : On peut de nir des droites, des segments de droites, des cercles, des arcs de cercle, des coniques, des arcs de coniques, des courbes de Bezier, des splines, des Bsplines, des NURBS 3 ... La liste est plus ou moins riche suivant les logiciels.
1. Mais le vocabulaire et la syntaxe des commandes varient d'un logiciel a l'autre. 2. par exemple le milieu d'un segment de droite existant, ou par l'intersection de deux lignes existantes, etc. 3. Non Uniform Rational B-Spline
J. Garrigues
17
CHAPITRE 3. LE MAILLAGE
de ni les limites. On peut de nir des plans, des cylindres, des spheres, des tores, des c^ ones, des surfaces de Bezier, des surfaces splines, des surfaces Bsplines, des surfaces NURBS... notion de surface delimitee : Il s'agit d'une surface porteuse associee a un ou plusieurs contours. L'un des contours (generalement le premier) est le contour exterieur, les autres sont les contours interieurs (surfaces a ((trous))) Il va de soi que ces contours ne doivent pas s'intersecter ! Une surface peut aussi ^ etre delimitee par son intersection avec d'autres surfaces. Elle peut aussi ^ etre engendree par le deplacement d'une ligne suivant un certain parcours (c'est a dire une autre ligne), ou suivant un certain mouvement (translation, rotation ou autre). notion de volume : C'est une portion limitee de l'espace. Il peut ^ etre de ni par les surfaces qui le delimitent, ou par le deplacement d'une surface suivant un parcours, ou par des instructions particulieres pour des volumes simples ( parallelepipedes, spheres, cylindres, tores ...) . Suivant les logiciels, la palette d'objets geometriques et la palette d'outils (intersections, prolongements, conges automatiques...) est plus ou moins riche et la de nition de domaines complexes est plus ou moins aisee 4 .
notion de surface ((porteuse)) : Il s'agit de surfaces sans bords, c'est a dire dont on n'a pas
3.1.2 Importation du domaine d'un autre logiciel

La de nition d'un domaine de forme complexe est parfois mal commode dans les preprocesseurs de logiciels d'elements nis. Il est quelquefois judicieux de le faire dans un autre logiciel plus adapte. On utilise alors un logiciel de CAO. Le probleme est que chaque logiciel (d'elements nis ou de CAO) a sa propre methode pour decrire la geometrie d'un domaine dans un chier. Un chier issu d'un logiciel de CAO n'est generalement pas compris par un autre logiciel de CAO ou d'elements nis ! Il existe cependant des formats standards normalises de description de geometries. Si on a la chance que le logiciel de CAO et le logiciel d'elements nis soient capables tous les deux de lire et ecrire dans l'un de ces formats standards, on peut communiquer une geometrie de l'un a l'autre. Malheureusement, les normes de ces formats standards evoluent de version en version et les logiciels dont on dispose ne parlent pas toujours la m^ eme version. D'autre part, beaucoup de logiciels prennent des libertes avec la norme ou ne reconnaissent qu'une partie des entites normalisees. Le resultat est que le transfert d'une geometrie d'un logiciel de CAO vers un logiciel d'elements nis se fait rarement sans pertes ! Il faut generalement retoucher la geometrie recue avec les outils du preprocesseur du logiciel d'elements nis.
3.2 De nition des sous domaines

3.2.1 Quelques ((evidences))
Le domaine etant de ni, il faut maintenant le diviser en sous domaines i , appeles mailles. Cette operation est appelee maillage.
{ Si le domaine est un volume, les mailles sont des volumes

4. Les concepteurs de logiciels d'elements nis ont tendance a considerer que cette phase de travail est annexe, et ils n'y apportent pas toujours le soin qui conviendrait pour repondre aux besoins de l'utilisateur. Une des raisons est que les initiateurs d'un projet de logiciel d'elements nis sont souvent des scienti ques specialistes dans le domaine du calcul, et les questions relatives aux interfaces utilisateurs les interessent peu. Les logiciels de CAO sont souvent plus e caces pour la de nition de geometries complexes.
18
J. Garrigues
3.2. DEFINITION DES SOUS DOMAINES
{ Si le domaine est une surface bornee, les mailles sont des surfaces bornees portees par la surface . { Si le domaine est une portion de ligne, les mailles sont des portions de lignes portees par la ligne . Une certaine confusion regne dans le vocabulaire employe en ce qui concerne les problemes quali es de ((3D)), ((2D)) ou ((1D)): { Dans certains contextes, ((2D)) signi e que les points du domaine sont reperes par deux reels sur une surface quelconque plongee dans l'espace physique a 3 dimensions dans d'autres contextes, ((2D)) signi e que l'espace physique possede seulement deux dimensions 5, il s'agit alors de problemes plans. { Dans certains contextes, ((1D)) signi e que les points du domaine sont reperes par un reel sur une ligne quelconque plongee dans l'espace physique a 3 dimensions dans d'autres contextes, ((1D)) signi e que l'espace physique possede seulement une dimension, il s'agit alors de problemes sur une droite. { Dans la plupart des logiciels, il existe une option dite ((axisymetrique)) : Elle signi e que le domaine et les conditions aux limites sont invariants par rotation autour d'un axe xe. Dans ce cas, on fait des economies de maillage: { Pour de nir un volume axisymetrique, on ne de nit qu'une section meridienne du volume, c'est a dire une surface. Faut-il appeler ces problemes ((2D)) ou ((3D))? { Pour de nir une surface axisymetrique, on ne de nit qu'une section meridienne de la surface, c'est a dire une ligne. Faut-il appeler ces problemes ((1D)), ((2D)) ou ((3D))? Dans ce cours on evitera d'employer les quali catifs ((1D)), ((2D)) ou ((3D)).
3.2.2 Les mailles disponibles

Les logiciels d'elements nis disposent de 3 sortes de mailles: Les mailles lineiques : Elles servent a mailler une courbe qui peut representer le domaine lineique (plonge dans un espace physique a 2 ou 3 dimensions) ou une section meridienne d'un domaine surfacique axisymetrique. Certains logiciels n'o rent que des mailles droites. Les mailles surfaciques : Elles servent a mailler une surface qui peut representer le domaine surfacique (plonge dans un espace physique a 2 ou 3 dimensions 6) ou une section meridienne 7 d'un domaine volumique axisymetrique. La plupart des logiciels ne disposent que de triangles et de quadrangles curvilignes portes par la surface. Certains n'o rent que des triangles et des quadrangles a ar^ etes rectilignes 8 . Les mailles volumiques : Les mailles sont des volumes. La plupart des logiciels ne disposent que de tetraedres, pentaedres ou hexaedres curvilignes (voir gure 3.1). Certains n'o rent que des tetraedres, pentaedres ou hexaedres a ar^ etes rectilignes 9 .
5. et donc que toutes les grandeurs physiques associees au probleme sont des tenseurs de cet espace (vecteur deplacement, tenseur plans de deformation, etc). 6. Une surface plongee dans un espace physique de dimension 2 est une surface plane (en physique classique, l'espace physique est non courbe). 7. donc plane 8. Tous les points de la maille ne sont donc pas exactement sur la surface! 9. Tous les points de la maille ne sont donc pas necessairement dans le volume, et il peut exister des points de qui ne sont pas dans une maille.
J. Garrigues
19
Tetrahedre
Pentahedre
Hexaedre
Fig. 3.1 {
Mailles volumiques
On voit donc que, suivant les possibilites des mailleurs, la geometrie de approchee 10 .
est plus ou moins
3.2.3 Les methodes de maillage

La plupart des logiciels proposent di erents algorithmes pour generer plus ou moins automatiquement des maillages. La phase de maillage est sans doute celle ou on passe le plus de temps dans la de nition d'un probleme! Obtenir un ((bon maillage)) resulte d'une certaine experience et d'une certaine intuition sur le resultat du calcul. On peut cependant degager quelques regles generales a suivre: { Les mailles doivent ^ etre ((bien proportionnees)), c'est a dire que le rapport de leur plus grande dimension sur leur plus petite dimension doit ^ etre aussi voisin de 1 que possible. Dans la pratique, on ne devrait pas depasser 5. Ce rapport est appele distorsion 11 de la maille. Dans un maillage surfacique, les mailles ideales sont les triangles equilateraux, et les carres. Dans un maillage volumique, les mailles ideales sont les tetraedres reguliers et les cubes. { Le maillage ne doit pas ^ etre inutilement n. On verra plus loin que plus le maillage est n, plus le calcul est co^ uteux. On a a faire un compromis entre la nesse de la representation geometrique et le co^ ut du calcul. Quand on a une intuition du resultat nal on peut decider de mailler grossierement dans certaines regions et plus nement dans d'autres. Il peut arriver que l'examen des resultats d'un calcul amenent a le recommencer avec un maillage remanie. { Dans certains calculs 12, certains experts pensent qu'on obtient de meilleurs resultats avec des quadrangles pour les maillages surfaciques, et avec des hexaedres pour les maillages volumiques. Il n'est pas toujours evident de satisfaire a ces regles avec les algorithmes de maillage dont on dispose. Parmi les algorithmes qu'on trouve dans les logiciels, on peut citer : le maillage de Delaunay : Cet algorithme est assez general : il peut mailler en triangles toute surface de nie par son contour maille.L'algorithme garantit une faible distorsion des triangles
10. On verra plus loin, en 4.4.2 page 32 que m^ eme si le mailleur respecte exactement la forme des ar^ etes et des surfaces, on fait quand m^ eme une autre approximation geometrique! 11. Il existe d'autres de nitions de la distorsion: pour un maillage surfacique, c'est le rapport entre le rayon du cercle circonscrit et le rayon du cercle inscrit. Pour les maillages volumiques, c'est le rapport entre le rayon de la sphere circonscrite et le rayon de la sphere inscrite. On peut aussi considerer le plus petit angle entre deux ar^ etes ayant un sommet commun, ou encore le plus petit angle entre deux faces ayant une ar^ ete commune. 12. notamment en plasticite
20
J. Garrigues
3.2. DEFINITION DES SOUS DOMAINES
generes. Il existe aussi une version volumique de cet algorithme qui genere des tetraedres 13 . La nesse du maillage est contr^ olee par la nesse de maillage des contours (ou des surfaces 14 ) limites. Certaines versions de cet algorithme permettent aussi de piloter le ra nement dans certaines zones ou autour de certains points a l'interieur du domaine. le maillage ((topologique)) : Le contour d'une surface est divise en quatre lignes, chacune d'elles est maillee de telle maniere que deux lignes opposees aient le m^ eme nombre de mailles. La surface a mailler est alors topologiquement equivalente a un rectangle. Le maillage de la surface est obtenu par ((transport)) d'un maillage regulier du rectangle sur la surface. La methode de ((transport)) peut ^ etre une simple projection ou une methode plus elaboree. La m^ eme methode existe en considerant la surface topologiquement equivalente a un disque 15. Ces algorithmes existent aussi en version volumique. Le volume topologiquement equivalent peut ^ etre un parallelepipede, un cylindre ou une sphere. Ces methodes ne garantissent rien sur la distorsion des mailles obtenues, et il faut souvent proceder par t^ atonnements successifs. le maillage ((par deplacement)) : On impose a une ligne maillee de se ((deplacer)) le long d'un chemin maille ou suivant un mouvement simple (translation ou rotation). On genere ainsi une maillage de quadrangles. La m^ eme methode est employee pour mailler des volumes: on ( (deplace) ) une surface maillee. La encore, ces methodes ne garantissent rien sur la distorsion des mailles obtenues, et il faut souvent proceder par t^ atonnements successifs. le maillage par sous domaines : On divise le domaine (volume ou surface) a mailler en sous domaines de forme plus simple, faciles a mailler avec l'une des methodes ci-dessus. La principale di culte de cette methode est de s'assurer de la coherence du maillage aux interfaces entre les sous domaines: les maillages aux interfaces entre sous domaines doivent coincider. Ce n'est pas evident a obtenir quand on utilise des algorithmes di erents dans des sous domaines voisins. le maillage manuel : Si aucun des algorithmes dont on dispose ne donne satisfaction, c'est le dernier recours ! On de nit les mailles individuellement. C'est parfois necessaire dans certains sous domaines, ou pour ((raccorder)) des sous domaines. Comme on le voit, le maillage d'un domaine de forme complexe n'est pas trivial. On ne le reussit que rarement du premier coup. Il faut essayer de combiner les outils de maillage dont on dispose pour obtenir quelque chose de satisfaisant.
3.2.4 L'importation d'un maillage

Si les outils de maillage du preprocesseur sont insu sants, on peut tenter de le faire dans d'autres logiciels. Notamment, un certain nombre de logiciels de CAO ont des fonctionnalites de maillage. Les problemes de communication entre maillages issus de logiciels di erents sont les m^ emes que pour les echanges de geometries : il faut qu'il existe un format de chier connu par les deux logiciels, et qu'ils utilisent des versions voisines. On peut aussi developper soi-m^ eme un ((traducteur)) 16 .
13. Actuellement, on ne conna^t pas d'algorithme general qui genere des quadrangles ou des hexaedres garantis bien proportionnes. 14. qui peuvent elle-m^ emes ^ etre un maillage de Delaunay. 15. Il faut alors choisir un point dans la surface comme ((centre)) de la surface. Cette methode est parfois appelee maillage rayonnant, car le disque de reference est maille par des rayons et des cercles concentriques. 16. C'est plus facile pour les maillages que pour les geometries, car la description d'un maillage est plus simple que celle d'une geometrie: l'information se reduit a la liste des coordonnees des sommets, et les mailles sont une liste de numeros de sommet (dans un certain ordre, qui varie, bien sur, d'un logiciel a l'autre).
J. Garrigues
21
22
J. Garrigues
Chapitre 4
Les elements et leur espace de fonctions d'interpolation

e de Les solutions approchees trouvees par la methode des elements nis sont une juxtaposition F e champs locaux fi de nis dans chaque maille. L'objet de ce chapitre est de monter comment on de nit dans les mailles une famille de champs locaux dans laquelle la methode des elements nis pourra choisir pour s'approcher au mieux de la solution. Pour qu'une maille devienne un element, il faut: { Choisir arbitrairement des n uds dans la maille. Comme on le verra plus loin, la resolution d'un probleme par la methode des elements nis se ramene a calculer les valeurs de la solution approchee aux n uds du maillage. Le nombre d'inconnues scalaires par n ud varie selon la nature tensorielle du champ inconnu et selon la dimension de l'espace physique. Ce nombre est appele nombre de degres de liberte 1 (abreviation courante : ddl). Par exemple, si le champ cherche est un champ vectoriel dans un espace physique a 3 dimensions, on a 3 ddl par n uds (1 ddl par composante du champ) { Choisir arbitrairement une famille de champs locaux, destines a donner une valeur approchee de la solution en tout point de la maille en ne connaissant que les valeurs d'une solution approchee aux n uds de l'element. Cette famille de fonctions s'appelle l'espace des fonctions d'interpolation de la maille. S'il y a plusieurs inconnues par n ud (c'est a dire plusieurs degres de liberte), on construit des interpolations pour chacun des degres de liberte. En general, les interpolations de chaque degre de liberte appartiennent au m^ eme espace de fonctions d'interpolation 2. Dans ce chapitre, on va montrer comment on construit les interpolations. On le fera en supposant qu'il n'y a qu'un seul ddl par n ud. S'il y en a plusieurs on procede de m^ eme pour chaque ddl. Pour construire les espaces de fonctions d'interpolation sur un element, on va proceder en deux temps : on construit un espace de fonctions d'interpolation fer sur une maille de reference standard topologiquement equivalente a la maille reelle, puis on le transforme pour qu'il devienne un espace de fonctions d'interpolation fe sur les mailles reelles. Ce procede a l'avantage de faire gagner du
1. Cette denomination est un heritage des precurseurs de la methode des elements nis: c'etaient des mecaniciens et les inconnues etaient donc les deplacements aux n uds. Ce n'est que plus tard que cette methode a ete utilisee pour d'autres problemes de physique. 2. mais ce n'est pas obligatoire: certains calculs en mecanique des uides utilisent des familles d'interpolation di erentes pour le champ des vitesses et le champ des pressions.
J. Garrigues
23
CHAPITRE 4. LES ELEMENTS ET LEUR ESPACE DE FONCTIONS D'INTERPOLATION

temps calcul. Les logiciels proposent des bibliotheques d'elements dans lesquels les interpolations dans les mailles de reference sont deja de nies et n'ont plus a ^ etre recalculees.
4.1 De nition des mailles de reference

Dans un maillage, toutes les mailles ont des formes et des dimensions di erentes. On trouve des mailles lineiques, des mailles surfaciques et des mailles volumiques de toutes formes et de toutes tailles. Dans le but d'uniformiser et d'automatiser les calculs, on introduit la notion de maille de reference. On ne donne ici que les mailles les plus classiques. Par convention generalement admise (voir gure 4.1):
La maille de reference lineique est le segment x1 2 ;1 1]. La maille de reference surfacique triangulaire est le triangle
x1 0 x2 0 x1 + x2 1.
La maille de reference surfacique quadrangulaire est le carre x1 2 ;1 1] x2 2 ;1 1] La maille de reference volumique tetraedrique est le tetraedre
x1 0 x2 0 x3 0 x1 + x2 + x3 1.
La maille de reference volumique pentaedrique est le pentaedre (prisme a base triangulaire) x1 0 x2 0 x1 + x2 1 x3 2 ;1 1]. La maille de reference volumique hexaedrique est le cube x1 2 ;1 1] x2 2 ;1 1] x3 2 ;1 1].
ou x1, x2 et x3 sont les coordonnees d'un point courant m de la maille de reference.
x2 1 -1 0 1 x1 0 x3 1 0 1 x1 x1 1 x2 1 x1 -1 x3 1 0 1 -1 1 x2 -1 x3 1 -1 1 0 -1 x1 1 x2 -1 1 x1 x2 1
Fig. 4.1 { Mailles
de reference
J. Garrigues
24
^ 4.2. INTERPOLATION POLYNOMIALE SUR LES MAILLES DE REFERENCE
4.2 Interpolation polyn^ omiale sur les mailles de reference

Les fonctions d'interpolations utilisees dans les logiciels sont pratiquement toujours des polyn^ omes. Si la fonction a interpoler est vectorielle ou tensorielle, on interpole de la m^ eme maniere chacune er des composantes. On est donc ramene a des interpolations de champs scalaires. Dans la suite, f represente donc soit l'interpolation d'un champ scalaire soit l'interpolation d'une composante de champ vectoriel ou tensoriel. { Pour les mailles lineiques fer est un polyn^ ome de x1, { Pour les mailles surfaciques fer est un polyn^ ome de x1 et x2, { Pour les mailles volumiques fer est un polyn^ ome de x1 , x2 et x3.
4.2.1 Rappels sur les bases de polyn^ omes

L'espace des polyn^ omes de degre d est un espace vectoriel dont la dimension depend du degre des polyn^ omes et du nombre de variables. Cet espace possede une base canonique constituee de tous les mon^ omes de degre non negatif inferieur ou egal a d. Par exemple : { La base canonique des polyn^ omes a une variable x1 de degre 3 est constituee des mon^ omes 3 . 1 x1 x2 x 1 1 { La base canonique des polyn^ omes a deux variables x1 et x2 de degre 2 est constituee des 2 mon^ omes 1 x1 x2 x2 1 x2 x1x2 . { La base canonique des polyn^ omes a trois variables x1, x2 et x3 de degre 2 est constituee des 2 2 mon^ omes 1 x1 x2 x3 x2 1 x2 x3 x1x2 x2x3 x3 x1 . Le tableau qui suit donne les dimensions des espaces de polyn^ omes pour des degres de 1 a 5 pour les polyn^ omes a 1, 2 ou 3 variables. degre 1 variable 2 variables 3 variables 1 2 3 4 2 3 6 10 3 4 10 20 4 5 15 35 5 6 21 56 A partir de la base canonique, on peut engendrer une in nite de bases : Si n est la dimension de l'espace de polyn^ omes, toute matrice reguliere n n de nit une autre base. Par exemple les 6 polyn^ omes P1 P2 P3 P4 P5 P6 sont une base de l'espace des polyn^ omes de degre 2 a 2 variables:
2P 1 6 P2 6 6 P3 6 6 P 6 4 P4
P6
5
3 2p 11 7 6 p21 7 6 7 6 p31 7 =6 7 6 p41 7 5 6 4p

51
p61
p12 p22 p32 p42 p52 p62
p13 p23 p33 p43 p53 p63
p14 p24 p34 p44 p54 p64
p15 p25 p35 p45 p55 p65
p16 p26 p36 p46 p56 p66
32 7 6 7 6 7 6 7 6 7 6 7 56 4
1 x1 x2 x2 1 x2 2 x1 x2
3 7 7 7 7 7 7 5
si la matrice des pij est reguliere 3 . Si on extrait d'une base une sous base de m polyn^ omes de base, on engendre un sous espace vectoriel de dimension m.
3. Certaines bases ou sous bases avec des proprietes particulieres ont recu des noms: Gegenbauer, Hermite, Jacobi, Laguerre, Legendre, Tchebiche ...
J. Garrigues
25
4.2.2 Interpolation polyn^ omiale sur une maille de reference

Soit Er l'espace de reference, et soit m un point courant de Er . Suivant le type de maille (lineique, surfacique ou volumique), le point m a 1, 2 ou 3 coordonnees. On note d la dimension de l'espace de reference 4 . Pour de nir une interpolation, on choisit arbitrairement n points (n > d) dans l'espace de reference appeles n uds de reference qu'on note m(j ) j = 1 : : : n. Une interpolation polyn^ omiale fer sur l'element de reference Er est de nie par :
er (m) = f er (m(j)) = telle que f
n X
i=1 n X i=1
u(i)P (i)(m) u(i)P (i)(m(j ) ) = u(j ) 8j 2 1 n]
(4.1) (4.2)
ou les P (i)(m) sont des polyn^ omes de dr variables, et ou les u(j ) sont n valeurs donnees aux n uds ( j ) m . Autrement dit, une interpolation fer attribue une valeur fer (m) a tout point m, et sa valeur aux n n uds m(j ) est u(j ). La de nition (4.1) montre que fer (m) appartient a un espace de polyn^ omes engendre par les n ( i ) er est de dimension polyn^ omes P . Si ces polyn^ omes sont une base, l'espace des interpolations f n. C'est generalement ce qu'on souhaite : a chaque n-uplet de valeurs aux n uds donnees u(i) , on veut une interpolation di erente. Dans une interpolation polyn^ omiale a n n uds, les P (i) doivent ^ etre une n-base de polyn^ omes. Soit fBk (m)g une n-base de polyn^ omes connue. On peut donc ecrire : P (i)(m) =
n X k=1 i) a( k Bk (m)
i) ou les a( k sont le terme general d'une matrice reguliere. La condition (4.2) implique les n2 egalites:
P (i)(m(j ) ) = ij soit encore :

n X k=1 i) (j ) a( k Bk (m ) = ij
(4.3)
ou encore sous forme matricielle:
2 6 6 6 6 6 6 6 4
(1) a(1) : : : a(1) 1 k : : : an .. .. .. . . . ( i) (i) (i) a1 : : : ak : : : an .. .. .. . . . (n) (n) (n) a1 : : : ak : : : ak
32 3 2 1 ::: 0 ::: 0 3 B1 (m(1) ) : : : B1 (m(j ) ) : : : B1 (m(n) ) 7 7 6 6 .. .. .. .. .. .. 7 7 7 6 6 . . . . . . 7 7 7 6 7 6 7 (1) ( j ) ( n ) 7 6 7 6 = 0 : : : 1 : : : 0 B ( m ) : : : B ( m ) : : : B ( m ) 7 k k k 7 6 7 6 7 7 6 7 6 . . . . . . 7 . . . . . . 5 4 4 . 5 . . . . . 5 (1) (j ) (n)

Bn (m ) : : : Bn (m ) : : : Bn (m ) 0 ::: 0 ::: 1
4. C'est aussi la dimension de la maille reelle.
26
J. Garrigues

ce qui montre que la matrice P ] des coe cients des polyn^ omes P (i) sur la base fBk g est l'inverse de la matrice C ] des valeurs des polyn^ omes de la base Bk aux n uds. La matrice C ] doit donc ^ etre reguliere 5 et les coe cients des polyn^ omes de base de l'interpolation sont donnes par: P ] = C ];1 Ainsi, les coe cients des polyn^ omes de base de l'interpolation dans l'espace de reference sont determines par les coordonnees des n uds dans la maille de reference 6 et par le choix de la base fBk (m)g. Le calcul de P ] est tres simple: il su t de construire C ] de terme general Ckj = Bk (m(j ) ) et de l'inverser. er sur une maille de reference est caracterisee par : En resume, une interpolation f { la dimension de l'espace de reference (lineique, surfacique, volumique), c'est a dire le nombre de coordonnees dans l'espace de reference { le choix de n n uds par leurs coordonnees dans l'espace de reference { le choix de la base fBk (m)g On en deduit la matrice P ] des coe cients des polyn^ omes de base de l'interpolation. ( i ) Les polyn^ omes P (polyn^ omes de base de l'interpolation) sont appeles fonctions de forme.
4.2.3 Exemple 1 (standard)

On prend une maille surfacique (les variables sont donc x1 et x2 ), et on veut construire une interpolation a 4 n uds (voir gure 4.2) dont les coordonnees sont:
(1) m(1) = (x(1) 1 x2 ) (2) m(2) = (x(2) 1 x2 ) (3) m(3) = (x(3) 1 x2 ) (4) m(4) = (x(4) 1 x2 )
x2 m
(4)
(3)
x1
(1) m (2) m
Fig. 4.2 { Interpolation
a 4 n uds dans un element de reference surfacique

4 X
L'interpolation est donc de la forme:
i=1 5. ce qui signi e geometriquement que les k points doivent ^ etre distincts, non alignes dans une maille surfacique non colineaires dans une maille surfacique et non coplanaires dans une maille volumique. 6. En general ils sont standardises. Ce ne sont donc pas de vraies variables.
J. Garrigues
er (x1 x2) = f
u(i)P (i) (x1 x2)
27

Les fonctions de forme P (i) sont une base de polyn^ omes a deux variables de dimension 4. Le tableau page 25 montre qu'on peut travailler dans un sous espace des polyn^ omes de degre 2. On peut choisir comme sous base canonique fBk g les quatre polyn^ omes (B1 = 1 B2 = x1 B3 = x2 B4 = x1x2) de la base canonique 7. Les conditions (4.2) s'ecrivent : P (1)(x(1) 1 (2) (1) P (x1 P (3)(x(1) 1 (4) P (x(1) 1 x(1) 2 )=1 x(1) 2 )=0 x(1) 2 )=0 x(1) 2 )=0 P (1)(x(2) 1 (2) (2) P (x1 P (3)(x(2) 1 (4) P (x(2) 1 x(2) 2 )=0 x(2) 2 )=1 x(2) 2 )=0 x(2) 2 )=0 P (1)(x(3) 1 (2) (3) P (x1 P (3)(x(3) 1 (4) P (x(3) 1 x(3) 2 )=0 x(3) 2 )=0 x(3) 2 )=1 (3) x2 )=0 P (1)(x(4) 1 (2) (4) P (x1 P (3)(x(4) 1 (4) (4) P (x1 x(4) 2 )=0 x(4) 2 )=0 x(4) 2 )=0 (4) x2 ) = 1
(i) (i) i) (i) omes P (i) si la En posant P (i) = a( 1 + a2 x1 + a3 x2 + a4 x1x2, on construit une 4-base de polyn^ (i) matrice des aj est reguliere. Les conditions (4.2) s'ecrivent alors:
2 (1) (1) (1) (1) 3 2 1 a1 a2 a3 a4 6 7 6 (1) a(2) a(2) a(2) a(2) 6 1 2 3 4 7 6 x1 6 7 6 (1) 4 a(3) a(3) a(3) a(3) 1 2 3 4 5 4 x2
a(4) a(4) a(4) a(4) 1 2 3 4
(1) x(1) 1 x2
1 x(2) 1 x(2) 2 (2) x(2) 1 x2
1 x(3) 1 x(3) 2 (3) x(3) 1 x2
1 x(4) 1 x(4) 2 (4) x(4) 1 x2
3 2 3 1 0 0 0 7 6 0 1 0 07 7 6 7 = 7 4 0 0 1 05 5
0 0 0 1
On en deduit les coe cients des fonctions de forme par inversion
2 (1) (1) (1) (1) 3 2 1 a1 a2 a3 a4 6 7 6 x(1) a(2) a(2) a(2) a(2) 6 6 1 1 2 3 4 7 = 6 7 6 4 a(3) 4 x(1) a(3) a(3) a(3) 2 1 2 3 4 5
a(4) a(4) a(4) a(4) 1 2 3 4
(1) x(1) 1 x2
1 x(2) 1 x(2) 2 (2) x(2) 1 x2
1 x(3) 1 x(3) 2 (3) x(3) 1 x2
3;1 1 7 x(4) 7 1 7 (4) x2 5 (4) x(4) x 1 2
Par exemple, si on prend comme n uds les points (;1 ;1) (1 ;1) (1 1) (;1 1) dans l'espace de reference :
2 (1) (1) (1) (1) 3 2 3 2 1 ;1 ;1 1 3 a1 a2 a3 a4 1 1 1 1 ;1 6 (2) (2) (2) (2) 7 6 16 a1 a2 a3 a4 7 ;1 1 1 ;1 7 1 1 ;1 ;1 7 6 6 7 6 7 = = 6 7 (3) (3) (3) (3) 4 5 4 4 a1 a2 a3 a4 5 ;1 ;1 1 1 4 1 1 1 15
a(4) a(4) a(4) a(4) 1 2 3 4 1 ; 1 1 ;1 1 ;1 1 ;1
Les quatre fonctions de forme de cet element sont donc :
2 P (1) 3 2 1 ;1 ;1 1 3 2 1 3 6 1 1 ;1 ;1 7 6 x1 7 P (2) 7 = 1 6 6 4 P (3) 7 5 46 41 1 1 17 56 4 x2 7 5

P (4) 1 ;1 1 ;1 x1 x2
7. Il s'agit bien d'un choix! Quatre polyn^ omes de degre 2 quelconques mais independants peuvent convenir. On s'arrange generalement pour garder des mon^ omes de la base canonique tels que les r^ oles des variables aient une certaine ((symetrie)). L'interpolation obtenue presente alors une certaine ((isotropie)). On pouvait prendre aussi une sous base de l'espace des polyn^ omes de degre 3.
28
J. Garrigues

1 (1 ; x ; x + x x ) P (1)(x1 x2) = 4 1 2 1 2 1 P (2)(x1 x2) = 4 (1 + x1 ; x2 ; x1x2 ) P (3)(x1 x2) = 1 4 (1 + x1 + x2 + x1x2 ) P (4)(x1 x2) = 1 4 (1 ; x1 + x2 ; x1x2 ) et l'interpolation en fonction des valeurs aux n uds est : er (x1 x2) = u(1)P (1)(x1 x2) + u(2)P (2)(x1 x2) + u(3)P (3)(x1 x2) + u(4)P (4)(x1 x2) f soit encore :
4.2.4 Exemple 2 (non standard)

On construit une interpolation a 7 n uds dans une maille surfacique (voir gure 4.3).
x2
2/3 1/3 (5) (4) m m (6) (7) (3) m m m
m
0 1/3
(1)
m
2/3
(2)
Fig. 4.3 { Interpolation
a 7 n uds dans un element de reference surfacique
L'espace d'interpolation est de dimension 7. Le tableau page 25 montre qu'il faut prendre au minimum des polyn^ omes de degre 3. On peut prendre, par exemple, comme sous base fBk g de la base canonique la suite de mon^ omes suivante : 2 2 2 B1 = 1 B2 = x1 B3 = x2 B4 = x2 1 B5 = x2 B6 = x1 x2 B7 = x1x2 . La matrice C ] de terme general Ckj = Bk (m(j ) ) est : 21 1 1 1 1 1 13 1 2 2 1 7 6 0 0 1 3 3 3 3 3 7 6 1 2 2 1 1 7 6 0 0 3 3 3 3 7 61 4 3 4 1 7 0 0 1 C] = 6 9 9 9 9 9 7 6 1 4 4 1 1 7 6 0 0 9 9 9 9 9 7 6 4 2 1 5 4 0 0 27 0 0 27 27 2 4 1 0 0 27 0 0 27 27 dont l'inverse est :
2 5 4 6 ;1 4 6 1 6 6 4 1 P] = 6 6 4 6 ;1 6 4 4 5 4
;3 2
; ; ; ;
Ce qui donne les coe cients des fonctions de forme.

J. Garrigues
3 8 3 8 9 8 9 8 9 8 45 8 27 4
; 45 ; 27 8 8 ; ; ;
9 8 9 8 9 8 3 8 3 8 27 4
; ;
27 8 9 8 9 8 9 8 45 8 27 4
; ; ;
45 8 9 8 9 8 9 8 27 8 27 8 27 4
; ; ; ;
81 8 81 8 81 8 27 8 27 8 27 8 27 4
3 ; 27 8 7 7 ; 7 7 7 7 ; 7 7 5
;
27 8 27 8 81 8 81 8 81 8 27 4
29
4.3 Continuite C0 inter-elements : Les elements conformes

On souhaite generalement que la solution approchee trouvee presente une continuite C0 interelement. Les elements qui ont cette propriete sont appeles elements conformes. Pour les mailles lineiques, la continuite C0 est facile a assurer : Il su t de mettre des n uds aux extremites de la mailles. On est ainsi assure que les fonctions d'interpolation prennent la m^ eme valeur au n ud commun. Dans la suite on ne parlera que de la conformite des elements surfaciques et volumiques. On remplit partiellement cette condition en choisissant de mettre des n uds sur la frontiere des elements 8. On garantit ainsi que la solution approchee pour deux elements voisins est continue au moins aux n uds communs sur la frontiere commune. Si on choisit les fonctions d'interpolation telles que : la valeur de l'interpolation sur une frontiere ne depend que des n uds de la frontiere , on garantit que la solution approchee est continue C0 a la traversee de la frontiere 9 . Autrement dit, sur une ar^ ete d'element conforme surfacique, l'interpolation doit se reduire a une interpolation de maille lineique sur les n uds de l'ar^ ete et sur une face d'element conforme volumique, l'interpolation doit se reduire a une interpolation de maille surfacique sur les n uds de la face. Considerons un element de reference et deux frontieres voisines de cet element, notees Fi et Fj . Si l'element est conforme, la valeur de l'interpolation sur Fi ne depend que des n uds de Fi . De m^ eme pour Fj . Sur l'intersection des frontieres, l'interpolation ne doit dependre que des n uds communs aux deux frontieres. Il faut donc qu'elles aient des n uds communs. On en deduit qu'un element de reference conforme a necessairement des n uds a ses sommets. On voit par exemple que si les n uds de l'exemple 2 page 29 sont ceux d'un element triangulaire tel que celui de la gure 4.1 page 24, l'interpolation de cet element n'est pas conforme. Mais cette condition n'est pas su sante. Il faut donc veri er la conformite des interpolations, c'est a dire : Soit F une frontiere de l'element de reference. On note IF l'ensemble des n uds appartenant a la frontiere F , et on note IF le reste des n uds. Soit mF un point appartenant a la frontiere. La conformite de l'element s'ecrit (pour chaque frontiere): P (i) (mF ) = 0 8i 2 IF 8mF 2 F Pour construire des elements conformes, il faut donc adjoindre les conditions ci-dessus 10 aux conditions d'interpolation 4.3 page 26. Par exemple, si les quatre n uds de l'exemple 1 precedent sont les sommets d'un element de reference carre, on veri e facilement que l'interpolation sur cet element est bien conforme. Les principaux elements de reference conformes qu'on trouve dans les logiciels sont : { les elements lineiques dont deux de leurs n uds sont les extremites { les triangles a 3 n uds (sommets) { les triangles a 6 n uds (sommets + milieux des ar^ etes)
8. avec une regle de position des n uds telle que deux elements voisins aient des n uds confondus sur leur frontiere commune. 9. La juxtaposition des interpolations des elements conformes est alors une fonction C0 de nie sur , et qui appartient donc a l'espace de Sobolev H1 ( ). C'est pour des solutions approchees appartenant a cet espace que les theoremes de convergence de la methode sont etablis. 10. La condition m 2 F est generalement une relation lineaire entre les coordonnees, car les frontieres des elements de reference sont des droites ou des plans. Le polyn^ ome P (i) ((m ) est donc un polyn^ ome a une ou deux variables qui doit ^ etre identiquement nul, c'est a dire que les coe cients de tous ses mom^ omes sont nuls.
F F
30
J. Garrigues
4.4. INTERPOLATION DANS L'ELEMENT REEL

{ les quadrangles a 4 n uds (sommets) { les quadrangles a 8 n uds (sommets + milieux des ar^ etes) { les tetraedres a 4 n uds (sommets) { les tetraedres a 10 n uds (sommets + milieux des ar^ etes) { les tetraedres a 20 n uds (sommets + n uds au tiers des ar^ etes) { les pentaedres a 6 n uds (sommets) { les pentaedres a 15 n uds (sommets + milieux des ar^ etes) { les pentaedres a 24 n uds (sommets + n uds au tiers des ar^ etes) { les hexaedres a 8 n uds (sommets) { les hexaedres a 20 n uds (sommets + milieux des ar^ etes) { les hexaedres a 27 n uds (sommets + milieux des ar^ etes + milieux des faces + n ud central) { des hexaedres a 32 n uds (sommets + n uds au tiers des ar^ etes) On trouve dans la litterature specialisee un grand nombre d'elements tout calcules. Le choix dans les logiciels est plus ou moins riche. L'utilisateur a seulement a se preoccuper de conna^tre le degre d'interpolation utilise dans les elements qu'il choisit, pour estimer sa ((souplesse)) pour approcher une solution exacte : Pour un maillage donne11 il est plus facile d'approcher une solution exacte avec des polyn^ omes de degre 2 qu'avec des interpolations lineaires. Pour bien approcher la solution avec des interpolations lineaires, il faut un maillage plus n.
4.4 Interpolation dans l'element reel

Les coordonnees des n uds de l'element reel sont di erentes de celles de l'element de reference. Il convient donc de construire des fonctions d'interpolation fe sur l'element reel en utilisant les interpolations fer construites sur les elements de reference.
4.4.1 Principe de construction

Soient m un point de l'element de reference, et M un point de l'element reel. er (m) la fonction d'interpolation connue sur l'element de reference, et fe(M ) la fonction Soient f d'interpolation cherchee sur l'element reel. Considerons une transformation inversible qui a tout point m de l'element de reference associe un point M de l'element reel : : m 2 element de reference telle que les n uds se correspondent : M (i) = (m(i) )
! M = (m) 2 element reel
e de nie par : Considerons l'application f
e : M 2 element reel ;! fe(M ) = fer (m) 8M = (m) f

11. C'est a dire pour une division en sous domaines donnee.

J. Garrigues
31

On a donc : soit encore :
e(M ) = fe( (m)) = fer (m) 8m 2 element de reference f e f er () fe = fer =f

n X
;1
e est une interpolation sur l'element reel car : L'application f e(M ) = fer (m) = f
= et
n X i=1 i=1 u(i)P (i)( ;1 (M ))
u(i) P (i)(m)
e(M (i)) = fer (m(i) ) = u(i) f er est une interpolation polyn^ e n'en est pas une en general, car On doit noter que si f omiale, f P (i)( ;1 (M )) n'est pas un polyn^ ome en M en general.
4.4.2 Choix de la transformation : l'approximation geometrique

On cherche idealement une transformation qui fait correspondre les points m de la maille de reference aux points M de la maille reelle, telle que 1. la transformation est inversible 2. la transformation fait correspondre les frontieres 3. la transformation fait correspondre les n uds : M (i) = (m(i) ) Malgre le nombre de conditions imposees a , il y a une in nite de solutions a ce probleme 12, mais ces transformations sont compliquees a trouver, surtout lorsque la geometrie de l'element reel est compliquee. On va donc chercher des transformations plus simples. e soit une interpolation est la correspondance des n uds. On La condition imperative pour que f va choisir des transformations qui ne respectent que cette condition. La consequence, est que si on appelle M 0 le transforme de m, en general M 0 n'appartient pas a l'element reel, sauf si m est un n ud. Cela signi e que : { Pour les elements lineiques, la courbe des points M 0 = (m) est di erente de la courbe reelle de l'element. { Pour les elements surfaciques, la surface des points M 0 = (m) est di erente de la surface reelle de l'element, ainsi que la courbe des ar^ etes. { Pour les elements volumiques, les ar^ etes et les faces transformees sont di erentes des ar^ etes et des faces reelles.
( (
12. En e et, on peut imaginer une in nite de transformations qui respectent les frontieres et les n uds et qui deforment)) l'interieur.
J. Garrigues
32
4.4. INTERPOLATION DANS L'ELEMENT REEL

La transformation transforme l'element de reference en une approximation geometrique de l'element reel. Les n uds de l'approximation geometrique et de l'element reel sont confondus, mais les autres points sont distincts. e = fer ;1 est une interpolation sur une approximation geometrique de la maille L'interpolation f reelle. Autrement dit, dans un calcul par elements nis, les mailles reelles sont remplacees par leur approximation geometrique, et cette approximation geometrique est determinee par le choix de . On cherche donc une fonction inversible qui fasse correspondre les n uds, c'est a dire telle que : M (i) = (m(i) ) 8n ud i
Ce probleme est un probleme d'interpolation sur l'element de reference : Choisir , c'est choisir l'interpolation les points M 0 entre les points M (i) . Le probleme est donc : e1 , X e2 et X e3 telles que : pour les mailles lineiques , trouver trois fonctions d'interpolation X
(i) e1(x(1i)) X1 =X (i) e2(x(1i)) X2 =X (i) e3(x(1i)) X3 =X
e1 , X e2 et X e3 telles que : pour les mailles surfaciques , trouver trois fonctions d'interpolation X
(i) e1(x(1i) x(2i)) X1 =X (i) e2(x(1i) x(2i)) X2 =X (i) e3(x(1i) x(2i)) X3 =X
e1, X e2 et X e3 telles que : pour les mailles volumiques , trouver trois fonctions d'interpolation X
(i) e1(x(1i) x(2i) x(3i)) X1 =X (i) e2(x(1i) x(2i) x(3i)) X2 =X (i) e3(x(1i) x(2i) x(3i)) X3 =X
i) eme ou Xj(i) est la j eme coordonnee dans l'espace physique du n ud reel i , et ou x( j est la j coordonnee dans l'espace de reference du n ud de refference i. Si on prend des interpolations polyn^ omiales, les fonctions d'interpolation des coordonnees sont de la forme 13 :
e1(x1 x2 x3) = X1 = X e2(x1 x2 x3) = X2 = X e3(x1 x2 x3) = X3 = X
n X i=1 n X
(i) (i) X1 Q (x1 x2 x3) (i) (i) Q (x1 x2 x3) X2 (i) (i) X3 Q (x1 x2 x3)
(4.4) (4.5) (4.6)
i=1 n X i=1
ou les Q(i) sont les polyn^ omes de base de l'interpolation geometrique, i 2 1 ::: n], et ou n est la dimension de l'espace des fonctions d'interpolation (c'est a dire le nombre de n uds). Cette interpolation doit ^ etre conforme, pour que les frontieres de deux elements voisins soient les m^ emes 14. Pour trouver ces interpolations on emploie les m^ emes techniques que precedemment, puisque ce ne sont que des interpolations sur l'element de reference. Quelques de nitions : { Si la methode d'interpolation des coordonnees est la m^ eme que celle des valeurs aux n uds on dit que l'element est isoparametrique. On a alors P ] = Q]
13. On les donne ici pour les elements volumiques, le lecteur les simpli era aisement pour les elements surfaciques ou lineiques. 14. En e et, les sous domaines de l'approximation geometrique doivent respecter les m^ emes regles que les sous domaines de
0
J. Garrigues
33

{ Si l'interpolation des coordonnees est de degre inferieur au degre d'interpolation des valeurs aux n uds on dit que l'element est subparametrique. { Si l'interpolation des coordonnees est de degre superieur au degre d'interpolation des valeurs aux n uds on dit que l'element est superparametrique. Il reste un dernier probleme : la transformation M = (m) doit ^ etre inversible. Il faut donc s'assurer que le determinant de la matrice jacobienne 15 de la transformation soit non nul en tout point de la maille de reference :
2 det J ] = det 6 4
e1 dX dx1 e2 dX dx1 e3 dX dx1
3 7 5 6= 0
Les equations (4.4) montrent que son calcul se ramene a des derivees des polyn^ omes de base Q(i) . (i) Le determinant est donc une fonction des coordonnees reelles des n uds Xj et des variables xi. Il peut arriver qu'il s'annule pour un certain ensemble de coordonnees nodales ou en certains points non nodaux 16. Tout bon logiciel se doit de faire cette veri cation avant de lancer un calcul 17. Lorsque le jacobien est non nul partout, on peut alors trouver ;1 en inversant le systeme (4.4). On trouve trois fonctions x e1, x e2 et x e3 x1 = x e1(X1 X2 X3)
e2(X1 X2 X3) x2 = x
x3 = x e3(X1 X2 X3)
(4.7)
Le calcul de la transformation inverse ;1 est di cile, sauf quand l'interpolation des coordonnees est strictement lineaire 18, car dans ce dernier cas, l'inversion de se ramene a des calculs algebriques matriciels. On verra plus loin que ce calcul n'est pas necessaire et qu'il su t de s'assurer que la transformation est inversible. Un grand nombre des elements proposes dans les logiciels ont une interpolation lineaire des coordonnees, mais on trouve aussi des interpolations de coordonnees de degre superieur, qu'on appelle souvent elements courbes 19. Dans le cas d'interpolation lineaire des coordonnees : etes sont approximees par des segments de droite joignant les { Les mailles lineiques et les ar^ n uds. S'il y a des n uds interieurs, c'est une ligne brisee . { Les mailles triangulaires et les faces triangulaires a trois n uds sont approximees par des triangles plans { Les mailles triangulaires et les faces triangulaires a plus de trois n uds sont approximees par des surfaces polyn^ omiales compliquees a contour polygonal gauche.
15. Les mecaniciens feront aisement le parallele avec la transformation d'un milieu continu, d'un domaine r de reference a un domaine t actuel. Dans ce cas, la matrice jacobienne n'est autre que la matrice des composantes du gradient de la transformation F dans un systeme de coordonnees cartesiennes, et on doit avoir detF 6= 0. 16. par exemple, si deux n uds d'un element reel sont confondus (ou trop proches), ou encore si les trois n uds d'un triangle reel sont alignes, etc. Il se peut aussi que le jacobien soit nul en certains points (xr yr zr ) non nodaux, par exemple si l'approximation geometrique de la maille (lineique, surfacique ou volumique) a des points doubles. 17. Ce n'est malheureusement pas toujours le cas! Certains logiciels se contentent de veri er seulement que la position des n uds reels ne conduit pas a une degenerescence de la geometrie de la maille (volume, surface ou longueur nuls), sans se preoccuper de savoir si l'approximation geometrique de la maille a des points doubles. Ce cas peut notamment se produire lorqu'on fait des calculs par pas : la reactualisation des n uds du maillage entre chaque pas peut rendre non inversible. 18. C'est a dire uniquement dans les elements lineiques a 2 n uds, les elements triangulaires a 3 n uds et les tetraedres a 4 n uds. Dans les autres, l'interpolation ne peut pas ^ etre strictement lineaire. Ce sont les seuls cas ou l'interpolation sur l'element reel reste polyn^ omiale. 19. Il convient de ne pas se laisser abuser par la representation graphique souvent simpli ee des elements: un element lineique a 3 n uds peut parfois ^ etre represente par un segment de droite avec les trois n uds alignes. Il faut conna^tre les possibilites du mailleur pour savoir si la gure a chee a l'ecran represente un element courbe ou si l'element est e ectivement rectiligne.
34
J. Garrigues
4.5. DERIVEES ET INTEGRALES DANS L'ELEMENT REEL

{ les mailles quadrangulaires et les faces quadrangulaires a quatre n uds sont approximees par des parabolo des hyperboliques. { les mailles quadrangulaires et les faces quadrangulaires a plus de quatre n uds sont approximees par des surfaces polyn^ omiales compliquees a contour polygonal gauche. On peut remarquer que dans le cas des elements lineiques a 2 n uds, des elements triangulaires a 3 n uds et des elements tetraedriques a 4 n uds, les approximations geometriques des lignes et des surfaces n'ont jamais de points doubles. Le jacobien ne peut s'annuler que pour des degenerescences geometriques (longueurs nulles, surfaces nulles ou volumes nuls). Pour les autres, il convient de veri er que le jacobien de n'est jamais nul non seulement aux n uds, mais aussi entre les n uds, c'est a dire que l'approximation geometrique ne contient pas de point double. Finalement , l'interpolation sur l'approximation geometrique de l'element reel est 20:
e(X1 X2 X3) = f
n X i=1
u(i) P (i) (x e1 (X1 X2 X3) x e2 (X1 X2 X3) x e3 (X1 X2 X3))
Elle est determinee quand on conna^t les coordonnees des n uds Xj(i) et les valeurs aux n uds u(i) . Quand on a choisi l'interpolation geometrique , on peut calculer a l'avance la matrice Q] des ei pour chaque element de coe cients des polyn^ omes de base des interpolations des coordonnees X reference. Par contre, le jacobien depend des coordonnees reelles, et sa non nullite ne peut ^ etre veri ee que quand on conna^t les coordonnees reelles des n uds. On peut toutefois preparer a l'avance la matrice Q0 ] des coe cients des derivees des polyn^ omes de base.
4.4.3 En resume ...

Un element est caracterise par : { la forme de sa maille de reference { la position des n uds dans la maille de reference (sommets obligatoires pour la conformite 21) er (c'est a dire la matrice P ] des coe cients des polyn^ { l'interpolation de reference f omes de base de l'interpolation) { l'interpolation des coordonnees des points non nodaux M 0 (c'est a dire la matrice Q] des coe cients des polyn^ omes de base de l'interpolation, et la matrice Q0 ] des coe cients des derivees de polyn^ omes de base pour le calcul de J ]). L'interpolation des coordonnees determine l'approximation geometrique de l'element. e sur l'element reel. { On en deduit une interpolation (en general non polyn^ omiale22) f
4.5 Derivees et integrales dans l'element reel

On verra plus loin que les derivees et les integrales dans l'element reel 23 des fonctions d'interpolation apparaissent dans la discretisation de la formulation variationnelle. On va montrer ici
20. Sauf dans le cas exceptionnel d'interpolations geometriques lineaires, les fonctions x ei (X1 X2 X3 ) (qui repree n'est donc pas polyn^ sentent 1 ) ne sont pas polyn^ omiales. L'interpolation f omiale en general. Cependant, elle appartient bien a l'espace de Sobolev H1 ( ). 21. mais ce n'est pas su sant, la conformite doit ^ etre veri ee. 22. sauf si est strictement lineaire. Voir note 18 23. En fait les integrales sont sur l'approximation geometrique de l'element
;
J. Garrigues
35

que les derivees de l'interpolation fe dans l'element reel s'expriment en fonction des derivees de l'interpolation fr dans l'element de reference, et qu'on peut ramener les integrales sur l'element reel a des integrales sur l'element de reference.
4.5.1 Mailles volumiques

Gradient gradfe = gradfer (grad );1 e et gradfer sont des tenseurs du premier ordre et grad un tenseur du second ordre 24. ou gradf
En termes de composantes dans un systeme de coordonnees cartesiennes : ; gradfe i = gradfer k grad ;1 ki er ;1 @ fe = @ f @Xi @xk J ] ki avec sommation sur l'indice k, et ou on a pose
er Puisque fe = f
;1, le gradient de l'interpolation f e est :
2 J] = 6 4
ou sous forme matricielle:
e1 @X @x1 e2 @X @x1 e3 @X @x1
e1 @X @x2 e2 @X @x2 e3 @X @x2
e1 @X @x3 e2 @X @x3 e3 @X @x3
3 7 5
e @f @X1
e @f @X2
e @f @X3
i h
e e e @f r @f r @f r = @x @x2 @x3 1
2 i6 4
e necessite donc l'inversion de J ] pour chaque point de element reel. Le calcul des derivees de f
Second gradient
Le second gradient de l'interpolation fe est : gradgradfe = grad gradfer (grad );1
e1 @X @x1 e2 @X @x1 e3 @X @x1
e1 @X @x2 e2 @X @x2 e3 @X @x2
e1 @X @x3 e2 @X @x3 e3 @X @x3
3;1 7 5
= (grad );T gradgradfer + gradfer grad (grad );1
En termes de composantes dans un systeme de coordonnees cartesiennes : er @fer @ ;J ;1 ki @ 2 f = ;J ;1 @2f ki @xk @xj + @xk @Xj @Xi @Xj avec sommation sur l'indice k. er peuvent ^ Les derivees secondes de f etre preparees d'avance en rangeant dans une matrice P 00 ] les coe cients des derivees secondes des polyn^ omes de base de l'interpolation fer . La matrice J ] 0 est calculee avec la matrice Q ], mais son inversion doit ^ etre faite pour chaque element reel car elle depend des coordonnees reelles des n uds.
24. Bien noter que grad 1 = (grad ) 1 . Il n'est donc pas necessaire de calculer explicitement les fonctions x ei (Mj ) pour calculer la matrice jacobienne de la transformation inverse.
; ;
36
J. Garrigues
Integrales
Si on note V 0 l'approximation geometrique de l'element reel et si on note Vr l'element de reference, l'element de volume de l'approximation geometrique de l'element reel est : dv = det J ] dx1dx2dx3 Une integrale sur l'element reel se ramene a une integrale sur l'element de reference par changement de variables :
G(X1 X2 X3)dv =
Vr
e1(x1 x2 x3) X e2(x1 x2 x3) X e3(x1 x2 x3)) det J ] dx1dx2dx3 G( X
ou det J ] est une fonction (en general non polyn^ omiale 25) de x1 x2 x3.
4.5.2 Mailles surfaciques

L'espace de reference est de dimension 2, On note S 0 l'approximation geometrique engendree par la transformation . La surface S 0 est une variete de dimension 2 plongee dans l'espace physique de dimension 3. L'interpolation fe n'est de nie que sur S 0, ses derivees par rapport aux 3 coordonnees Xi0 de M 0 n'ont aucun sens. Les derivees et les integrales sur les surfaces demandent un travail de preparation. La surface S 0 de point courant M 0 = (m) est naturellement parametree par x1 et x2:
ei sont la de nition parametrique de S 0. Les trois fonctions d'interpolation X La base naturelle tangente de S 0 pour ces parametres est : 0 @M 0 e e1 = @M 2= @x1 @x2 En general, cette base naturelle n'est ni orthogonale ni normee. Le plan de ni par e1 et e2 est le plan tangent a S 0 en M 0. Les composantes de la base naturelle de S 0 dans une base globale fE k g orthonormee sont : 2 e1]E = 6 4
e1 @X @x1 e2 @X @x1 e3 @X @x1
0 =X e1(x1 x2) X1
0 =X e3(x1 x2) X2
0 =X e3(x1 x2) X3
3 7 5
2 e2]E = 6 4
Une partie du travail peut ^ etre prepare a l'avance: la matrice Q0 ] des derivees des polyn^ omes de ei. Mais les vecteurs de la base naturelle sont fonction des coordonnees reelles et doivent base des X ^ etre calculees pour chaque element reel. On aura aussi besoin de la matrice 2 2 symetrique
1 e1 e1 e2 a ]= e e2 e1 e2 e2
e1 @X @x2 e2 @X @x2 e3 @X @x2
3 7 5
et de son inverse
25. voir note 18 page 34.
J. Garrigues
a ] = a ];1
37
Gradient de surface
Les composantes covariantes du gradient de surface de fe dans la base naturelle sont (voir sectionB.3 page 75) :
e(M ) = fer (m) = fe(x1 x2) le gradient de fe est encore : Puisque f

@ fr gradfe = @x Par changement de base (voir section B.2 page 75), on obtient les composantes de ce vecteur dans la base cartesienne orthonormee fEk g :
@f gradfe = @x
gradfe E = e1]E e2]E
]i
e e @f r @f r a ] @x @x2 1
La derivee dans une direction unitaire tangente u est donnee par :
e = gradfe u = gradfe u]E Du f E

Second gradient
]i
e dans la base naturelle Les composantes deux fois covariantes du second gradient surfacique de F sont rangees dans la matrice symetrique 2 2 (voir section B.3 page 75) :
grad gradfe
ou
]]
2f e @fe = @x@ @x ; @x ; er @ 2 fer ; @ f = @x @x @x ;
1a ; =2
@a + @a ; @a @x @x @x
Pour chaque element, il faut donc calculer les 8 ; . En fait on n'en calcule que 6 car ; = ; Remarque : Dans le cas d'elements triangulaires a trois n uds avec interpolation lineaire des coordonnees, la surface S 0 est un plan et les tous ; sont tous nuls. Dans tous les autres cas, ils ne sont pas tous nuls 26. Les composantes de grad gradfe dans la base fE k g se rangent dans la matrice 3 3 symetrique 27 (voir section B.2 page 75) :
i ] ] g rad g rad e
]]
0
f E=
e1]E e2]E
er @fer @2f a ] @x @x ; @x ;
"
T 1 ]E a ] e e2]T E
"
26. Dans beaucoup de logiciels (pour ne pas dire dans tous!), les approximations geometriques S des mailles reelles sont considerees comme quasi planes, m^ eme quand l'interpolation des coordonnees n'est pas lineaire. Cette hypothese simpli catrice permet d'eviter le calcul des ; , mais elle peut amener des erreurs importantes sur la valeur des derivees secondes si S est tres eloignee d'un plan. 27. on ne calcule donc que 6 termes.
0
38
J. Garrigues

Si on se donne deux directions tangentes unitaires u et v, la derivee seconde dans les directions u et v est : e = u grad gradfe v Duv f Si on calcule cette derivee avec les composantes dans fE k g : h ei u]E Duv fe = u]T E grad gradf E
]] ]]
Integrales
On montre en analyse des surfaces que l'element de surface de S 0 est : ds0 = (e1 ^ e2) n dx1dx2 On pose K (m) = (e1 ^ e2 ) n. Cette quantite doit ^ etre evaluee pour chaque element reel. On a donc : Z Z G(M 0) ds0 = G ( (m)) K (m) dx1dx2
S
0
Sr
4.5.3 Mailles lineiques

L'espace de reference est de dimension 1. On note C 0 l'approximation geometrique de la courbe e n'est de nie que sur C0 et les derivees reelle C engendree par la transformation . L'interpolation f 0 e de f par rapport aux 3 coordonnees Xi de M n'ont aucun sens. La courbe C 0 est naturellement parametree par x1 : 0 =X e1(x1) X20 = X e2(x1) X30 = X e3(x1) X1 ei sont la de nition parametrique de C0. Les trois fonctions d'interpolation X La base naturelle tangente de C 0 est le vecteur tangent non unitaire 0 dX e1 e2 e3 dX dX = E E e1 = dM 1+ 2+ dx1 dx1 dx1 dx1 E 3
Gradient
Le gradient lineique est dF e1 = dfr e1 gradfe = dx dx
1 1
La direction de derivation unitaire u est : 1 u = ke ek=r

1
1
e1 dX dx1
2
e2 2 e3 X dX + d dx1 + dx1
e1
La derivee le long de l'approximation geometrique de la courbe (c'est a dire la derivee par rapport a l'abscisse curviligne) est donc :
e df e e dl = Du f = gradf u = r dX e1
dx1
J. Garrigues
e df r dx1 e2 dX dx1
e3 X + ddx 1
39
second gradient
On deduit de ce qui precede :
0 e dB d2F Br 2 dl = dl @ dX e dx 0 d B B = dx @r 0 B =B B @r
=
1
1 1
+
2
e dF dx1 e2 dX dx1 e df r dx1 e2 dX dx1
e3 X + d dx1
1 C C A
2
e1 dX dx1
+
2e r 2 1
e3 X + d dx1
1 C dx1 C A
dl
2 2 2
e1 dX dx1
df dx e2 dX dx1
e1 d X e1 e2 d X e2 e3 d X e3 X X X +d + ddx er ddx d f dx1 dx2 1 dx2 1 dx2 1 1 1 ; dx 3 2 2 2 2 2 1 e3 2 e e d X e dX1 + dX2 + dX3 + dx1 dx1 dx1 dx1 e dX
2
1 C dx1 C C A dl
e1 dX dx1
e d2 f r dx2 1 e2 dX dx1
X dX2 d X2 dX3 d X3 1 2 + dx1 dx2 + dx1 dx2 dfer dx11 ddx 1 1 1 ; 2 2 2 2 2 2 e dx 1 X3 e2 e3 e1 dX dX dX + d + + dx1 dx dx dx
e e
2
Remarque : Si l'interpolation des coordonnees est lineaire (l'approximation geometrique de la 2X ei maille est une ligne brisee), alors les derivees secondes ddx 2 sont nulles. 1 Integrales
L'element de longueur est
v u u e1 !2 @X e2 !2 @X e3 !2 @ X t dl = + + dx1
dx1 dx1 dx1
On a donc :
Zl
l1
G(l) dl =
Z1
;1
v !2 !2 !2 u u e e e @ X @ X @ X 1 2 3 G (l(x1 )) t + + dx1
dx1 dx1 dx1
4.5.4 Remarque
On voit parfois proposer des elements de reference avec des continuites superieures a C0. Les fomules precedents donnant les gradients des interpolations reelles en fonction des gradients de l'interpolation de reference montrent que si les gradients sur l'element reel sont continus, ils ne le sont plus sur l'interpolation de reference (et inversement) car le jacobien de la transformation geometrique change d'un element a l'autre. La continuite C1 n'est donc qu'approximative sur l'element reel.
40
J. Garrigues
4.6. CONCLUSION
4.6 Conclusion
Dans les logiciels, on dispose d'une bibliotheque d'elements predetermines 28 . Les polyn^ omes de base des interpolations de reference sont calcules a l'avance et preprogrammes. L'utilisateur ordinaire n'a donc plus a s'en preoccuper. L'objectif de ce chapitre n'etait que de faire comprendre comment on les construit, et surtout de permettre a l'utilisateur d'un logiciel de les choisir en connaissance de cause. Il su t de retenir que choisir un element dans une bibliotheque c'est a la fois : { Choisir une forme de maille { Choisir la place des n uds dans la maille { Choisir une interpolation sur l'element de reference { Choisir une approximation geometrique de la forme de l'element reel Le choix des elements est delicat : c'est un compromis entre la qualite de la solution approchee et le co^ ut du calcul. L'exemple traite dans la section 2.2 page 9 montre clairement que la qualite de la solution approchee augmente avec le degre des polyn^ omes d'interpolation: Il est plus facile d'approcher une courbe avec des polyn^ omes de degre 2 qu'avec des segments de droite et plus le nombre de mailles est grand, plus la courbe est facile a approcher. Par contre, le nombre de n uds (et donc le nombre d'inconnues 29) augmente avec le degre d'interpolation et le nombre de mailles. La taille du systeme d'equations a resoudre (et donc le temps de resolution) peut devenir tres grande. Pour xer les idees, on peut retenir que pour un algorithme de resolution par une methode directe 30 , le temps de resolution est proportionnel au cube du nombre d'inconnues. Sans qu'on puisse le chi rer exactement, on peut dire que les elements de faible degre d'interpolation demandent un maillage plus n que les elements de degre plus eleve pour obtenir une ( (precision equivalente) ). Il faut donc se preoccuper de cette question des la phase de maillage, et l'examen d'une solution approchee peut parfois amener a recommencer un calcul avec un maillage et/ou des interpolations di erents. Une fois qu'on a choisi un maillage et des elements conformes, on a construit un espace F de fonctions F (M ) de nies sur tout le domaine , continue C0 , obtenue par juxtaposition des interpolations sur les elements reels 31 , a valeurs reelles et de carre integrable 32. Cet espace est donc inclus dans l'espace fonctionnel de Sobolev H1 ( ). S'il y a nddl inconnues par n ud, on a construit un espace de fonctions inclus dans H1 ( )nddl . Cette propriete est essentielle dans les theoremes de convergence de la methode des elements nis. Chaque fonction F de F est une interpolation sur tout le domaine , completement determinee par les Nddl = N nddl valeurs de ddl aux n uds. On peut donc ecrire : F (M ) =
N ddl X i=1
Les fonctions d'interpolation sur l'approximation geometrique, leurs derivees, et leurs integrales sur l'approximation geometrique sont alors determinees par les valeurs aux n uds.
u(i) F (i)(M )
(4.8)
ou u(i) est la valeur au ddl numero i (i 2 1 2 : : : Nddl ]).

28. et on a parfois la possibilite d'en de nir soi m^ eme. 29. le nombre d'inconnues est egal au nombre de n uds multiplie par le nombre d'inconnues par n ud. 30. Pour les systemes de grande taille, on prefere les methodes iteratives qui sont plus economiques dans ce cas. 31. ou plut^ ot, sur leur approximation geometrique. 32. Sous reserve que le jacobien de la transformation ne s'annule pas en un point.
J. Garrigues
41

Les fonctions F (i)(M ) sont les fonctions de base 33 de l'espace F . Ce sont les Nddl fonctions F (i)(M ) obtenues pour les Nddl -uplets de valeurs u(1) = i1 u(2) = i2 : : : u(j ) = ij : : : u(Nddl ) = iNddl . Les fonctions de base F (i)(M ) sont donc nulles partout sauf dans les elements qui contiennent le n ud i. On a en particulier : F (i) (M (j )) = ij La quasi totalite des logiciels d'elements nis utilise des elements de reference pour construire les interpolations sur l'approximation de l'element reel. Ce choix a des avantages informatiques evidents puisque les resultats d'une partie des calculs de l'interpolation sur l'element reel peuvent ^ etre faits a priori et sont codes dans le programme. Cependant, on a vu que cette methode presente quelques inconvenients : { la geometrie des mailles reelles est approximee { l'interpolation sur l'element reel n'est pas polyn^ omiale en general 34, ce qui complique les calculs de derivees et d'integrales sur l'element reel 35 . { la veri cation de l'inversibilite en tout point de la transformation entre l'element de reference et l'element reel est di cile et bien souvent ignoree. On pourrait construire directement des interpolations polyn^ omiales conformes 36 sur l'element reel, ce qui eliminerait les inconvenients precedents. Ce serait au prix d'une plus grande utilisation de la memoire 37, et le temps de preparation (pretraitement) est augmente 38. On est en droit de se demander, compte tenu du materiel actuel, si les contraintes informatiques qui justi aient l'utilisation d'elements de reference sont encore justi ees aujourd'hui...
Remarque :
33. en general non polyn^ omiales, voir note 18 page 34. 34. sauf pour les mailles lineiques a deux n uds, les mailles surfaciques triangulaires a trois n uds et les mailles volumiques tetrahedriques a quatre n uds. 35. Voir section 4.5 page 35 et suivantes. On verra plus loin que pour calculer les integrales sur les elements reels, on utilise des methodes approchees, le plus souvent la methode des points de Gauss. Ces integrales pourraient ^ etre calculees exactement si les interpolations sur l'element reel etaient polyn^ omiales. 36. comme on le fait habituellement sur l'element de reference 37. La description de chaque element doit contenir, en plus de sa geometrie, la description de ses propres fonctions de forme 38. Le mailleur doit non seulement de nir la geometrie des elements, mais aussi determiner les fonctions de forme de chaque element
42
J. Garrigues
Chapitre 5
Les formulations variationnelles

Pour resoudre un systeme di erentiel 1 , modelisant un systeme physique, il faut le mettre sous une forme integrale, appelee aussi forme variationnelle ou encore forme faible. Ces formulations peuvent ^ etre deduites du systeme di erentiel par des raisonnements mathematiques ou par des raisonnements physiques. On traitera a titre d'illustration des exemples en n de chapitre. Pour ^ etre su samment general, on utilise les notations suivantes : { On note u(M ) le champ inconnu 2. { On symbolise le systeme di erentiel par un operateur di ferentiel D devant ^ etre nul en tout point d'un domaine . etre nul 3 sur { On symbolise les conditions aux limites sur le bord par un operateur C devant ^ la frontiere @ . Par exemple, dans la section 2.2 page 9, on a : 2f (x) df (x) D : f (x) ;! d dx 2 + 2 dx ; 3 x " f (0) # C : f (x) ;! h df (x) i On verra d'autres exemples en n de chapitre. Tout probleme peut donc s'enoncer ainsi : Trouver le champ u(M ) de ni sur tel que { D(u(M )) = 0 8M 2 { C (u(N )) = 0 8N 2 @
dx x=1
5.1 Rappel de quelques resultats d'analyse fonctionnelle

Dans des espaces fonctionnels de fonctions de nies sur un domaine , avec des conditions sur les fonctions de cet espace et des conditions sur qu'on supposera satisfaites, on peut de nir un
1. aux derivees partielles si est une surface ou un volume. 2. Ici u est note en caracteres gras par souci de generalite. Le champ inconnu peut ^ etre scalaire vectoriel ou tensoriel. 3. l'operateur C peut ^ etre di erentiel ou non, selon que les conditions aux limites portent sur les derivees de u ou non.
J. Garrigues
43
CHAPITRE 5. LES FORMULATIONS VARIATIONNELLES

produit scalaire entre deux fonctions f et g : < f g >= ou
f (M ) g(M ) d
f (M ) g(M ) = f (M ) g(M ) si f et g sont a valeurs scalaires f (M ) g(M ) = f (M ) g(M ) si f et g sont a valeurs vectorielles f (M ) g(M ) = f (M ) g(M ) si f et g sont a valeurs tensorielles du second ordre
D'autre part, un produit scalaire a la propriete suivante: < f g > = 0 8g () f = 0 En utilisant cette propriete, il vient :
f (M ) g(M ) d = 0 8g(M ) () f (M ) = 0
(5.1)
5.2 Formes variationnelles

En utilisant (5.1), une nouvelle formulation du probleme est donc : { Trouver le champ u tel que
(M ) D(u(M )) = 0 8 (M )
{ avec C (u(N )) = 0 8N 2 @ Cette formulation est une formulation variationnelle ((equivalente)) 4 au systeme di erentiel initial. On obtient d'autres formulations variationnelles en tranformant les integrales. En e et, l'operateur D fait intervenir des operateurs tels que gradient, divergence, rotationnel, laplacien, et on donne dans la section suivante des formules pour transformer les integrales. Ces transformations d'integrales font appara^tre des integrales sur les frontieres de , qu'on appelle integrales de bord, dans lesquelles on peut tenir compte de tout ou partie des conditions aux limites: le nombre de conditions aux limites diminue. D'autre part, apres chaque transformation d'integrale, l'inconnue u appara^t avec un ordre de derivation diminue de 1. L'inter^ et de cette particularite appara^tra plus loin. Finalement, apres d'eventuelles transformations d'integrales, le probleme peut se mettre sous la forme variationnelle generale suivante : { Trouver le champ u tel que A( (M ) u(M )) = B( (M )) 8 (M ) (5.2) { avec C 0 (u(N )) = 0 8N 2 @ ou A et B sont des operateurs produisant des integrales sur et sur @ portant sur (M ) et u(M ) et leurs derivees, et ou C0 est un operateur produisant les conditions aux limites restantes.
0 0 0 0
4. dans un certain sens precise en analyse fonctionnelle: la fonction D(u(M )) et les fonctions doivent appartenir a certains espaces fonctionnels, eventuellement di erents. En fait, si 2 F , alors D(u(M )) 2 F ou F est l'espace dual de F . Si les fonctions de F sont tres regulieres, alors les fonctions de F le sont peu, et inversement. On montre en analyse fonctionnelle que si F est l'espace des fonctions de carre integrable sur , alors F = F .
44
J. Garrigues
5.3. TRANSFORMATIONS D'INTEGRALES
5.3 Transformations d'integrales

La validite des formules qui suivent est soumise a des conditions de regularite des champs (scalaires, vectoriels ou tensoriels) de nis sur que nous supposerons satisfaites.
5.3.1 Transformations d'integrales sur un volume

On note V le volume , S la frontiere du volume, et n la normale unitaire exterieurea S . Dans les formules qui suivent : { g est un champ scalaire quelconque de ni sur V , { V est un champ vectoriel quelconque de ni sur V , { T est un champ tensoriel du second ordre quelconque de ni sur V , { S est un champ tensoriel du troisieme ordre quelconque de ni sur V . La demonstration de ces formules est en annexe A page 71.
si u est un champ scalaire: gradu est un vecteur, u est un scalaire.
Z Z V grad u dv = ; u div V dv + u V n ds V Z S ZV Z
V
(5.3) (5.4)
g u dv = ;
grad g grad u dv + g grad u n ds

S
si u est un champ vectoriel: gradu est un tenseur du second ordre, div u est un scalaire, rot u est un vecteur, u est un vecteur.
Z Z T grad u dv = ; u div T dv + u T n ds V Z S ZV Z g div u dv = ; u grad g dv + g u n ds Z V ZV ZS V rot u dv = + u rot V dv ; (V ^ u) n ds V S ZV Z Z V u dv = ; grad V grad u dv + V grad u n ds

V V S
(5.5) (5.6) (5.7) (5.8)
si u est un champ tensoriel du second ordre : gradu est un tenseur du troisieme ordre, div u est un vecteur, u est un tenseur du second ordre.
V Z
S grad u dv = ; u div S dv + u S n ds SZ ZV V div u dv = ; u grad V dv + V u n ds V S Z ZV Z T u dv = ; grad T grad u dv + T grad u n ds

V V S
(5.9) (5.10) (5.11)
J. Garrigues
45

n n t
nt
Fig. 5.1 { Domaine
surfacique
5.3.2 Transformations d'integrales sur une surface

On note S la surface , C la courbe frontiere de la surface, n la normale unitaire a S et nt la normale exterieure a la courbe C et tangente a S . (voir gure 5.1). Dans les formules qui suivent : { g est un champ scalaire de ni sur S , { V est un champ vectoriel tridimensionnel de ni sur S , { V est la partie tangente de V , { T est un champ tensoriel du second ordre tridimensionnel de ni sur S , { T est la partie tangente de T { S est un champ tensoriel du troisieme ordre tridimensionnel de ni sur S , { S est la partie tangente de S f d g g Les operateurs di erentiels de surface grad, div, iv, r ot, e et e sont de nis en annexe B page 73. La demonstration des formules qui suivent est en section B.5 page 78.
si u est un champ scalaire: grad u est un vecteur tangent, e u est un scalaire. Z Z Z f V grad u ds = ; u d iv V ds + u V nt dl
si u
est un scalaire,
S
] ] ] ] est un champ vectoriel tangent: g rad u eu Z Z Z g ] T grad u ; u div T u T n Z Z Z fu ] ; u g rad u n Z Z Z ] ] ] rad u V g V eu ; g rad V g rad u n
S
]
S
g e u ds = ;
ZS
grad g grad u ds + u grad g nt dl

C
(5.12) (5.13)
est un vecteur tangent. ds =

S
f est un tenseur tangent du second ordre, d iv u

C
t dl
ds +
(5.14) (5.15)
t
g div ds =
ds =
g ds +
t dl
ds +
dl
(5.16)
J. Garrigues
46
5.4. EXEMPLE 1: LE PROBLEME THERMIQUE STATIONNAIRE
si u est un champ tensoriel tangent du second ordre : grad u est un tenseur tangent du g troisieme ordre, d iv u est un vecteur tangent, e u est un tenseur tangent du second ordre.
Z Z g V d iv u ds = ; grad V u ds + V (u nt) dl
S C
(5.17)
5.3.3 Transformations d'integrales sur une courbe

On note C la courbe , A et B sont les points frontiere de la courbe. Dans les formules qui suivent : { g est un champ scalaire de ni sur C , { V est un champ vectoriel tridimensionnel de ni sur C .
si u est un champ scalaire:
Z du Z dg g dl dl = ; u dl + u g]B A
C C
(5.18)
(C'est la formule classique d'integration par parties.)
si u est un champ vectoriel tridimensionnel:
Z dV du V dl dl = ; u dl dl + V u]B A C C
(5.19)
5.4 Exemple 1 : Le probleme thermique stationnaire

On suppose que le domaine est un volume V de frontiere S et de normale exterieure n. Le champ inconnu T est un scalaire. Le probleme s'enonce ainsi: { Trouver le champ de temperatures T tel que : T (M ) = f (M ) 8M 2 V ou f (M ) est une source de chaleur volumique connue. { avec les conditions aux limites: T (N ) = Tf (N ) 8N 2 ST gradT (N ) n(N ) = f (N ) 8N 2 Sf (On impose un champ de temperature stationnaire sur la partie de frontiere ST et un champ de ux thermique stationnaire sur la partie de frontiere Sf ) En utilisant (5.1) page 44, on obtient une premiere forme variationnelle: { Trouver le champ de temperatures T tel que :
( T ; f ) dv = 0 8 (M ) de ni sur V
{ avec les m^ emes conditions aux limites.

J. Garrigues
47

En utilisant la formule (5.4) page 45 , on obtient une seconde forme variationnelle: { Trouver le champ de temperatures T tel que :
; grad
V
gradT dv +
gradT n ds ;
f dv = 0 8
{ avec les m^ emes conditions aux limites. Or, l'integrale de bord peut se decomposer :
gradT n ds =
ST
gradT n ds +
Z
Sf
gradT n ds
En tenant compte de la condition aux limites de ux sur Sf , il vient:
gradT n ds =
ST
gradT n ds +
Sf
ds
Finalement, la seconde forme variationnelle du probleme est : { Trouver le champ de temperatures T tel que :
; grad
gradT dv +
A( T )
{z
ST
gradT n ds = ;
Z
Sf
} |
f ds +
Z
V
B( )
{z
f dv 8
{ avec la seule condition aux limites: T (N ) = Tf (N ) 8N 2 ST Le probleme thermique stationnaire se met bien sous la forme donnee en (5.2) page 44. On peut obtenir une troisieme forme variationnelle en utilisant la formule (5.3) page 45 : { Trouver le champ de temperatures T tel que :
dv ;
T grad
n ds +
Sf
f
ds +
ST
gradT n ds =
V
f dv = 0 8
{ avec la seule condition aux limites: T (N ) = Tf (N ) 8N 2 ST En tenant compte de la condition aux limites de temperature sur ST , il vient :
T grad
n ds =
=
ZSf
Sf
T grad T grad
n ds + n ds +
ZST
ST
T grad Tf grad
n ds n ds
J. Garrigues
48
5.5. EXEMPLE 2: STATIQUE DU SOLIDE DEFORMABLE EN PETITES DEFORMATIONS. La troisieme formulation variationnelle s'ecrit donc encore : { Trouver le champ de temperatures T tel que : Z Z Z gradT n ds = T dv ; T grad n ds +
|V
Sf
Z | ST
Tf grad
A( TZ ) n ds ;
{z
ST
f
B( )
{zSf
ds +
Z
V
f dv 8
(5.20)
{ sans conditions supplementaires. Bien que pouvant sembler plaisante (il n'y a plus de conditions a respecter sur les frontieres), cette derniere formulation presente des inconvenients qui appara^tront dans le prochain chapitre.
5.5 Exemple 2 : Statique du solide deformable en petites deformations.

Le probleme est : { Trouver le champ de deplacements tel que : div | {z+ f} = 0 D( ) 1 grad + gradT 5 , si le solide deformable est elastique ou = 2 + Tr G, et = 2 isotope lineaire en petites deformations. Cette equation est une equation vectorielle. { avec les conditions aux limites C ( ) : { des deplacements imposes 6 0(N ) sur la partie de frontiere Sd : (N ) ; 0(N ) = 0 8N 2 Sd { des contraintes imposees 7 q0 (N ) sur la partie de frontiere Sf : (N ) n(N ) ; q0 (N ) = 0 8N 2 Sf
5.5.1 Premiere forme variationnelle

En utilisant (5.1) page 44, on obtient une premiere formulation equivalente : { Trouver le champ de deplacements tel que :
5. On ne remplace pas et par leurs valeurs en fonction de pour alleger l'ecriture. Si on le faisait, on obtiendrait l'equation de Navier. On peut tres bien obtenir les m^ emes formes variationnelles a partir des equations de Navier. L'inconvenient serait que les equations de Navier ne sont valables que pour les solides elastiques lineaires isotropes en petites deformations alors que la formulation variationnelle qu'on va ecrire ici est valable pour toutes les lois de comportement : il su t de changer l'expression de en fonction d'un tenseur de deformation et de preciser la de nition du tenseur de deformation qu'on utilise. 6. eventuellement nuls 7. eventuellement nulles
J. Garrigues
(div + f ) dv = 0 8 (M ) champ vectoriel dans V .
(5.21)
49

T . ou = 2 + Tr G, et = 1 2 grad + grad Cette equation est une equation scalaire. { avec les conditions aux limites: (N ) ; 0 (N ) = 0 8N 2 Sd (N ) n(N ) ; q0 (N ) = 0 8N 2 Sf
5.5.2 Seconde forme variationnelle

En utilisant la formule (5.10) page 45, on obtient le probleme equivalent suivant : { Trouver le champ de deplacements tel que :
grad dv +
n ds +
f dv = 0 8 (M )
{ avec les conditions aux limites:
(N ) ; 0 (N ) = 0 8N 2 Sd (N ) n(N ) ; q0 (N ) = 0 8N 2 Sf
L'integrale de bord peut ^ etre decomposee en deux integrales :
On peut alors eliminer la condition aux limites de contrainte imposee en la reportant dans l'integrale sur Sf :
n ds =
Sd
n ds + n ds =
Sf
n ds
Sf
Sf
q0 ds
Le probleme devient : { Trouver le champ de deplacements tel que :
grad dv +
f dv +
Sf
q0 ds +
Z
Sd
n ds = 0 8 (M )
{ avec la condition aux limites:
(N ) ; 0 (N ) = 0 8N 2 Sd
Le tenseur des contraintes etant symetrique, on peut ecrire : grad = Sym(grad ) Le probleme devient : { Trouver le champ de deplacements tel que :
Sym(grad ) dv +
f dv +
Z
Sf
q0 ds +
Z
Sd
n ds = 0 8 (M )
(5.22)
{ avec la condition aux limites:
(N ) ; 0 (N ) = 0 8N 2 Sd
J. Garrigues
50
5.5. EXEMPLE 2: STATIQUE DU SOLIDE DEFORMABLE EN PETITES DEFORMATIONS.
Interpretation mecanique :
Appelons ((champ de vitesses virtuelles )) 8 le champ vectoriel arbitraire . Si est une vitesse virtuelle, Sym(grad ) s'interprete comme un tenseur des vitesses de deformations virtuel. R L'integrale V ; Sym(grad ) dv s'interprete alors comme une puissance virtuelle des e orts interieurs dans le champ de vitesses virtuelles . R L'integrale Sf q0 ds s'interprete comme une puissance virtuelle des e orts de surface dans le champ de vitesses virtuelles . R L'integrale Sd n ds s'interprete comme une puissance virtuelle des e orts de liaison dans le champ de vitesses virtuelles . R f dv s'interprete comme une puissance virtuelle des e orts de volume dans le L'integrale V champ de vitesses virtuelles . L'equation (5.22) exprime donc que la puissance virtuelle des e orts interieurs et exterieurs est nulle dans tout champ de vitesses virtuelles applique a la solution, ce qui est l'expression du theoreme des puissances virtuelles 9 . On pouvait donc deduire cette formulation par des raisonnements mecaniques. On reecrit le probleme (5.22) en isolant les termes contenant le champ inconnu : { Trouver le champ de deplacements tel que :
grad dv +
A(
Z
)
{z
Sd
n ds = ; } |
f dv ; {z
T ou = 2 + Tr G, et = 1 2 grad + grad { avec la condition aux limites: 0 (N ) = 0 8N 2 Sd | (N ) ; {z } 0 C( )
B( )
Sf
q0 ds 8 (M ) (5.23) }
Le probleme du solide elastique se met bien sous la forme donnee en (5.2) page 44.
Remarques :
On peut noter que la loi de comportement n'a pas ete explicitement utilisee dans les transformations d'integrales. Cette formulation variationnelle est en fait valable pour toute loi de comportement. Si la loi de comportement elastique lineaire isotrope en petites deformations, l'integrale R on admet grad dv peut encore ^ etre transformee. On laisse le soin au lecteur de trouver une troisieme V forme variationnelle du probleme 10, dans laquelle on pourra faire dispara^tre les conditions aux limites de deplacement impose sur le bord.
8. On peut aussi bien l'appeler deplacement virtuel. Dans ce cas, les puissances virtuelles s'appellent travaux virtuels. A vrai dire, est un champ arbitraire qu'on peut baptiser comme on veut. Ces denominations sont censees aider a xer les idees. Elles peuvent toutefois induire en erreur: on peut montrer en mecanique que quand on remplace les ((vitesses virtuelles)) par les vitesses reelles , les di erentes integrales donnent bien des puissances reelles, alors qu'il n'en est rien si on remplace les ((deplacements virtuels)) par les deplacements reels! On n'obtient que des integrales homogenes a des travaux, mais qui ne sont pas les travaux reels. 9. Si est appele ((deplacement virtuel)), le theoreme s'appelle theoreme des travaux virtuels. 10. C'est la forme variationnelle qu'on obtiendrait directement a partir des equations de Navier.
J. Garrigues
51
5.6 Exemple 3 : Statique des coques

Le domaine est une surface S de normale unitaire n ayant pour frontiere une courbe C dont le triedre de Darboux-Ribeaucourt est fnt t ng (voir gure 5.1 page 46). Le probleme s'enonce ainsi: { Trouver le champ de deplacements (M ), M 2 S satisfaisant en tout point M de S aux equations di erentielles11 :
f g f ( ) = 0 (2 equations d'equilibre tangent) p+d ivT ( ) ; B ( ) d ivM fd g f ( ) = 0 (1 equation d'equilibre normal) ivM p3 + T ( ) B ( ) + div
{ p = p + p3 n est le chargement surfacique { T ( ) est le tenseur tension (tenseur tangent). Cette fonction de est la loi de comportement membranaire de la coque 12. f ( ) = E M ( ), ou M ( ) est le tenseur des moments lineiques (tenseur tangent). {M Cette fonction de est la loi de comportement de exion de la coque 13. { B ( ) est le tenseur courbure de la surface deformee 14 { avec les conditions aux limites: { des deplacements imposes 0 (N ) sur une partie de frontiere Cd : (N ) ; 0(N ) = 0 8N 2 Cd { des rotations r0 (N ) (autour de t) imposees sur une partie de frontiere Cr :
ou
t grad nt ; r0(N ) = 0 8N 2 Cr
{ des tensions tangentes T 0 (N ) imposees sur une partie de frontiere CT :
T ( (N )) nt ; t0(N ) = 0 8N 2 CT
{ des e orts tranchants q0 (N ) (suivant n) imposes sur une partie de frontiere Cq :
g f ( (N )) nt ; q0(N ) = 0 8N 2 Cq d ivM
{ des moments lineiques tangents m0 imposes sur une partie de frontiere Cm :
f ( (N )) nt ; m0 = 0 8N 2 Cm M ( (N )) nt ; m0 = ;E M f ( (N )) nt ; E m0 = 0 8N 2 Cm () M
11. Voir un cours de theorie des coques. 12. En petites deformations, T ( ) s'exprime en fonction du tenseur des petites deformations de surface ( ). En deformations nies, les lois de comportement sont plus compliquees. f ( )) s'exprime en fonction du tenseur de variations de courbures K ( ). 13. En petites deformations M ( ) (ou M En deformations nies, les lois de comportement sont plus compliquees. 14. En petites deformations de coques on a B = B 0 + K ( ), ou B 0 est la courbure initiale et K ( ) est le tenseur de variation de courbure (tenseur tangent) linearise en . En deformations nies, l'expression de la courbure actuelle B en fonction de B 0 et est une expression non lineaire plus compliquee.
52
J. Garrigues
5.6. EXEMPLE 3 : STATIQUE DES COQUES
5.6.1 Premiere forme variationnelle

Trouver le champ de deplacements (M ) tel que :
SZ S
f T ( ) ; B( ) d g f ( ) ds = 0 8 (M ) champ vectoriel tangent p + div ivM f g f ( ) ds = 0 8 (M ) champ scalaire p3 + T ( ) B ( ) + d ivd ivM
(5.24) (5.25)
avec les m^ emes conditions aux limites que pour la forme di erentielle du probleme.
5.6.2 Seconde forme variationnelle

En utilisant (5.17) page 47, on transforme les integrales de (5.24):
Z f divT ds = ; grad ZS
=
S S
Z
S
g f ds = B d ivM
=
] ] ; g rad Z ] g rad Z ] g rad Z

S C;Cm
T nt dl Z T ds + t0 dl + T nt dl CT Z C;CT f ds ; f nt dl ( B) M B M C Z f ds ; ( B) M B E m0 dl Cm f nt dl B M
T ds +
ZC
En utilisant (5.15) page 47, on transforme les integrales de (5.25):
Z Z f g f g f g f nt dl divdivM ds = ; divM grad ds + d ivM S C Z Z Z g f grad ds + g f nt dl =; d ivM q0 dl + d ivM

S Cq C;Cq
] ]
;
On peut noter qu'avec ces transformations d'integrales, les conditions aux limites en e ort (tensions, moments lineiques et e orts tranchants) peuvent ^ etre prises en compte. Finalement, (5.24) devient :
; grad
S
T ds +
C;CT
Z T nt dl + grad ( S Z
Cm
f ds + B) M
B E m0 dl = ;
C;C Z m S
f nt dl B M Z p; t0 dl
CT
et (5.25) devient :
Z Z Z Z g f grad ds + g f nt dl = ; p3 ds ; q0 dl (5.27) ivM d ivM T B; d S C;Cq S Cq S ou B = B ( ) est la courbure apres deformation. Elle est fonction (en general non lineaire) du champ de deplacements (et aussi de la courbure initiale) f=M f ( ) ) sont les lois de comportement de la coque qui et T = T ( ) et M = M ( ) (ou M varient suivant les theories utilisees. R On peut encore noter que l'integrale de bord Cm B E m0 dl reste une fonction de a Z
cause de la variation de courbure 15 .
S
(5.26)
15. Certains auteurs negligent la variation de courbure dans ce terme (B ' B 0 ) pour pouvoir le faire passer dans R ; R f ds est parfois simpli e en f ds le second membre. De m^ eme, le terme grad B M grad B M
S
J. Garrigues
53
5.6.3 Une variante de la formulation precedente

En utilisant (5.17) page 47, on peut encore transformer les integrales de (5.27) page 53:
g f grad ds = ; grad grad d ivM

=
S S
L'equation (5.27) page 53 devient :
T B+
Z f f n ] ] ] M g rad M Z Z f ] ] ] ; grad grad M g rad E m Z f n ] g rad M Z Z Z f ; f n g f n ] ] ] g rad g rad M g rad M d ivM Z Z Z ] ; ; g rad E m Z
ds + ds +
C
t dl
0
Cm
dl
C;Cm
t dl
ds
C;Cm S
t dl +
C;Cq
t dl =
0
p3 ds
Cq
q0 dl +
Cm
dl
5.6.4 Remarques
L'interpretation de ces integrales en termes de ((puissances virtuelles)) (ou de ((travaux virtuels))) est beaucoup plus delicate que pour les problemes de solides tridimensionnels en raison de la complexite geometrique due a la variation de la courbure des coques pendant la deformations. Les auteurs qui tentent d'etablir des formulations variationnelles a partir d'interpretations mecaniques des champs arbitraires et sont amenes a faire des hypotheses simpli catrices plus ou moins justi ees pour garder des interpretations simples. D'autre part, bien que la deformation et les changement de courbure d'un milieu continu bidimensionnel soit completement determinee par le champ de deplacement, la plupart des elements nis utilises dans les logiciels utilisent des degres de liberte ((de rotation)) 16. Ces degres de liberte de rotation posent deux problemes : { Dans les problemes modelisant a la fois des coques et des milieux massifs, on a du mal a ecrire les conditions de raccord entre les deux types de milieux continus { Les seules rotations engendrant une ((puissance virtuelle)) sont les rotations tangentes a la coques. Que faire de la partie normale de la rotation? En n, les formulations variationnellles donnees ici sont valables pour toutes les lois de comportement (T = T ( ) et M = M ( )). Point de vue personnel : Plut^ ot que de chercher a ((faire de la mecanique avec les elements nis)) 17, il me semble plus sain d'etablir un theorie de coques 18 avec une loi de comportement etablie sur des hypotheses physiques precises 19 , et de ne passer qu'ensuite a des formulations variationnelles du systeme di erentiel pour la resolution numerique du probleme.
5.7 Formulations a plusieurs champs

Bien qu'introduisant un plus grand nombre de d'inconnues par n uds, les formulations qui suivent presentent des avantages quand les grandeurs pertinentes pour juger des resultats sont des gran16. Ils sont introduits par des interpretations mecaniques de ((puissances virtuelles de exion ou de torsion)) 17. C'est a dire considerer les elements nis comme de petits bouts de matiere avec une loi de comportement. 18. Les equations di erentielles d'equilibre des coques donnees ici sont exemptes de toute hypothese 19. petites ou grandes deformations, hypotheses cinematiques diverses sur les normales materielles, developpements asymptotiques, etc
54
J. Garrigues
5.7. FORMULATIONS A PLUSIEURS CHAMPS

deurs derivees de l'inconnue principale. En e et, on verra en 7.1.2 page 68 que la continuite interelements des grandeurs derivees n'est pas assuree 20 , ce qui rend malcommode leur representation et leur exploitation. Il peut ^ etre judicieux de considerer les grandeurs auxquelles on s'interesse comme des inconnues (aux n uds) et non comme des grandeurs derivees de la solution. Exemple : Formulation a deux champs d'un probleme de statique des solides : Le probleme de statique des solides peut ^ etre vu comme la solution du systeme di erentiel suivant:
div + f = 0 ; G( ) = 0
ou la loi de comportement G est un operateur di erentiel. Le probleme est complete par des conditions aux limites. Ce qui conduit a la formulation:
ZV
V
(div + f ) dv = 0 8 champ vectoriel ( ; G ( )) dv = 0 8 champ tensoriel symetrique
(5.28) (5.29)
On en deduit d'autres par transformations d'integrales : l'equation (5.28) se transforme comme precedemment (voir (5.23) page 51)
grad dv +
A(
Z
)
{z
Sd
n ds = ; } |
f dv ; {z
B( )
Sf
q0 ds 8 (M ) champ vectoriel }
La transformation de l'equation (5.29) depend de la loi de comportement G . Par exemple, en elasticite lineaire isotrope petites deformations,
G ( ) = 2 Sym(grad ) + div G
et
Z
V
G ( ) dv =
=
Z ZV Z
VZ
(2 Sym(grad ) + div G) dv 2
V
grad dv +
2
=;
div dv + 2 gradTr dv +
ZV
Tr div dv
;
l'equation (5.29) devient :
ZS
S
n ds
Tr
n ds
dv +
Z
V
div dv +
S
Tr
gradTr dv ; 2 n ds V S n ds = 0 8 champ tensoriel symetrique 55
20. Parfois, on n'evalue ces grandeurs derivees qu'en certains points de l'element (les points de Gauss)
J. Garrigues

dans laquelle on peut faire intervenir des conditions aux limites de deplacement dans les integrales de bord. Dans cette formulation, les inconnues sont (6 ddl) et (3 ddl), mais on obtiendra des valeurs de contraintes aux n uds. L'equation d'equilibre et la loi de comportement ont une formulation ( (faible) ), c'est a dire qu'elles seront ne seront respectees que globalement par la solution approchee 21 . On laisse le soin au lecteur d'ecrire une formulation a trois champs en prenant comme inconnues , et 22. Le choix entre ces deux formulations resulte d'un compromis: si on utilise des elements conformes, { soit on a 9 inconnues par n ud et une solution approchee continue en deplacement et en contraintes, { soit on a seulement trois inconnues par n ud avec une solution approchee continue en deplacements mais discontinue en contraintes 23 . On peut employer le m^ eme procede pour obtenir des contraintes continues dans les problemes de coques. Il faut constater que dans les logiciels industriels courants actuels, on a generalement choisi des formulations avec le plus faible nombre d'inconnues possibles.
5.8 Conclusion
On peut donner beaucoup de formes variationnelles a un probleme de physique donne sous forme di erentielle : { en choisissant les inconnues principales du probleme { en utilisant plus ou moins de transformations d'integrales Le choix de l'une d'entre elle pour obtenir une approximation de la solution par elements nis est guide par : { le nombre de champs inconnus (co^ ut du calcul) { l'ordre de derivation des inconnues et des champs arbitraires (On verra dans le chapitre suivant que cet ordre de derivation a des consequences sur le choix des fonctions d'interpolation). Dans les logiciels d'elements nis orientes ((metier)), ce choix est rarement laisse a l'utilisateur. Ce sont souvent des ((bo^tes noires)) dans lesquelles la formulation variationnelle est codee.
21. alors que dans la formulation a un champ, seule l'equation d'equilibre avait une formulation faible. 22. qui presente peu d'inter^ et si la relation entre et n'est pas di erentielle (elle conserve donc la continuite interelement). 23. La taille du probleme (et sa place en memoire) est divisee par 3 , et le temps de calcul est environ divise par 33 = 27. La continuite des contraintes co^ ute cher en ressources informatiques, mais le materiel evolue...
56
J. Garrigues
Chapitre 6
La discretisation
Que le probleme a resoudre soit formule sous forme di erentielle ou sous forme integrale, sa solution exacte est toujours inaccessible. L'objet de ce chapitre est de montrer comment choisir une solution approchee dans la famille F H1 ( )nddl de fonctions construite par les operations de maillage et de choix d'elements (voir (4.8) page 41). On rappelle que chaque fonction F (M ) de F peut s'ecrire :
N ddl X j =1
F (M ) =
u(j ) F (j ) (M )
et que les fonctions de base F (j )(M ) de l'espace F sont des fonctions 1 nulles partout, sauf dans les elements qui contiennent le n ud j .
6.1 Principe de la discretisation

On a vu dans le chapitre precedent que la formulation variationnelle exacte d'un probleme peut se mettre sous la forme: { Trouver u(M ) tel que
A( u) = B( ) 8 (M ) de ni dans C 0 (u) = 0
{ avec les conditions (eventuellement vides) Si on cherche une solution appartenant a l'espace F , on obtient le probleme approche suivant: { Trouver F (M ) tel que
A( F ) = B( ) 8 (M ) de ni dans C 0 (F ) = 0
(6.1)
{ avec les conditions aux limites:

1. On rappelle que ces fonctions de forme ne sont pas polyn^ omiales en general, sauf si l'interpolation geometrique est strictement lineaire (Voir note 18 page 34).
J. Garrigues
57
CHAPITRE 6. LA DISCRETISATION
Ce probleme n'a en general pas de solution 2. On va donc remplacer la condition ((8 (M ))) par une condition moins contraignante : On ne cherche a satisfaire l'equation (6.1) que pour certains (M ). D'autre part, il faut respecter les conditions aux limites C 0 (F ) = 0 si elles existent 3 . On restreint donc l'espace F aux seules fonctions qui respectent ces conditions. On appelle Fad l'espace des fonctions admissibles, c'est a dire l'espace des fonctions qui satisfont aux conditions aux limites C 0 (F ) = 0. Fad est un sous espace de F . En e et, imposer les conditions C 0 (F ) = 0 revient a imposer NC relations sur les valeurs des u(i). La dimension de l'espace des fonctions admissibles est donc Ninc = Nddl ; NC . Pour determiner une fonction F (M ) 2 Fad , il faut determiner un Ninc-uplet de valeurs :
u(1) u(2) : : : u(i) : : : u(Ninc )
Si on a choisi correctement Ninc fonctions (i) (M ), pour chacune d'elles, l'equation (6.1) va fournir une equation scalaire :
A(
(i)
F ) ; B(
(i)
) = 0 ou F =
N ddl X j =1
u(j ) F (j )
(6.2)
La construction de ce systeme d'equations s'appelle l'assemblage. Les fonctions (i) choisies sont souvent appelees fonctions test. Pour trouver une solution approchee, on a donc un systeme de Ninc equations a Ninc inconnues u(1) u(2) : : : u(i) : : : u(Ninc ) a resoudre. { Si A( u) est une application lineaire en u, les termes A( (i) F ) conduisent a un systeme lineaire en u(i). Cette classe de problemes est appelee problemes lineaires 4 . { Si A( u) n'est pas lineaire en u, la solution du systeme (6.2) est plus complexe 5. Il existe cependant des methodes numeriques pour le resoudre qui se ramenent a une successsion de problemes lineaires. Tout le probleme reside donc dans le choix judicieux des fonctions (i)(M ) : { Une condition imperative est que le choix des (i) (M ) doit ^ etre tel que le systeme d'equations (6.2) soit regulier. En particulier, le choix des (i) (M ) et le choix de l'espace des fonctions d'interpolation F (M ) ne sont pas independants du choix de la formulation variationnelle: En e et, a chaque transformation d'integrale sur , on diminue l'ordre de derivation de u et on augmente l'ordre de derivation de . Il faut donc que les fonctions F et aient des derivees d'un ordre su sant non nulles partout pour que le systeme soit regulier 6. { On dispose d'algorithmes e caces pour la resolution de systemes lineaires symetriques. On preferera donc un choix de (i) (M ) qui conduit a un systeme symetrique 7.
2. sauf si par chance la solution exacte appartenait a F . 3. En e et, on a vu dans le chapitre precedent que suivant le nombre de transformations d'integrales qu'on e ectue, le nombre de conditions a imposer a la solution diminue, et peut m^ eme ^ etre reduit a 0. On peut aussi se referer a l'exemple 2.2 page 9. 4. Le lecteur est invite a veri er que les problemes de mecanique traites en exemple dans le chapitre precedent qu'une condition necessaire (mais pas su sante) est que la loi de comportement soit lineaire. 5. On peut noter que l'unicite et la convergence de la solution par elements nis n'est classiquement etablie que pour les problemes lineaires. Dans presque tous les autres cas, la methode des elements nis reste heuristique. 6. Par exemple, dans la formulation (5.20) page 49, la fonction appara^t sous la forme de son laplacien. Il faut choisir des fonctions a laplacien non nul partout. Une interpolation lineaire ne conviendrait pas. 7. Le choix fait dans l'exemple 2.2 page 9 n'est donc pas judicieux de ce point de vue.
58
J. Garrigues
6.2. LA METHODE DE GALERKINE

Il existe plusieurs methodes pour choisir les fonctions (i) (M ). Elles peuvent conduire a des systemes symetriques ou non, et la precision de la solution est plus ou moins sensible au choix des . On ne developpe ici que la methode la plus utilisee dans les logiciels d'elements nis : la methode de Galerkine.
6.2 La methode de Galerkine

On choisit comme fonctions (i) les fonctions de base de l'espace Fad . Ce choix a les consequences suivantes : { Le nombre de fonctions (i) (et donc le nombre d'equations scalaires) est automatiquement egal au nombre d'inconnues Ninc. { On montre en analyse numerique que le systeme d'equations (6.2) page 58 obtenu pour ce choix de (i) est necessairement regulier pour une classe de problemes appeles problemes elliptiques : dans ce cas, A( u) et B( ) doivent ^ etre lineaires en u et en ,
A(
(i)
(i)
F ) = A(F
(i)
N ddl X j =1
u(j ) F (j ))
Pour les autres problemes, soit le systeme d'equations n'est pas lineaire en u(j ), soit il est lineaire mais non garanti regulier, et la methode de Galerkine peut se reveler inadequate. { Le fait de choisir des fonctions (i) appartenant a Fad peut simpli er certaines integrales de bord. { La methode de Galerkine impose une contrainte sur le choix de la forme variationnelle: les fonctions (i) (M ) et les fonctions d'interpolation F (M ) appartenant au m^ eme espace de fonctions, il faut choisir une forme variationnelle dans laquelle les derivees des fonctions et u soient d'un ordre tel qu'elles ne s'annulent pas partout. La recherche de la solution approchee se ramene donc a la resolution du systeme d'equations :
PNddl u(j) F (j)) doit ^ et A(F etre une forme quadratique de nie positive en F (j ) . j =1
0 Nddl 1 X A @F (i) u(j ) F (j )A ; B F (i) = 0

j =1
Les u(j ) trouves determinent la solution approchee F . Ce systeme est lineaire et regulier pour les problemes elliptiques.
6.3 Application de la methode de Galerkine au probleme thermique stationnaire

On a trouve trois formes variationnelles pour le probleme thermique stationnaire.
6.3.1 forme 1
{ Trouver le champ de temperatures T tel que :
|V
J. Garrigues
( T ; f ) dv = 0 8 (M ) de ni sur V
A( T )
{z
59
{ avec les conditions aux limites: T (N ) = Tf (N ) 8N 2 ST gradT (N ) n(N ) = f (N ) 8N 2 Sf La fonction inconnue T (M ) appara^t par son laplacien. Il faut donc choisir un espace d'interpolation F tel que le laplacien de ses fonctions de base soit non nul, ce qui interdit les interpolations lineaires. D'autre part, la construction de l'espace Fad respectant la condition de ux est di cile: Il faut choisir des elements qui contr^ olent la derivee normale au bord, ce qui est di cile 8. Par contre, l'espace des fonctions test (qui est aussi Fad ) para^t trop luxueux, car n'est pas derive dans cette formulation. La forme A( T ) est bien lineaire en et T . Le systeme en T (i) est donc lineaire et regulier. Bien que la methode de Galerkine appliquee a cette formulation paraisse ((desequilibree)) du point de vue des derivees de et T , elle est pratiquable si on dispose d'elements adequats.
6.3.2 forme 2
{ Trouver le champ de temperatures T tel que :
; grad
gradT dv +
A( T )
{z
ST
gradT n ds = ;
Z
Sf
} |
ds +
Z
V
B( )
{z
f dv 8
{ avec la seule condition aux limites: T (N ) = Tf (N ) 8N 2 ST Les fonctions et T n'apparaissent que sous la forme de leur gradient. On peut donc choisir un espace d'interpolation F plus simple: il su t que les fonctions de base de F aient un gradient non nul. On peut donc prendre des elements plus simples: des interpolations lineaires 9 sont su santes. L'espace Fad est plus facile a construire 10 : il su t d'imposer des valeurs aux T (j ) sur la frontiere ST . D'autre part, puisque 2 Fad , on a
ST
gradT n ds =
ST
Tf gradT n ds
Cette integrale est plus facile a calculer. La forme A( T ) est toujours lineaire en et T . Le systeme algebrique
; gradF (i)
V
N inc X j =1
T (j ) gradF (j ) dv +
ST
Tf
N inc X
Zj=1
Sf
T (j ) gradF (j ) n ds =
F (i) f ds +
F (i) f dv
8. On verra plus loin (Remarque page 63, qu'en fait il n'est pas absolument necessaire de construire Fad . Cette objection peut ^ etre contournee. 9. Il n'est pas interdit de prendre des interpolations de degre supereur pour ameliorer la qualite de la solution approchee. 10. voir note 8
60
J. Garrigues
6.4. APPLICATION DE LA METHODE DE GALERKINE AU PROBLEME DE STATIQUE DE SOLIDE DEFORMABLE est lineaire en T (i). Cette formulation variationnelle est bien adaptee a la methode de Galerkine. Les integrales sur V se ramenent a une somme d'integrales sur les mailles. Par exemple:
F (i) f dv =
Nm Z X
m=1
F (i) f dv
ou Nm est le nombre de mailles. Il est a noter que les integrales de volume et de surface de l'equation i ne portent nalement que sur les fonctions de base de l'espace d'interpolation11. Ces integrales sur les mailles ne sont non nulles que sur les mailles qui contiennent le n ud i. On ne les calcule donc que sur les mailles qui contiennent ce n ud. Les inconnues qui interviennent dans l'equation i sont donc uniquement les T (k) tels que le n ud k appartient a une maille contenant le n ud i. Il en est de m^ eme pour les integrales de surface. Il reste donc :
F (i) f dv =
XZ
m2I
F (i) f dv
ou I est l'ensemble des mailles contenant le n ud i.
6.3.3 forme 3
La troisieme formulation variationnelle du probleme thermique stationnaire a les m^ emes inconvenients que la premiere : les fonctions tests apparaissent sous la forme de laplacien, et demandent donc un espace de fonctions d'interpolation avec des fonctions de base a laplacien non nul partout. La methode de Galerkine est peu adequate pour cette formulation du probleme.
6.4 Application de la methode de Galerkine au probleme de statique de solide deformable

Pour les m^ emes raisons que dans l'exemple precedent, on choisit la formulation variationnelle: { Trouver le champ de deplacements tel que :
Sym grad dv +
A(
{z
Sd
n ds = ; } |
f dv ; {z
B( )
Sf
q0 ds 8 (M ) }
(6.3)
1 grad + gradT ou = 2 + Tr G, et = 2 { avec la condition aux limites:
) = 0 8N 2 Sd | (N ) ; {z 0(N } 0
C( )
Dans cette formulation, l'inconnue et les fonctions de leur gradient.
apparaissent tous les deux sous la forme
11. On rappelle aussi que, si on utilise des elements de reference, les integrales sur un element reel se ramenent a une integrale sur l'element de reference, et qu'elles peuvent ^ etre partiellement preparees a l'avance. Voir section 4.5 page 35.
J. Garrigues
61
D'autre part, appartenant a Fad , on a
Sd
n ds =
Z
Sd
0
n ds
Le probleme discretise est donc : Trouver les Ninc valeurs (i) tel que :
; e Sym gradF (i) dv +

V
Sd
(i) 0
Z Z e n ds = ; F (i) f dv ; F (i) q0 ds V Sf (6.4) (i) 8F base de Fad
ou e = 2
e + Tre G, et e =
N inc X j =1
(j )
Sym gradF (i)
La encore, les integrales de volume et de surface de l'equation i se reduisent a des integrales sur les seules mailles qui contiennent le n ud i.
6.5 L'assemblage dans la methode de Galerkine

L'assemblage est l'operation qui consiste a construire le systeme d'equations a resoudre. Le fait de choisir comme fonctions les fonctions de base de l'espace Fad simpli e beaucoup le calcul des integrales. Comme on l'a vu dans les exemples precedents, l'equation relative a la fonction de base F (i) ne fait intervenir en fait que les integrales sur les mailles qui contiennent le n ud i. On en deduit une methode pour construire la matrice K ] du systeme a resoudre :
2 u(1) 3 2 b(1) 3 K] 6 4 ... 7 5=6 4 ... 7 5

u(Ninc ) b(Ninc )
Pour chaque maille, on calcule les integrales de la formulation. Pour chaque F (i) dans la maille, on porte le coe cient de u(j ) en ligne i colonne j de K ] en l'ajoutant a ce qui s'y trouve deja, et on porte en ligne i de b] le resultat des integrales second membre toujours en l'ajoutant a ce qui s'y trouve deja. Informatiquement, on ameliore encore cette methode en utilisant les particularites de la matrice K] : { Si la matrice K ] est symetrique, il su t de n'en calculer (et de ne stocker) que le triangle superieur (ou inferieur). { La matrice K ] est creuse : sur la ligne i n'interviennent que les n uds des mailles contenant le n ud i. On utilise di erentes methodes de stockage pour eviter de remplir la memoire avec des zeros. Principalement, on trouve : { le stockage ((morse)) : On ne stocke que les termes non nuls ainsi que leur place dans K ], c'est a dire la ligne et la colonne. { le stockage ((ligne de ciel)) : On ne stocke que la partie de ligne de la matrice (ou du triangle, si elle est symetrique), depuis la diagonale jusqu'au dernier terme non nul. { le stockage ((bande)) : On ne stocke que les diagonales non nulles du triangle superieur (ou inferieur). Le nombre de diagonales est appele largeur de bande.
62
J. Garrigues
6.6. QUELQUES INDICATIONS SUR L'INTEGRATION NUMERIQUE

Il existe des algorithmes de renumerotation des n uds pour reduire la largeur de bande ou la longueur des lignes de ciel. Le plus utilise est l'algorithme de Gibbs 12. Il consiste a renumeroter les n uds pour que les di erences de numeros entre deux n uds voisins soient aussi petites que possible. Les termes non nuls de K ] sont donc groupes autour de la diagonale principale. Cette disposition ameliore tres sensiblement le temps de resolution du systeme lineaire pour certains solveurs 13. Ces choix d'options (type de stockage, renumerotation, et algorithme de resolution) etant souvent proposees dans les logiciels, il est indispensable de les comprendre pour les choisir a bon escient 14 .
Remarque :
Les conditions aux limites sont souvent plus complexes que celles qui ont ete donnees dans les exemples proposes precedemment. On peut par exemple imposer des relations entre certains degres de liberte 15. Le sous-espace Fad est alors plus di cile a construire. En fait, on ne le construit pas 16 et on rajoute des equations pour ((forcer)) les conditions aux limites. Le nombre d'inconnues est donc Nddl et non Ninc. Il existe en analyse numerique des methodes pour resoudre un systeme algebrique avec des contraintes sur la solution: on peut citer (entre autres) la methode des multiplicateurs de Lagrange et la methode de penalisation 17. Certains logiciels o rent ces choix en option.
6.6 Quelques indications sur l'integration numerique

Comme on l'a vu precedemment, la construction des termes de la matrice K ] consiste essentiellement en des calculs d'integrales sur des elements. Par exemple, dans le probleme thermique, les equations (6.3.2) page 60 montrent qu'on a a calculer les integrales
gradF (i)
N inc X j =1
T (j )gradF (j ) dv
ou les n uds i et j appartiennent a la maille k . On a donc a calculer des integrales du type
T (j ) gradF (i) gradF (j ) dv
On peut les ramener a une integrale sur un element de reference :
12. C'est un algorithme heuristique. Il diminue la largeur de bande mais il ne fournit pas necessairement la largeur minimale. Dans les maillages simples, une numerotation manuelle astucieuse peut parfois fournir de meilleurs resultats. L'e cacite de la renumerotation est sans in uence sur la precision du resultat, mais pour certains algorithmes de resolution, elle a une forte in uence sur le temps de calcul. 13. Le lecteur interesse est invite a se reporter a un cours d'analyse numerique 14. par exemple, il est inutile de renumeroter si le stockage est de type morse et si le solveur est iteratif. 15. Par exemple, pour simuler des n uds lies a un solide indeformable, ou pour simuler des liaisons complexes. 16. Fad n' a ete introduit dans ce chapitre que pour des commodites de raisonnement. 17. Pour plus de precisions sur la resolution de systemes algebriques ((sous contraintes)), se reporter a un cours d'analyse numerique. En fait, on a a resoudre un systeme algebrique de Nddl equations, comprenant celles issues de la formulation variationnelle, et celles issues des conditions aux limites. Si on utilise la methode de Galerkine, le sous-systeme issu de la formulationvariationnelle est symetrique,mais l'ajout des conditions aux limites detruit la symetrie du systeme complet. C'est cette perte de symetrie qui justi e le traitement special des conditions aux limites sous forme de ( (contraintes) ). Pour garder la symetrie, on pourrait satisfaire aux conditions aux limites sous forme variationnelle sur des elements frontiere, mais apparemment ce n'est pas l'usage. Cette maniere de voir est probablement un heritage des interpretations mecaniques des fonctions test.
J. Garrigues
T (j ) gradfe(i) grad ;1
gradfe(j) grad ;1 det J ] dx1dx2dx3
63
e(i)(x1 x2 x3) est la fonction d'interpolation sur l'element de reference obtenue pour une valeur ou f de 1 au n ud i, et 0 pour les autres 18 : e(i) (x1 x2 x3) = P (i)(x1 x2 x3) f
Les integrales a calculer ont donc des integrandes compliques (generalement une fraction rationnelle) 19. Plut^ ot que de les calculer exactement, on prefere les calculer numeriquement. La methode la plus souvent utilisee dans les logiciels est la methode des points de Gauss et ses variantes. Elle consiste a evaluer une approximation de l'integrale a calculer sous la forme
h(m) d =
n X q=1
wq h(mq )
ou mq sont des points bien choisis dans l'element de reference, appeles points de Gauss, et wq des coe cients numeriques bien choisis, tels que la formule ci dessus integre exactement un polyn^ ome de degre n donne 20. On ne donnera pas ici le detail de la maniere dont on choisit les points de Gauss mq et les valeurs des coe cients wq . Le lecteur est invite a se referer a un cours d'analyse numerique pour plus de precisions. La plupart des logiciels proposent des options de choix du nombre de points d'integration dans les elements, et ce qui precede n'a pour but que de donner la signi cation de ces options. Les positions des points de Gauss sont predeterminees dans les elements de reference 21 . On doit toutefois signaler un abus de langage couramment fait dans les logiciels: lorsque le nombre de points de Gauss est su sant pour integrer exactement la fonction d'interpolation de reference, l'integration est quali ee d'((exacte))22. Lorqu'on utilise un nombre de points de Gauss inferieur, on dit que l'element est sous integre 23.
6.7 Quelques indications sur les solveurs

Un probleme traite par elements nis se ramene toujours a la resolution de systemes d'equations algebriques (lineaires ou non, symetriques ou non). Ces systemes sont generalement de grande taille. Pour les problemes discretises par la methode de Galerkine, le systeme d'equations est symetrique. On ne trouvera pas ici un expose detaille des methodes de resolution, mais comme beaucoup de logiciels proposent des choix sur l'algorithme de resolution, on donne quelques indications pour guider ce choix. Le lecteur insatisfait est invite a se reporter a un cours d'analyse numerique.
18. ATTENTION ! dans la formule qui suit, les numeros de n uds sont des numeros dans l'element de reference. C'est le logiciel qui assure cette correspondance. 19. On y trouve par exemple det J ] 1 et grad 1 20. et evidemment tous les polyn^ omes de degre inferieur a n. 21. La position des points de Gauss dans l'element reel ou dans l'element de reference n'est pas toujours bien precisee dans les notices de logiciels ou dans les chiers de resultats. Ceci peut rendre di cile l'exploitation des resultats. 22. Mais on a vu dans ce qui precede que les integrales a calculer sont des integrales des interpolations et de leurs gradients sur l'element reel, dont les integrandes ne sont pas toujours polyn^ omiaux a cause de la transformation . Cette inexactitude disparaitrait si on se passait d'elements de reference, et si on construisait les fonctions d'interpolation et les points de Gauss sur l'element reel. On pourrait m^ eme se passer de points de Gauss et calculer ces integrales de polyn^ omes exactement. 23. On utilise souvent cette option pour masquer les anisotropies induites par les fonctions d'interpolation d'un maillage trop regulier avec des fonctions de forme trop simples, ou pour ((ameliorer la convergence)) de problemes di ciles.
; ;
64
J. Garrigues
6.7. QUELQUES INDICATIONS SUR LES SOLVEURS
6.7.1 Resolution des systemes lineaires

On peut distinguer deux classes de methodes 24 : { les methodes directes : Elles aboutissent a la solution en un nombre ni d'operations. { les methodes iteratives : La solution est atteinte en un nombre in ni d'operations convergeant vers la solution. On arr^ ete les iterations lorque qu'on estime qu'on est su samment pres de la solution. Si on utilise une telle methode de resolution, le logiciel demande generalement de preciser le critere d'arr^ et des iterations.
Les methodes directes:

Si la matrice n'est pas symetrique, on utilise la methode d'elimination de Gauss. Si la matrice est symetrique, on utilise generalement une decomposition de Choleski. Di erentes versions de ces methodes sont adaptees sont adaptees aux di erents modes de stockage de la matrice K ]. Les methodes directes sont employees pour des systemes de taille raisonnable (quelques centaines ou quelques milliers d'inconnues 25). Le temps de calcul est sensiblement proportionnel au cube du nombre d'inconnues.
Les methodes iteratives:

Que la matrice soit symetrique ou non, on part d'une solution initiale 26 et on l'ameliore par etapes successives. Si le systeme a resoudre est K ] U ] = B] on cherche a diminuer la quantite 27 kRk = k K ] Ui] ; B ] k, Ui] etant la solution approchee a l'iteration i. Suivant les methodes, la norme k : k peut ^ etre la norme euclidienne ou la norme du max. La plupart des methodes sont basees sur une descente de gradient : on fait evoluer les Ui ] dans le sens de ;gradkRk. Les di erentes variantes se distinguent par la vitesse a laquelle on descend dans cette direction. Le critere d'arr^ et est soit kRk < 1, soit kUj ; Uj ;1 k < 2. Le choix de 1 ou 2 depend de ce qu'on considere comme petit 28. Les methodes iteratives sont employees pour les systemes de grande taille, ou le temps de calcul serait prohibitif avec une methode directe. C'est un compromis entre l'exactitude (determinee par le critere d'arr^ et des iterations) et le co^ ut du calcul.
6.7.2 Resolution des systemes non lineaires

Un systeme non lineaire peut toujours s'ecrire sous la forme K (U )] U ] = B (U )] On le resout par approximations successives : Ui+1 ] est solution de K (Ui )] Ui+1] = B (Ui )]
24. Cette distinction dichotomique est a nuancer: La methode du gradient conjugue aboutit en un nombre ni d'operations, mais on s'arr^ ete souvent avant la n. 25. Ces limites sont seulement indicatives. Elles dependent de la machine dont on dispose. 26. le plus souvent c'est le vectur nul. 27. souvent appelee residu. 28. par exemple, dans un probleme de mecanique, 1 , est une force negligeable, 2 est un deplacement negligeable. Il est donc absurde de xer un sans considerer des unites qui ont ete utilisees pour de nir le probleme numerique!
J. Garrigues
65
On voit donc que la resolution d'un systeme non lineaire se ramene a la resolution d'une suite de problemes lineaires. Di erentes variantes existent pour ameliorer la convergence 29 ou diminuer le volume des calculs 30. La convergence de ces methodes est heuristique et peut poser probleme 31: il faut parfois ajuster des parametres algorithmiques pour obtenir une solution. D'autre part, la solution d'un systeme non lineaire n'est pas toujours unique 32 et il n'est pas toujours evident de savoir vers laquelle des solutions l'algorithme a converge.
6.8 La discretisation du temps

Tous les exemples traites precedemment etaient des problemes stationnaires. Dans la plupart des logiciels d'elements nis, le temps est traite comme un parametre33. On resout donc un probleme a chaque pas de temps, en partant d'un instant initial ou les derivees premieres et secondes sont supposees connues : ce sont les conditions initiales. Ce sont donc des methodes de di erences nies en temps. On peut distinguer deux grandes classes de methodes : { Les methodes explicites : les derivees premieres et secondes deduites du pas precedent sont utilisees pour calculer la valeur suivante du champ inconnu. La stabilite 34 de la solution demande generalement des pas de temps tres petits. Ces methodes sont le plus souvent utilisees dans les logiciels de modelisation de crash. { Les methodes implicites: les derivees premieres et secondes en n de pas 35 sont utilisees pour calculer la valeur suivante. La stabilite est meilleure 36, mais la resolution de chaque pas de temps exige un calcul iteratif, avec les problemes associes a tout calcul iteratif : les parametres algorithmiques et le choix de criteres d'arr^ et.
29. Methode de Newton-Raphson: on construit une matrice tangente KT ] pour remplacer K ] 30. Par exemple, on peut ne pas mettre a jour la matrice K ] ou KT ] a toutes les iterations. On evite ainsi des assemblages. 31. Beaucoup d'exposes sur la methode des elements nis orientes ((mecanique)) cherchent a justi er les iterations des algorithmes de resolution de systemes non lineaires par des interpretations mecaniques: ((on charge la structure par increments de charge)), au cours d'un ((temps ctif)), ou on considere des ((lois de comportement tangentes)). 32. Elle peut notamment comporter des bifurcations. 33. On pourrait avoir un point de vue di erent, mais ce n'est pas l'usage: un probleme d'evolution est solution d'une equationdi erentielledont les variables sont les variablesd'espace et le temps. On peut donc en tirer des formes variationnelles en traitant le temps comme les autres variables. Il reste que pour les problemes tridimensionnels evolutifs, on ne dispose pas (ou pas encore ...) de mailleur dans un espace a 4 dimensions, et que le nombre d'inconnues des systemes a resoudre augmente fortement. 34. Ne pas confondre avec l'exactitude de la solution! La stabilite ne traduit que la sensibilite des solutions approchees a la valeur des pas de temps. 35. ou en milieu de pas, ou en tout autre point, d'ou les di erentes variantes. 36. Pour certaines variantes, la solution est toujours stable.
66
J. Garrigues
Chapitre 7
Le post-traitement
La solution approchee d'un probleme d'elements nis est donc un Ninc -uplet de valeurs de ddl aux n uds du maillage. Pour un probleme d'evolution, on a un Ninc -uplet pour chaque pas de temps. Le nombre de valeurs dans un resultat de calcul est donc tres grand ! Par exemple, dans un probleme de mecanique de solide deformable tridimensionnel evolutif, si on a 1000 n uds, chaque n ud a 3 ddl (les composantes du champ de deplacement), et si on a calcule pour 100 pas de temps, la solution est constituee de 300 000 nombres! D'autre part, le champ solution n'est pas toujours la grandeur physique pertinente a observer. L'ingenieur a besoin d'observer des grandeurs qui se deduisent de la solution. Par exemple, le champ des deplacements d'un solide deformable n'est pas une information pertinente pour savoir si le solide respecte les criteres de resistance. On a donc besoin : { de calculer des grandeurs physiques deduites de la solution: ce sont les grandeurs derivees. { d'outils pour analyser e cacement tous ces nombres.
7.1 Les grandeurs derivees

7.1.1 Quelques considerations generales
Les grandeurs derivees dont on a besoin dependent fortement de la physique du probleme traite : Un thermicien est sans doute interesse par les ux thermiques dans le domaine et a travers ses frontieres un mecanicien des solides est sans doute interesse par le champ de tenseur des deformations 1 , par le champ de tenseurs des contraintes 2 (pour veri er un critere de resistance par exemple), et par les actions de liaison sur la frontiere 3. Il faut donc calculer ces grandeurs derivees a partir du champ solution approche primaire fourni par les valeurs aux n uds. Si le logiciel d'elements nis est specialise dans un certain type de physique, (mecanique des solides, mecanique des uides, thermique, electromagnetisme ...) il ne propose que les grandeurs derivees susceptibles d'interesser les specialistes de ce domaine, avec un vocabulaire adapte, qui evite a l'utilisateur de se souvenir des de nitions des grandeurs qui l'interessent. C'est le plus souvent le cas.
1. qui se calcule a partir du gradient du champ de deplacement. 2. qui se calcule algebriquement a partir du champ des deformations 3. qui ne sont que le ux du tenseur des contraintes a travers la frontiere:
J. Garrigues
67
CHAPITRE 7. LE POST-TRAITEMENT
Cependant, le post-traitement se ramene le plus souvent a des calculs de gradient, de divergence, de rotationnel et de laplacien du champ solution (ou de grandeurs derivees). On peut donc concevoir un logiciel avec un post-traitement assez general, et ce serait a l'utilisateur de fabriquer les grandeurs derivees qu'il desire observer. Certains logiciels ont cette pretention, avec plus ou moins de souplesse et un vocabulaire plus ou moins heureux. Tous les operateurs di erentiels usuels se ramenent a des calculs de gradient. On n'examinera donc que ce cas. Le gradient d'un champ vectoriel ou tensoriel se ramene au calcul du gradient des composantes quand celles-ci sont donnees dans un repere orthonorme xe 4 . On se limite donc au cas d'un champ scalaire. La generalisation se fait sans di culte.
7.1.2 Evaluation du gradient d'un champ scalaire

La solution approchee primaire est une suite de valeurs aux n uds u(1) : : : u(Nddl) qui permet de reconstruire un champ solution approche dans chaque element. On considere un element et u(1) : : : u(n) la suite de ses valeurs aux n n uds 5 de l'element. La valeur du champ sur l'approximation geometrique de l'element est :
e(M 0) = fe(m) = F
n X i=1
u(i) P (i)(m)
ou les P (i)(m) sont les polyn^ omes de base de l'interpolation sur l'element de reference. Les formules donnant le gradient ont ete donnees en section 4.5 page 35 et suivantes. Ces formules montrent que { la repartition du gradient a l'interieur de l'element n'a pas beaucoup de signi cation physique : elle ne depend que du choix des polyn^ omes d'interpolation dans l'element de reference et de la transformation , qui sont tous les deux arbitraires. L'equation di erentielle initiale n'est pas strictement respectee a l'interieur de l'element. { le gradient sur l'element reel est di erent du gradient sur l'element de reference. D'autre part, si l'utilisation d'elements conformes garantit bien la continuite C0 de la solution F , la continuite du gradient n'est pas assuree. M^ eme si les ecarts entre la solution approchee F et la solution exacte u sont faibles, on ne peut mathematiquement rien a rmer sur les ecarts entre les derivees. Toutefois, si on postule que la solution exacte est ((assez reguliere)), on peut esperer que ces ecarts soient faibles. On cherche donc souvent a donner un gradient plus signi ant: { Si les integrales sur les elements ont ete calculees avec des points de Gauss, on peut donner les valeurs du gradient aux points de Gauss 6. { On peut aussi donner un gradient moyen 7. Le gradient d'un champ de valeurs d'un ddl est alors decrit par une fonction vectorielle constante par element.
4. Ce qui est toujours le cas dans les logiciels d'elements nis, sauf dans les problemes axisymetriques. 5. Attention, il s'agit maintenant de la numerotation interne a l'element. 6. Ceci revient a calculer le gradient en des points non nodaux, ce qui pose un probleme pour leur representation graphique, notamment si on veut presenter des isocouleurs. D'autre part, la position des points de Gauss dans l'element reel n'est pas toujours donnee dans les chiers resultats! R 7. On peut le de nir par V (1 i ) gradFds, ou V ( i ) est le volume, l'aire ou la longueur de l'element , ou par la moyenne des valeurs aux points de Gauss.
68
J. Garrigues
7.2. LES OUTILS DE VISUALISATION

Dans tous les cas, le gradient etant discontinu aux frontieres des elements, on a du mal a de nir une valeur du gradient aux n uds 8 . Si le logiciel donne de telles valeurs, il s'agit d'une moyenne des gradients moyens dans les elements qui contiennent ce n ud 9. Toutes ces erreurs restent moderees si les sauts interelements sont faibles. Dans le cas contraire, il faut considerer ces valeurs avec circonspection et eventuellement refaire un calcul avec un maillage ra ne dans cette region.
7.2 Les outils de visualisation

Etant donne le grand nombre de valeurs a examiner pour analyser une solution et ses ses grandeurs derivees, la plupart des logiciels proposent des outils de visualisation graphique plus ou moins riches pour representer les resultats. Principalement on trouve : des outils de visualisation de champs scalaires : Ce sont normalement des valeurs scalaires (temperature, energies, normes de vecteurs, invariants ou valeurs propres de tenseurs). On peut aussi representer des composantes de vecteurs ou de tenseurs 10 . { Si est un domaine lineique, ces champs scalaires sont normalement representes par des courbes. { Si est un domaine surfacique, on utilise un reseau de courbes isocouleurs { Si est un domaine volumique, on utilise des surfaces isocouleurs. Pour observer l'interieur d'un champ volumique, beaucoup de logiciels permettent de faire des coupes de maniere plus ou moins souple, et certains visualisateurs permettent de jouer sur la transparence pour observer l'interieur. des outils de visualisation de champ vectoriel : On represente des vecteurs colineaires a la grandeur et de longueur proportionnelle a sa norme. Cette representation n'est vraiment lisible que pour des domaines lineiques ou surfaciques plans. Dans les autres cas, les vecteurs sont vus a l'ecran en perspective, ce qui rend di cile l'estimation de leur intensite. Pour paillier partiellement a ce probleme, l'intensite du vecteur est a la fois representee par la longueur (faussee par la projection) et une couleur. des outils de visualisation de champ tensoriel : On represente les tenseurs symetriques par trois vecteurs dans les directions propres. L'intensite des vecteurs est la valeur propre associee. A ces outils de visualisation, qui permettent une vue globale des champs, sont generalement associees des recherches de maxima ou de minima pour aider a la recherche de valeurs critiques. Il est parfois indispensable de faire son propre programme d'exploitation des resultats pour evaluer des grandeurs derivees qui ne sont pas prevues dans le logiciel.
8. Sauf si on a utilise une formulation a plusieurs champs dans laquelle ce gradient est une inconnue. 9. Cette moyenne peut ^ etre ponderee par le volume, l'aire ou la longueur, des elements qui contiennent le n ud. 10. Bien que des valeurs de composantes dans un systeme d'axes arbitraire ne presentent en general aucun inter^ et!
J. Garrigues
69
70
J. Garrigues
Annexe A
Formules de transformation d'integrales de volume

On suppose admis le theoreme d'Ostrogradski :
div V dv =
Z
S
V n ds
Si f est un champ scalaire et V un champ vectoriel, on veri e facilement l'identite: div(f V ) = V grad f + f div V et donc
En appliquant le theoreme d'Ostrogradski,
f div V dv = ; f div V dv = ;
V grad f dv + div(f V ) dv
V
Pour V = grad g on obtient :
Z V grad f dv + f V n ds
S
f g dv = ;
Z
V
grad g grad f dv + f grad g n ds

S
Si V et W sont deux champs vectoriels, on veri e facilement l'identite: div(V ^ W ) = W rot V ; V rot W et donc
V rot W dv = W rot V dv ; div(V ^ W ) dv

V V
V
J. Garrigues
V rot W dv = W rot V dv ; (V ^ W ) n ds
V S
71
ANNEXE A. FORMULES DE TRANSFORMATION D'INTEGRALES DE VOLUME

Si V est un champ vectoriel et T un champ tensoriel du second ordre, on veri e facilement l'identite: div(V T ) = V div T + grad V T et donc
Z
V
V div T dv = ; grad V T dv + div(V T ) dv

V V
V div T dv = ; grad V T dv + V T n ds
V S
Pour T = gradW ,
W dv = ; grad V grad W dv + V grad W n ds

V S
Si T est un champ tensoriel du second ordre et S un champ tensoriel du troisieme ordre, on veri e facilement l'identite : div(T S ) = gradT S + T div S et donc
Z
V
T div S dv = ; gradT S dv + div(T S) dv

V V
T div S dv = ; grad T S dv + T S n ds
V S
Pour S = gradR,
Z Z R dv = ; grad T grad R dv + T grad R n ds

V S
72
J. Garrigues
Annexe B
Rappels sur les surfaces

B.1 Caracterisation des surfaces
Une surface S est un ensemble de points a deux parametres :
S = M 2 E3
OM = OM (x1 x2)
Si on repere les points de l'espace physique tridimensionnel par un systeme ce coordonnees cartesiennes orthonormees, d'origine O et de base fE k g les coordonnee de M sont : ~1 (x1 x2) X1 = X ~2 (x1 x2) X2 = X ~3 (x1 x2) X3 = X ce qui est la de nition parametrique de S . Les deux vecteurs de la base naturelle de S en M sont : OM e1 = @@x 1 OM e2 = @@x 2
Ils de nissent le plan tangent au point M . En general ces deux vecteurs ne sont ni orthogonaux ni normes. Leurs composantes dans fE k g sont :
2 e1]E = 6 4
~1 @X 1 @x ~2 @X 1 @x ~3 @X @x1
3 7 5
2 e2]E = 6 4
~1 @X 2 @x ~2 @X 2 @x ~3 @X @x2
3 7 5
Dans (M ), on peut de nir des vecteurs tangents et des tenseurs tangents. Si V = v e Les v sont appeles composantes contravariantes de V . Les v = V e sont appeles composantes covariantes de V . Si T est un tenseur tangent du second ordre c'est a dire
1 ^ e2 La normale unitaire a S en M est n = ke e ^ e k. 1 2
T (V W ) = t v w = t v w = t v w = t v w 2 8V W 2
J. Garrigues
73
ANNEXE B. RAPPELS SUR LES SURFACES

ses 4 composantes deux fois covariantes sont : t = T (e e ) etc. On de nit le tenseur metrique de surface A, du second ordre symetrique : A(V W ) = V W 8V W 2 Ses composantes deux fois covariantes dans la base naturelle sont :
1 e1 e1 e2 a ]= e e2 e1 e2 e2 Ses composantes deux fois contravariantes dans la base naturelle sont 1 : a ] = a ];1 Ses composantes mixtes sont la matrice unite. On de nit aussi le tenseur d'orientation de surface E , du second ordre antisymetrique: E(V W ) = (V ^ W ) n 8V W 2 Ses composantes covariantes dans la base naturelle sont : pa e ] = 0 1 ;1 1 et ses composantes contravariantes dans la base naturelle sont : 1 e ]= 0 1 p ;1 1 a ou a = det a ]. Le produit vectoriel de deux vecteurs tangents est : V ^ W = (V E W ) n On a aussi n ^ V = ;E V
En etudiant les variations de la normale unitaire n quand M varie sur S , on met en evidence l'existence du tenseur de courbure normale B . Ce tenseur est du second ordre, symetrique de ni par dn = ;B dM Ses composantes covariantes dans la base naturelle sont : @ 2 OM b = n @x @x Les derivees des vecteurs de base sont : @e @x = ; e + b n ou 1 a @a + @a ; @a ; =2 @x @x @x Propriete : ; = ; .
1. Attention: Cette propriete n'est vraie que pour le tenseur metrique.
74
J. Garrigues
B.2. CHANGEMENTS DE BASE
B.2 Changements de base

Certains vecteurs ou tenseurs de ce chapitre sont de nis par leurs composantes dans la base naturelle en M , qui n'est en general ni orthogonale ni normee. On peut donner leurs composantes dans une base xe orthonormee 2. Pour les vecteurs et les tenseurs dont la de nition est le resultat d'operations tensorielles, on peut e ectuer ces operations dans toute base. Pour les autres 3 , on donne les formules de changement de base : Si V un vecteur tangent connu par ses composantes contravariantes sur la base naturelle : ses composantes dans la base E k sont :
V =v e
V ]E = e1]E e2]E
v =a v soit encore sous forme matricielle: v1 = a ] v1 v2 v2
v1 v2
Si V est connu par ses composantes covariantes, ses composantes contravariantes sont :
Si T est un tenseur du second ordre tangent connu par ses composantes contravariantes, ses composantes dans la base E k sont la matrice 3 3:
T ]E = e1]E e2]E t ] e1]E e2]E

t =a t a soit encore sous forme matricielle: t ]= a ] t ] a ]
Si T est connu par ses composantes covariantes, ses composantes contravariantes sont :
Pour les tenseurs d'ordre superieur, on ne peut plus ecrire de matrice pour representer les sommations 4. On laisse les resultats sous forme de sommations.
B.3 Operateurs di erentiels tangents sur une surface

B.3.1 Gradient surfacique
Si f (M ) = f (x1 x2) est un champ scalaire de ni sur S , le gradient surfacique de f est le vecteur de ni par : df = gradf dM
2. On l'a deja fait plus haut pour les vecteurs de la base naturelle. 3. On se limite ici a ce qui est strictement utile. 4. Il faudrait des ((matrices cubiques)) pour les tenseurs d'ordre 3 !
J. Garrigues
75

Ses composantes covariantes dans la base naturelle sont :
1 (gradf ) = @x @f
@f
@x2
Si V (M ) = V (x1 x2) est un champ de vecteurs tangents, de ni sur S , le gradient surfacique de V est le tenseur du second ordre tangent de ni par : Tang(dV ) = gradV dM Ses composantes covariantes dans la base naturelle sont : @v ; v ; gradV = @x ATTENTION ! Le gradient surfacique ne donne que la partie tangente de la variation du champ tangent V (M ). Quand la surface a une courbure non nulle, les variations d'un champ de vecteurs tangents ont une partie normale. Si X est un champ de vecteurs non tangent de ni sur S , sa derivee dans la direction unitaire tangente t est : Dt X = gradV t + (V B t) n
Si T (M ) = T (x1 x2) est un champ de tenseurs tangents, de ni sur S , le gradient surfacique de T est le tenseur du troisieme ordre tangent de ni par : Tang(dT ) = gradT dM Ses composantes covariantes dans la base naturelle sont : gradT = @t @x ; t ; ; t ; On generalise au gradient de champs tensoriels tangents de tous ordres.
B.3.2 Divergence surfacique

La divergence surfacique d'un champ vectoriel tangent V (M ) est le champ scalaire : f d ivV = gradV A La divergence surfacique d'un champ tensoriel tangent du second ordre T (M ) est le champ vectoriel tangent: g A d ivT = On generalise aux tenseurs tangents de tous ordres.
] ] g radT ]
B.3.3 Rotationnel surfacique

Le rotationnel surfacique d'un champ vectoriel tangent V (M ) est le champ scalaire : rf otV = ;gradV E Si X (M 0 ) est un prolongement di erentiable de V (M ) hors de la surface S , on a la relation: rf otV (M ) = rotX (M ) n On l'identite suivante fV rf ot(n ^ V ) = div (B.1)
76
J. Garrigues
B.4. FORMULES D'ANALYSE SUR LES SURFACES
B.3.4 Laplacien surfacique

Le laplacien d'un champ scalaire f (M ) est le champ scalaire de ni par :
Le laplacien d'un champ vectoriel tangent (M ) est le champ vectoriel tangent de ni par : = Le laplacien d'un champ tensoriel tangent du second ordre T (M ) est le champ tensoriel tangent du second ordre de ni par :
] V g eV d ] iv g radV
f ef = d iv gradf eT = d g iv gradT
On generalise aux tenseurs tangents de tous ordres.
B.4 Formules d'analyse sur les surfaces

Soit C une courbe fermee de la surface, entourant la portion de surface S . On de nit sur tout point N de C le triedre de Darboux-Ribeaucourt. On note :
n n t
nt
Fig. B.1 {
Domaine surfacique
n la normale unitaire a S , nt la normale unitaire exterieure a C et et tangente a S , t = n ^ nt la tangente unitaire a C. Le triedre de Darboux-Ribeaucourt ft n nt g est orthonorme.
L'element de surface est : ds = E (e1 e2 ) dx1 dx2 = (e1 ^ e2) n dx1 dx2 = e dx dx On suppose admise la formule de Stockes classique :
C 0 ou X (M ) est un champ vectoriel de ni dans l'espace.

J. Garrigues
rotX (M 0)
n ds = X (M 0) tdl
77
B.4.1 Formule de Stockes pour les surfaces

Il su t d'appliquer la formule de Stockes classique a un champ vectoriel X (M 0 ) qui est le prolongement di erentiable d'un champ vectoriel tangent V de ni sur la surface. On a alors :
S Z
rotX (M 0) n ds = X (M 0) t dl rotV (M ) n ds = V (M ) t dl
C S
ZC
rf otV (M ) ds =
Z
C
V (M ) t dl
B.4.2 Formule d'Ostrogradski pour les champs vectoriels tangents

Il su t d'appliquer la formule de Stockes au vecteur tangent W = n ^ V . on etablit facilement que
f rf otW (M ) = d ivV (M ) Z
C
On obtient alors:
Z
S
f V (M ) ds = div
V (M ) nt dl
B.4.3 Formule d'Ostrogradski pour les champs tensoriels tangents du second ordre
Soit T un tenseur R de surface du second ordre. On se propose de calculer le (( ux)) de T a travers la courbe C : C T nt dl. C'est un vecteur non tangent. Pour calculer cette integrale vectorielle on la projette sur la base xe E k .
Ek Ek
Z
C
T nt dl = Ek (T nt) dl
C
En developpant les calculs en composantes sur la base naturelle on trouve :
(La courbure appara^t car on a a deriver les vecteurs de base.) et donc :
T nt dl = Ek T nt dl =
f T + B T n ds div
f d ivT + B T n ds
B.5 Transformations d'integrales

Dans les formules qui suivent : { f et g sont des champs scalaires de nis sur S ,
78
J. Garrigues
B.5. TRANSFORMATIONS D'INTEGRALES

{ V et W sont des champs vectoriels tridimensionnels de nis sur S , { V et W sont les parties tangentes de V et W , { T et R sont des champs tensoriels tridimensionnels du second ordre de nis sur S , { T et R ont les parties tangentes de T et R, { S est un champ de tenseur de surface du troisieme ordre. 0n veri e facilement l'identite : f f V ) = V gradf + f div fV div( et donc
en appliquant le theoreme d'Ostrogradski pour les champs vectoriels tangents :
f V ds = f div
soit encore :
f fd iv V ds =
Z
S
Pour V = gradg, on obtient :
f V ds = f div Z
S
] Z ] ; V g rad Z ] ; V g rad Z ] ; V g rad

S S S
f ds + f ds + f ds +
f d iv(f V ) ds
f V nt dl f V nt dl
Z
C
f e g ds = ;
gradg gradf ds + f gradg nt dl

C
] ]
0n veri e facilement l'identite : g f V T) = V d iv T + gradV T div( et donc
Z Z g f V d iv T ds = ; grad V T ds + d iv(V T ) ds
S S
soit encore :
g V d iv T ds = g V d iv T ds =
Z
S
Pour T = grad W , on obtient :

S
] Z
] Z ] ; g rad V T Z ] ; g rad V T
S S
ds +
V T nt dl
ds +
Z
C
V (T nt) dl
Z Z V e W ds = ; grad V grad W ds + V grad W nt dl

S C
] ]
J. Garrigues
79

0n veri e facilement l'identite : et donc
g g d iv(T S) = T d iv S + grad T S
S
Z
S
Z Z g g T d iv S ds = ; grad T S ds + d iv(T S ) ds
S
Pour S = grad R, on obtient :

S
]Z
Z Z g T d iv S ds = ; grad T S ds + T S nt dl
S C
Z Z T e R ds = ; grad T grad R ds + T grad R nt dl

S C
] ]
80
J. Garrigues
Table des matieres

1 Introduction 2 Grandes lignes de la methode des elements nis
2.1 Expose de la demarche . . . . . . . . . . . . . . 2.2 Un exemple detaille . . . . . . . . . . . . . . . 2.2.1 Choix du maillage . . . . . . . . . . . . 2.2.2 Choix des n uds et des champs locaux 2.2.3 Formulation variationnelle du probleme 2.2.4 Discretisation du probleme . . . . . . . 2.2.5 Resolution . . . . . . . . . . . . . . . . . 2.2.6 Examen des resultats . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5
7 9 10 10 11 12 14 14
3 Le maillage
3.1 De nition du domaine . . . . . . . . . . . . . . . . . . . . . . 3.1.1 De nition du domaine dans le logiciel d'elements nis 3.1.2 Importation du domaine d'un autre logiciel . . . . . . 3.2 De nition des sous domaines i . . . . . . . . . . . . . . . . . . 3.2.1 Quelques ((evidences)) . . . . . . . . . . . . . . . . . . . 3.2.2 Les mailles disponibles . . . . . . . . . . . . . . . . . . . 3.2.3 Les methodes de maillage . . . . . . . . . . . . . . . . . 3.2.4 L'importation d'un maillage . . . . . . . . . . . . . . . . . . . . . .
17
17 17 18 18 18 19 20 21
4 Les elements et leur espace de fonctions d'interpolation
4.1 De nition des mailles de reference . . . . . . . . . . . . . . . 4.2 Interpolation polyn^ omiale sur les mailles de reference . . . . . 4.2.1 Rappels sur les bases de polyn^ omes . . . . . . . . . . . 4.2.2 Interpolation polyn^ omiale sur une maille de reference 4.2.3 Exemple 1 (standard) . . . . . . . . . . . . . . . . . . 4.2.4 Exemple 2 (non standard) . . . . . . . . . . . . . . . .
23
24 25 25 26 27 29
J. Garrigues
81
TABLE DES MATIERES

4.3 Continuite C0 inter-elements : Les elements conformes . . . . . . . 4.4 Interpolation dans l'element reel . . . . . . . . . . . . . . . . . . . 4.4.1 Principe de construction . . . . . . . . . . . . . . . . . . . . 4.4.2 Choix de la transformation : l'approximation geometrique 4.4.3 En resume ... . . . . . . . . . . . . . . . . . . . . . . . . . . 4.5 Derivees et integrales dans l'element reel . . . . . . . . . . . . . . 4.5.1 Mailles volumiques . . . . . . . . . . . . . . . . . . . . . . . 4.5.2 Mailles surfaciques . . . . . . . . . . . . . . . . . . . . . . . 4.5.3 Mailles lineiques . . . . . . . . . . . . . . . . . . . . . . . . 4.5.4 Remarque . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 31 31 32 35 35 36 37 39 40 41
5 Les formulations variationnelles
5.1 Rappel de quelques resultats d'analyse fonctionnelle . . . . . . . . 5.2 Formes variationnelles . . . . . . . . . . . . . . . . . . . . . . . . . 5.3 Transformations d'integrales . . . . . . . . . . . . . . . . . . . . . . 5.3.1 Transformations d'integrales sur un volume . . . . . . . . . 5.3.2 Transformations d'integrales sur une surface . . . . . . . . . 5.3.3 Transformations d'integrales sur une courbe . . . . . . . . . 5.4 Exemple 1: Le probleme thermique stationnaire . . . . . . . . . . 5.5 Exemple 2: Statique du solide deformable en petites deformations. 5.5.1 Premiere forme variationnelle . . . . . . . . . . . . . . . . . 5.5.2 Seconde forme variationnelle . . . . . . . . . . . . . . . . . 5.6 Exemple 3: Statique des coques . . . . . . . . . . . . . . . . . . . 5.6.1 Premiere forme variationnelle . . . . . . . . . . . . . . . . . 5.6.2 Seconde forme variationnelle . . . . . . . . . . . . . . . . . 5.6.3 Une variante de la formulation precedente . . . . . . . . . . 5.6.4 Remarques . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.7 Formulations a plusieurs champs . . . . . . . . . . . . . . . . . . . 5.8 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
43
43 44 45 45 46 47 47 49 49 50 52 53 53 54 54 54 56
6 La discretisation
6.1 Principe de la discretisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.2 La methode de Galerkine . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.3 Application de la methode de Galerkine au probleme thermique stationnaire . 6.3.1 forme 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.3.2 forme 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
57
57 59 59 59 60
82
J. Garrigues
TABLE DES MATIERES

6.4 6.5 6.6 6.7 6.8 6.3.3 forme 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Application de la methode de Galerkine au probleme de statique de solide deformable L'assemblage dans la methode de Galerkine . . . . . . . . . . . . . . . . . . . . . . Quelques indications sur l'integration numerique . . . . . . . . . . . . . . . . . . . Quelques indications sur les solveurs . . . . . . . . . . . . . . . . . . . . . . . . . . 6.7.1 Resolution des systemes lineaires . . . . . . . . . . . . . . . . . . . . . . . . 6.7.2 Resolution des systemes non lineaires . . . . . . . . . . . . . . . . . . . . . La discretisation du temps . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 61 62 63 64 65 65 66
7 Le post-traitement
7.1 Les grandeurs derivees . . . . . . . . . . . . . . . . 7.1.1 Quelques considerations generales . . . . . 7.1.2 Evaluation du gradient d'un champ scalaire 7.2 Les outils de visualisation . . . . . . . . . . . . . .
67
67 67 68 69
A Formules de transformation d'integrales de volume B Rappels sur les surfaces

B.1 Caracterisation des surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.2 Changements de base . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.3 Operateurs di erentiels tangents sur une surface . . . . . . . . . . . . . . . . . . . B.3.1 Gradient surfacique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.3.2 Divergence surfacique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.3.3 Rotationnel surfacique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.3.4 Laplacien surfacique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.4 Formules d'analyse sur les surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . . B.4.1 Formule de Stockes pour les surfaces . . . . . . . . . . . . . . . . . . . . . . B.4.2 Formule d'Ostrogradski pour les champs vectoriels tangents . . . . . . . . . B.4.3 Formule d'Ostrogradski pour les champs tensoriels tangents du second ordre B.5 Transformations d'integrales . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
71 73
73 75 75 75 76 76 77 77 78 78 78 78
J. Garrigues
83

La Méthode Des Éléments Finis, J. GARRIGUES

Transféré par

Droits d'auteur :

Formats disponibles

La Méthode Des Éléments Finis, J. GARRIGUES

Transféré par

Informations du document

Titre original

Copyright

Formats disponibles

Partager ce document

Partager ou intégrer le document

Options de partage

Avez-vous trouvé ce document utile ?

Ce contenu est-il inapproprié ?

Droits d'auteur :

Formats disponibles

La Méthode Des Éléments Finis, J. GARRIGUES

Transféré par

Droits d'auteur :

Formats disponibles

La methode des elements nis

Jean GARRIGUES [email protected] Janvier 2002

Ecole Superieure de Mecanique de Marseille

Ecole Superieure de Mecanique de Marseille

Exemple 1: Le probleme thermique stationnaire

Ecole Superieure de Mecanique de Marseille

Exemple 2: Le probleme elastique lineaire d'evolution

Ecole Superieure de Mecanique de Marseille

Grandes lignes de la methode des elements nis

Ecole Superieure de Mecanique de Marseille

CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS

Solution approchee discontinue

Solution approchee continue C0

Solution approchee continue C1

Ecole Superieure de Mecanique de Marseille

2.2. UN EXEMPLE DETAILLE

2.2 Un exemple detaille

Ecole Superieure de Mecanique de Marseille

CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS

2.2.1 Choix du maillage

2.2.2 Choix des n uds et des champs locaux

Ecole Superieure de Mecanique de Marseille

2.2. UN EXEMPLE DETAILLE

2.2.3 Formulation variationnelle du probleme

(x) h(x) dx = 0 8 (x) () h(x) = 0 8x 2

Ecole Superieure de Mecanique de Marseille

CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS

(x) dx = ; (x) dfdx

Z 1 d (x) df (x) df (x) dx + ( x ) dx dx dx Z01 d (x) 1

dx f (x) dx + (x) f (x)]0

Z 1 d (x) df (x) (x) ; dx + 0 ; (0) dfdx dx dx 0 Z 1 d (x)

dx f (x) dx + 2 (1) f (1) ; 0

(x) x dx = 0 8 (x) (2.9)

2.2.4 Discretisation du probleme

Ecole Superieure de Mecanique de Marseille

2.2. UN EXEMPLE DETAILLE

, 1 ,..., 6 de nies par la gure 2.5.

Fig. 2.5 { Les

sept fonctions i (x)

" ~ # Z 1 d i(x) df~(x) (x) dx ; i(0) df ; dx dx dx 0 x=0 Z 1 d i (x )

~ ~ dx f (x) dx + 2 i (1) f (1)

Ecole Superieure de Mecanique de Marseille

CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS

" ~ # Z d 1(x) df~(x) df (x) ; dx dx dx ; 1 (0) dx x=0 0 Z d 1(x)

~ ~ dx f (x) dx + 2 1(1) f (1)

2.2.6 Examen des resultats

2.2. UN EXEMPLE DETAILLE

Comparaison des solutions exacte et approchee

0.5 0 -0.5 -1 -1.5 -2 0.2 0.4

Approximation par 3 elements lineaires

CHAPITRE 2. GRANDES LIGNES DE LA METHODE DES ELEMENTS FINIS

Approximation par 3 elements lineaires,avec C.L. imposees exactes

0 -0.5 -1 -1.5 -2 -2.5

Fig. 2.9 { Approximation

par 3 elements d'ordre 2,avec C.L. imposees exactes

Ecole Superieure de Mecanique de Marseille

3.1 De nition du domaine

3.1.1 De nition du domaine dans le logiciel d'elements nis

Ecole Superieure de Mecanique de Marseille