Mathématiques des éléments finis
Mathématiques des éléments finis
Méthode des
éléments finis
Dans ce (de moins en moins court) document, plutôt à destination d’ingénieurs mécaniciens
connaissant déjà la méthode des éléments finis, nous allons essayer de faire une présentation un
peu plus théorique que ce qui leur est généralement proposé (et qui est quand même souvent de
type « preuve par les mains », ce qui occulte trop de points).
Nous ne ferons appel qu’à des notions mathématiques de bases généralement déjà vues pour la
plupart en taupe (ou en tout début de cycle d’ingénieur)... bien que des compléments que l’on peut
qualifier d’élémentaires nous aient été demandés et aient été inclus.
Nous espérons, grâce à cette présentation théorique montrer toute la souplesse et la puissance
de la méthode, afin de permettre au lecteur d’envisager d’autres simulations que celles qu’il a pu
déjà réaliser par le passé.
Pourquoi un ingénieur pratiquant déjà les éléments finis devrait-il s’intéresser plus en profondeur
aux mathématiques derrière la méthode ?
Tout d’abord parce que c’est beau et intéressant : deux raisons parfaitement licites et suffisantes.
Mais surtout parce que le monde change ! On souhaite des modélisations toujours plus proches
du réel, toujours plus détaillées, toujours plus complexes, toujours plus couplées... Par ailleurs un
constat s’impose : si la physique d’hier était essentiellement celle du continu, celle d’aujourd’hui
est le règne du discontinu. Ainsi, connaître l’intégrale de Riemann et savoir intégrer par parties
étaient autrefois des outils suffisants, alors qu’aujourd’hui il faut en passer par les dérivées au sens
des distributions, les espaces de Sobolev, l’intégrale de Lebesgue...
Rester sur les outils d’hier c’est se condamner à résoudre les problèmes d’hier ! C’est pourquoi
ce document a vu le jour, pour essayer de présenter et d’expliquer « simplement » (nous l’espérons)
les mathématiques derrière la méthode, mais sans demander au lecteur de se transformer en
mathématicien.
But du document
Le but initial était de présenter brièvement la théorie mathématique derrière les éléments finis afin
que les ingénieurs utilisant cette méthode puisse en envisager toutes les applications, ainsi que de
couvrir les aspects qui, selon nous, devraient être connus de tout ingénieur mécanicien impliqué ou
intéressé par le calcul numérique.
Toutefois, il s’envisage comme support de référence à plusieurs cours, cours qui ne portent pas
sur tous les aspects traités dans ce document, et pendant lesquels les aspects pratiques sont plus
développés (avec mise en situation sur machine).
Même si nous avons voulu rester le plus succinct possible, l’introduction de notions de proche
en proche à conduit à un document qui fait aujourd’hui une certaine taille (par exemple, nous avons
besoins des espaces de Sobolev, mais comment les introduire sans parler des espaces de Lebesgue,
mais comment les introduire sans parler...).
Aussi le document a-t-il finalement été découpé en plusieurs parties : un survol des notions
mathématiques, puis le traitement du problème continu constituent l’ossature théorique nécessaire
à assoir la méthode des éléments finis sur un socle solide. La discrétisation par éléments finis
à proprement parler n’est aborder qu’ensuite, et d’ailleurs un seul chapitre suffirait à en faire le
3 But du document
tour... sauf à entrer plus dans le détail concernant « ce qui fâche » : homogénéisation, non linéarité,
dynamique, ce qui est fait dans des chapitres séparés.
Enfin, d’autres méthodes sont abordées car également très employées aujourd’hui. Aussi est-il
indispensable selon nous d’en avoir entendu parlé et d’en connaître les principales notions (BEM,
FEEC...).
En annexes, se trouve un petit fourre-tout comprenant des choses censées être maîtrisées depuis
la taupe (mais qui parfois nous sont demandées) et les compléments qui alourdiraient encore les
propos précédents.
Certaines notions (essentiellement de topologie) ne sont pas présentées dans ce document. Il
nous a semblé que le lecteur devait avoir quelques souvenirs de ce qu’est un ouvert, un fermé,
l’adhérence, la densité... Par ailleurs, leur nom peut être suffisamment évocateur pour se passer
d’une définition formelle dans le contexte de ce document.
Attention, ce document n’est pas un document de mathématiques, il ne contient d’ailleurs
aucune preuve. C’est, dans ces deux premières parties, un document de vulgarisation de notions
mathématiques nécessaires à une bonne compréhension de la méthode des éléments finis.
Nous avons voulu réaliser un survol des notions importantes, mais malgré tout, afin de ne pas
être parfois trop laconique, nous avons un peu débordé.
En fin de document, un petit index des noms propres permettra au lecteur de replacer les divers
développements mentionnés dans l’histoire... Il se peut qu’il subsistent quelques erreurs, notamment
au niveau des nationalités mentionnées, car il n’est pas toujours aisé de déterminer rapidement cette
information (et nous ne connaissons pas toutes les biographies des personnes citées).
Ce document a été réalisé très rapidement, et de manière extrêmement hachée. Il comporte
forcément encore beaucoup de fautes : merci de m’en faire part.
Remerciements :
Nous remercions Mathias Legrand pour ses conseils avisés et sa relecture pertinente.
Notre collaboration a permis une très nette amélioration de la qualité typographique générale du
document, et a conduit à la coexistence de deux versions, issues du même code source : l’une que
nous appelons « version cours » correspond à ce que nous proposons en cours (couleurs, notations) ;
l’autre que nous nommons « version livre », plus classique et sage dans sa forme, est plus proche
d’un ouvrage.
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
But du document 3
Démarche de l’ingénieur numéricien 4
I SURVOL MATHÉMATIQUE
2 Applications et morphismes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.1 Fonction, application, injection, surjection, bijection 28
2.2 Morphismes 29
2.2.1 Présentation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.2.2 Cas des espaces vectoriels : application et forme linéaires . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.2.3 Endo, iso, auto -morphismes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.2.4 Espace dual d’un espace vectoriel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.2.5 Noyau et image . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
2.3 Opérateur 32
3 Continuité et dérivabilité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.1 Continuité et classe C0 33
3.2 Continuité de Hölder et Lipschitz 35
3.3 Dérivée 35
3.4 Fonctions de classe Ck 37
3.5 Différentielle 37
3.6 Dérivées partielles 37
3.7 Retour sur les classes Ck pour une fonction de plusieurs variables 39
4 Espaces de Lebesgue . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4.1 Présentation des espaces de Lebesgue Lp 48
4.2 Construction de Lp 49
4.3 Espace L0 50
4.4 Espace L1 et dualité avec L1 50
4.5 Espace L2 50
4.6 Compléments et retour sur les fonctions continues et différentiables 51
5 Espaces de Sobolev . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
5.1 Distributions 56
5.2 Dérivées au sens des distributions 58
5.3 Espaces Wm;p ./ 59
5.4 Espaces Hm ./, Hm
0 ./ et H
m ./ 60
5.5 Prise en compte du contour du domaine 61
5.5.1 Trace . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
5.5.2 Espace trace . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
5.6 Espaces H1 ./, H10 ./ et H 1 ./ 63
5.7 Espaces [Link]/ et [Link]/ 65
5.8 Inégalités utiles 66
II PROBLÈME CONTINU
14 Le maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 197
14.1 Maillage de Delaunay 198
14.1.1 Maillage simplexial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 198
14.1.2 Maillage de Delaunay-Voronoï . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 199
14.1.3 Remarques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 201
14.2 Maillage par avancement de fronts 202
14.3 Maillage par transformation 202
14.4 Remarques sur le maillage quadrangulaire et hexaédrique 204
15 Homogénéisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 205
15.1 Méthodes d’homogénéisation 206
15.1.1 Méthode de développement régulier . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 207
15.1.2 Méthode de la couche limite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 207
15.1.3 Méthode de développement asymptotique infini . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 208
15.1.4 Cas des coefficients discontinus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 208
18 L’acoustique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 239
18.1 Introduction à l’acoustique physique 240
18.1.1 Émission . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 240
18.1.2 Transmission . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 240
18.1.3 Réception . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 246
18.2 Calculs acoustiques par éléments finis 252
18.2.1 Modèles simplifiés pour valeurs de référence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 252
18.2.2 Constitution d’un modèle éléments finis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 254
18.2.3 Convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 258
18.2.4 Vers l’infini... . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 258
18.2.5 ... et au-delà . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 259
18.2.6 Post-Traitement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 259
18.3 Quelques illustrations avec F REE F EM ++ 260
18.3.1 Un exemple en acoustique des salles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 260
18.3.2 Un silencieux automobile . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 260
18.3.3 Deux mots de statistiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 261
18.3.4 Sur les conditions aux limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 262
18.3.5 Un peu d’amortissement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 264
18.3.6 Un obstacle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 265
18.3.7 Transmission entre deux milieux . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 266
IV ANNEXES
15 SURVOL MATHÉMATIQUE I
Chapitre 1
Résumé — Dans les problèmes que nous envisageons de traiter in fine, il n’est finalement
besoin que d’espaces vectoriels normés finis sur le corps des réels. On pourrait rapidement
en donner les principales propriétés qui sont sans doute encore en mémoire des lecteurs
tant elles sont « naturelles ».
Mais en faisant cela, nous ne respecterions pas nos engagements de présenter un peu
plus avant les fondements mathématiques.
Sans toutefois aller trop loin dans les notions topologiques générales, nous allons
essayer dans ce chapitre de passer en revue la liste des espaces depuis les espaces
topologiques jusqu’aux espaces de Hilbert, de manière succincte et didactique (nous
l’espérons) en relevant les points importants essentiellement du point de vue de leur
utilité.
Le mot « topologie » vient de la contraction des noms grecs « topos » (lieu) et « logos » (étude),
c’est donc l’étude du lieu. On a d’ailleurs commencé par l’appeler Analysis Situs, le terme « topolo-
gie » n’étant introduit qu’en 1847, en allemand, par Johann Benedict Listing dans Vorstudien zur
Topologie.
La topologie vise à définir ce qu’est un lieu (i.e. un espace) et quelles peuvent être ses propriétés
(je dirais uniquement en tant que tel, sans autre ajout). Elle s’intéresse plus précisément à ce que
l’on appelle aujourd’hui espaces topologiques et aux applications, dites continues, qui les lient, ainsi
qu’à leurs déformations (“A topologist is one who doesn’t know the difference between a doughnut
and a coffee cup”).
En analyse, elle fournit des informations sur l’espace
considéré permettant d’obtenir un certain nombre de résultats
(existence et/ou unicité de solutions d’équations au dérivées
partielles, notamment). Les espaces métriques ainsi que les
espaces vectoriels normés sont des exemples d’espaces topo-
logiques. L’origine de la topologie est l’étude de la géométrie
dans les cultures antiques. Le travail de Leonhard Euler da- Euler Poincaré Hausdorff
tant de 1736 sur le problème des sept ponts de Königsberg
est considéré comme l’un des premiers résultats de géométrie qui ne dépend d’aucune mesure, i.e.
l’un des premiers résultats topologiques.
Henri Poincaré publia Analysis Situs en 1895, introduisant les concepts d’homotopie et d’homo-
logie. Bien d’autres mathématiciens ont contribué au sujet parmi lesquels nous citerons : Cantor,
Volterra, Arzelà, Hadamard, Ascoli, Fréchet, Hausdorff...
Finalement, une dernière généralisation en 1922, par Kuratowski, donna le concept actuel
d’espace topologique.
Notons bien que rien n’est dit sur les éléments de l’ensemble E. Ce peuvent être des éléments
discrets, des scalaires, des vecteurs, des fonctions...
Unifiant les travaux de ses prédécesseurs sur les espaces de fonctions, c’est en 1906 que Maurice
Fréchet introduit le concept d’espace métrique.
La métrique qui nous est la plus usuelle est évidemment la métrique
euclidienne, qui est celle que nous utilisons en géométrique « classique »
(euclidienne) : la distance entre deux points est égale à la longueur du
segment les reliant. La structure métrique fournit beaucoup plus d’infor-
mation sur la forme géométrique des objets que la structure topologique.
Enfin nous redonnerons le si important critère de Cauchy (qui est va-
lable pour tout espace uniforme, dont notamment les espaces métriques) Fréchet Cauchy
qui permet de définir la toute aussi importante notion de complétude.
E, un ensemble
Une distance ou métrique d est une application de E E ! RC telle que :
— d.x; y/ D d.y; x/ (symétrie) ;
— d.x; y/ > 0 si x ¤ y, et d.x; x/ D 0 (positivité) ;
— d.x; z/ 6 d.x; y/ C d.y; z/ (inégalité triangulaire).
Espace métrique
On peut réexprimer les notions d’ouvert, fermé, adhérence... densité, continuité... avec la
métrique (les "...).
Deux normes sont équivalentes si elles définissent la même topologie.
Deux normes k k1 et k k2 sont équivalentes s’il existe deux constantes strictement positives k 0
et k 00 telles que 8x 2 E, kxk1 6 k 0 kxk1 et kxk2 6 k 00 kxk1 .
Critère de Cauchy :
Soit E un espace métrique et soit x0 , x1 , ..., xn , ... une suite d’éléments de E. Cette suite est
de Cauchy de E si :
8" > 0; 9N 2 N; .8n 2 N; 8m 2 N; n > N; m > N/ W [Link] ; xn / 6 " (1.1)
ou encore : [Link] ; xn / tend vers 0 quand m et n tendent vers l’infini.
Espace métrique complet
Toute suite convergente est de Cauchy.
Réciproque : si x0 , x1 , ..., xn , ... est de Cauchy sur R ou C, alors elle converge.
Cette propriété est fondamentale car elle permet de reconnaître si une suite converge sans
connaître sa limite.
Attention, la notion d’espace complet est une notion métrique et non topologique (elle peut donc
être vraie pour une métrique et fausse pour une autre).
L’importance de la complétude tient à ce que la solution d’un problème passe souvent par une
solution approchée. Si la suite des solutions approchées est de Cauchy d’un espace convenable,
et si cet espace est complet, alors la convergence vers une solution du problème est assurée.
E, un ensemble
Structure d’espace vectoriel :
C’est une structure comportant une loi de composition interne et une loi de composition
externe sur un corps K permettant d’effectuer des combinaisons linéaires (voir cours sur les
structures algébriques).
La loi de composition interne, notée C en fait un groupe abélien, la loi de composition
externe est la multiplication par un scalaire, scalaire pris sur le corps K considéré.
Espace vectoriel (topologique)
Sur un espace vectoriel de dimension finie sur R ou C, deux normes quelconques sont équiva-
lentes.
Une norme sur un espace vectoriel E est une fonction, x 7! kxk possédant les propriétés :
— positivité : kxk > 0 pour x ¤ 0, k0k D 0 ;
— transformation par les homothéties : kxk D jjkxk; 2 K ;
— inégalité de convexité : kx C yk 6 kxk C kyk.
La distance issue de la norme est la distance définie par d.x; y/ D kx yk.
Espace vectoriel normé
Tout espace vectoriel normé est automatiquement un espace métrique (avec la distance issue de
sa norme). De plus :
— la distance est invariante par translation : d.x a; y a/ D d.x; y/.
— une homothétie de rapport multiplie la distance par jj : d.x; y/ D jjd.x; y/.
Tout espace vectoriel normé de dimension finie est localement compact (i.e. tout point possède
au moins un voisinage compact), car toute boule fermée est compacte.
L’intérêt des distances issues d’une norme est qu’elles rendent continues les opérations de
l’espace vectoriel et qu’en particulier les translations sont des homéomorphismes.
Comme un espace vectoriel normé muni de la distance issue de sa norme est un espace mé-
trique, on peut se demander si cet espace métrique vérifie le critère de Cauchy (voir paragraphe
précédent) afin d’en faire un espace complet.
Espace de Banach
C’est donc finalement un espace vectoriel normé sur un sous-corps K de C (en général, K D R
ou C), complet pour la distance issue de sa norme.
Comme la topologie induite par sa distance est compatible avec sa structure d’espace vectoriel,
c’est un espace vectoriel topologique. p
Une norme k k découle d’un produit scalaire ou hermitien h; i si l’on a : kxk D hx; xi.
Espace de Hilbert
C’est donc un espace préhilbertien complet, i.e. un espace de Banach dont la norme découle
d’un produit scalaire ou hermitien unique.
Un espace de Hilbert est la généralisation en dimension quelconque d’un espace euclidien ou
hermitien (Un espace vectoriel euclidien est un espace vectoriel réel de dimension finie muni
d’un produit scalaire).
Théorème 1 Tout sous-espace vectoriel fermé d’un espace de Hilbert est un espace de Hilbert.
où xi désigne le conjugué de xi .
— L’espace des fonctions (dans R ou C) continues et bornées sur un espace topologique X,
muni de la norme kf k D supx2X .jf .x/j/. En particulier, l’espace des fonctions continues
sur un espace X compact, comme un intervalle réel Œa I b.
— Pour tout réel p > 1, l’espace Lp des classes de fonctions mesurables, dans R ou C, sur un
espace mesuré X.
Voici quelques espaces de Hilbert classiques :
— L’espace euclidien Rn muni du produit scalaire usuel est hilbertien,
— L’espace L2 des fonctions de carré intégrable (voir chapitre 4).
— Certains espaces de Sobolev (voir chapitre 5) sont des espaces de Hilbert (ceux qui nous
intéresserons en général, ça tombe bien).
Ce résultat subsiste si l’on suppose seulement que H est un espace préhilbertien et que F est un
sous-espace vectoriel complet de H.
Définition 1 — Algèbre. Une algèbre (de Boole) E sur un ensemble X est un ensemble non
vide de parties de X, stable par passage au complémentaire et par union finie (donc aussi par
intersection finie), i.e. :
1. E 6D ¿
2. 8A 2 E, c A 2 E, où c A désigne le complémentaire de A dans X. (donc X 2 A)
3. si A1 ; : : : ; An 2 E alors A1 [ ::: [ An 2 E
Définition 2 — Tribu. Une tribu A sur un ensemble X est un ensemble non vide de parties de X,
stable par passage au complémentaire et par union dénombrable (donc aussi par intersection
dénombrable), i.e. :
1. A 6D ¿
2. 8A 2 A, c A 2 A S
3. si 8n 2 N, An 2 A alors n2N An 2 A
Ces deux dernières définitions étant placées l’une sous l’autre, leur différence doit apparaître
clairement : dans le cas d’une algèbre, on a à faire à un nombre fini d’intersections ou de réunions,
dans celui d’une tribu, on prend en compte un ensemble dénombrable (donc fini ou non).
Il est donc évident que toute tribu est une algèbre, la réciproque étant fausse.
Définition 3 — Espace mesurable. Le couple .X; A/ est appelé espace mesurable ou espace
probabilisable en fonction du contexte. Sur les espaces mesurables on définit des mesures (voir
ci-après) ; sur les espaces probabilisables, on appelle ces mesures des probabilités.
Les parties de X qui appartiennent à la tribu A sont appelées ensembles mesurables. Dans
un contexte probabiliste, on les appelle événements (il suffit que .X/ D 1 pour que soit une
probabilité).
Une mesure sur un ensemble X est une fonction qui associe à chaque élément d’une tribu
d’un ensemble X un scalaire positif (une longueur, un volume, une probabilité...).
Définition 4 — Mesure. Soit .X; A/, un espace mesurable. Une application définie sur A, à
valeurs dans Œ0; C1 est appelée mesure lorsque les deux propriétés suivantes sont satisfaites :
1. L’ensemble vide a une mesure nulle : .¿/ D 0
2. L’application est -additive : si E1 , E2 , ... est une famille dénombrable de parties de X
appartenant à A et si ces parties sont deux à deux disjointes, alors la mesure .E/ de leur
réunion E est égale à la somme des mesures des parties :
1 1
!
[ X
Ek D .Ek /:1 (1.3)
kD1 kD1
On appelle espace mesuré un triplet .X; A; /, où X est un ensemble, A une tribu sur X et
une mesure sur A.
Complément sur le produit d’espaces mesurables. Soient .X; T/ et .Y; U/ deux espaces mesurables. On
appelle rectangle mesurable du produit D X Y, toute partie de de la forme A B où A et B sont
des éléments respectivement T et U-mesurables.
On appelle produit tensoriel des deux tribus T et U, la tribu engendrée par l’ensemble des rectangles
mesurables. Cette tribu est notée T ˝ U, et est la plus petite tribu de qui contient toutes les parties de la
forme A B, A 2 T, B 2 U.
Le produit des espaces mesurables est .X Y; T ˝ U/ et est un espace mesurable :
Z Z Z
f d. ˝ / D f .x; y/d.x/d.y/ (1.4)
On aura besoin de ces notions (dont la tribu borélienne) pour l’intégrale de Lebesgue et
par extension toute la théorie de l’intégration qui est à la base de l’analyse numérique (et oui,
formulation faible, quand tu nous tiens) ainsi notamment que pour le théorème de représentation de
Riesz fondamental en éléments finis.
Les masses de Dirac sont très importantes, notamment dans la pratique car elles permettent par
exemple de construire des mesures par approximations successives.
Définition 7 — Mesure complète. Soit .X; A; / un espace mesuré, on dit que est une
mesure complète lorsque tout ensemble négligeable pour appartient à la tribu A sur laquelle
est définie, i.e. :
La mesure de Lebesgue a permis de bâtir une théorie de l’intégration palliant les insuffisances
de l’intégrale de Riemann (il suffit de vouloir intégrer 1Q sur Œ0 I 1 avec Riemann pour être dans
l’impasse). Nous détaillerons cela un peu plus au chapitre 4 sur les espaces de Lebesgue.
Parmi les définitions de cette mesure, nous présentons la plus intuitive, celle qui consiste à
généraliser la notion de volume en gardant les mesures sur les pavés de Rn .
Définition 8 — Mesure de Lebesgue. Il existe une plus petite mesure définie sur une tribu de
parties de Rn qui soit complète et coïncide sur les pavés avec leur volume (i.e. avec le produit
des longueurs de leurs côtés).
Cette mesure est appelée la mesure de Lebesgue (notée n ) et sa tribu de définition la
tribu de Lebesgue (notée Ln et que nous ne définissons pas ici, et dont les éléments sont dits
Lebesgue-mesurables).
Cette restriction aux boréliens de la mesure de Lebesgue est parfois dénommée mesure de
Borel-Lebesgue.
Si pour la mesure de Lebesgue sur R, la notion de « presque partout » correspond bien à
l’intuition, ce n’est pas vrai en général. Par exemple, pour D ıa sur la tribu de l’ensemble des
parties de X, une propriété est vraie -pp simplement si elle est vérifiée en a.
Si la mesure est nulle, toute propriété est vérifiée pp (ainsi que sa négation).
Théorème 4 La mesure de Lebesgue est invariante sous toutes les isométries. Elle est en
particulier invariante sous les translations : c’est une mesure de Haar du groupe topologique Rn .
Oui, la mesure de Lebesgue peut être vue comme une mesure de Haar. Mais, historiquement, la
mesure de Haar est définie plus tard (on parle souvent de la mesure de Haar, alors que l’on devrait
parler d’une mesure de Haar). Elle généralise celle de Lebesgue.
Définition 9 — Mesure régulière. Une mesure (positive) définie sur une tribu contenant la
tribu borélienne d’un espace topologique X est dite régulière lorsque elle est à la fois intérieure-
ment régulière et extérieurement régulière :
1. 8A X de la tribu, .A/ D supf.K/ j K compact contenu dans Ag ;
2. 8A X de la tribu, .A/ D inff.O/ j O ouvert contenant Ag.
Théorème 5 La mesure de Lebesgue est finie sur tout compact ; chaque compact, qui est borné,
pouvant être enfermé dans un cube.
Elle est par voie de conséquence régulière, Rn étant métrisable, localement compact et sépa-
rable.
— Rn est métrisable : un espace topologique est un espace métrisable lorsqu’il est homéomorphe
à un espace métrique (un homéomorphisme est une application bijective continue entre deux
espaces topologiques dont la réciproque est continue) ;
— Rn est localement compact : un espace localement compact est un espace séparé qui, sans
être nécessairement compact lui-même, admet des voisinages compacts pour tous ses points
(De plus, on a la propriété : Tout espace localement compact est régulier) ;
— Rn est séparable : un espace séparable est un espace topologique contenant un sous-ensemble
fini ou dénombrable et dense, i.e. contenant un ensemble fini ou dénombrable de points dont
l’adhérence est égale à l’espace topologique tout entier.
Applications et morphismes
Résumé — Au chapitre précédent, des espaces ont été définis, mais on ne s’est pas
intéressé beaucoup aux relations entre eux ou au sein d’eux.
Des distances, normes... ont été introduites, sans utiliser plus que ça le vocable de
fonction. Le terme d’injection a été prononcé à la fin du chapitre précédent comme fil
conducteur pour introduire celui-ci...
Dans ce chapitre, nous ne présenterons pour le coup que des choses extrêmement
« rudimentaires » et toutes vues en taupe ou avant. Il s’agit uniquement d’un aide-mémoire.
L’univers mathématiques du début du XVIIIe siècle est dominé par Leonhard Euler et par ses apports
Histoire
tant sur les fonctions que sur la théorie des nombres, tandis que Joseph-Louis Lagrange éclairera la
seconde moitié de ce siècle.
Euler a introduit et popularisé plusieurs conventions de notation par
le biais de ses nombreux ouvrages largement diffusés. Plus particulière-
ment, il a introduit la notion de fonction (dans L’Introductio in analysin
infinitorum, premier traité dans lequel le concept de fonction est à la
base de la construction mathématique, et dont les premiers chapitres
lui sont consacrés) et a été le premier à écrire f .x/ pour désigner la
fonction f appliquée à l’argument x, en 1734 (bien que le terme de Euler Lagrange
« fonction » apparaisse pour la première fois dans un manuscrit d’août
1673 de Leibniz, resté inédit, et intitulé la Méthode inverse des tangentes ou à propos des fonctions).
Il a également introduit la notation moderne des fonctions trigonométriques, la lettre e pour la
base du logarithme naturel (également connue sous le nom de nombre d’Euler) en 1727, la lettre
grecque † pour désigner une somme en 1755 et la lettre i pour représenter l’unité imaginaire, en
1777. L’utilisation de la lettre grecque pour désigner le rapport de la circonférence d’un cercle à
son diamètre a également été popularisée par Euler, mais celui-ci n’est pas à l’origine de la notation.
Les différentes techniques mises au point (par exemple pour la résolution des équations différen-
tielles, le développement en séries entières ou asymptotiques, applications aux réels négatifs, aux
complexes...) conduisent à s’intéresser à la « fonction » en tant que sujet d’étude.
À la fin du XVIIIe siècle, les mathématiciens croient encore, pour
peu de temps, que la somme infinie de fonctions continues est continue,
et (pour plus longtemps) que toute fonction continue admet une dérivée...
(sur ces notions, voir chapitre suivant).
C’est Cauchy qui met un peu d’ordre dans tout cela en montrant que
la somme d’une série numérique n’est commutativement convergente que
Cauchy Abel
si la série est absolument convergente. Mais Cauchy, qui pourtant n’est
qu’à un doigt de la notion de convergence uniforme, énonce un théorème
faux de continuité d’une série de fonctions continues qu’Abel contredit par un contre-exemple
le 16 janvier 1826 : Cauchy affirme en 1821 que la somme d’une série de fonctions continues
est toujours continue. Cinq ans plus tard, Abel propose un contre-exemple en considérant la suite
.fn /n>1 d’applications continues de R dans R de terme général fn .x/ D . 1/n =n: [Link]/.
= élément n’ayant pas de relation ; = élément ayant 1 relation ; = élément ayant plus d’une
relation.
Tableau 2.1: Types de fonctions
Une étoile rouge dans E n’a pas de sens, cela voudrait dire qu’un élément de E peut avoir
plusieurs valeurs différentes par la relation considérée...
En fait, on sait donner un sens à cela. C’est ce que l’on appelle une fonction multivaluée
ou fonction multiforme ou fonction multivoque ou multifonction. L’exemple le plus simple d’un
fonction multiforme est la fonction réciproque d’une application non injective (penser simplement
aux fonctions circulaires).
On trouve les fonctions multiformes en analyse complexe : lorsque l’on veut utiliser le théorème
des résidus pour calculer une intégrale réelle, on peut être amené à considérer des restrictions
(déterminations) qui font de ces fonctions multiformes des fonctions (univoques), par exemple en
utilisant la théorie des revêtements qui considère des fonctions sur des surfaces de Riemann.
En restreignant une fonction à son domaine de définition, on en fait une application. En la
restreignant en plus à son ensemble d’arrivée on en fait une surjection (une surjection, c’est un
« truc » défini partout sur E et F).
Quand on a une surjection, on est sûr que tout élément de l’ensemble de départ à une image, et
2.2 Morphismes
Cette section est extraite du cours sur les structures algébriques. Nous l’avons toutefois sérieusement
amputée pour coller à l’objectif de ce document.
2.2.1 Présentation
Histoire
Mort au cours d’un duel à l’âge de vingt ans (ce qui en fait un héros romantique), il laisse un
manuscrit élaboré trois ans plus tôt, dans lequel il établit qu’une équation algébrique est résoluble
par radicaux si et seulement si le groupe de permutation de ses racines a une certaine structure,
qu’Emil Artin appellera justement résoluble.
Son Mémoire sur les conditions de résolubilité des équations par radicaux, publié
par Joseph Liouville quatorze ans après sa mort, a été considéré par ses successeurs,
en particulier Sophus Lie, comme le déclencheur du point de vue structural et
méthodologique des mathématiques modernes.
Toutefois, pour être tout à fait exact, Lagrange, reprenant une idée de d’Alembert
en vue de prouver qu’un polynôme de degré n possède n racines (i.e. que C est
algébriquement clos), utilise des résultats sur les fonctions semblables des racines, Galois
i.e. des fonctions rationnelles des racines qui restent invariantes par les mêmes
permutations. Ce résultat, qu’il a établi dans son fameux mémoire de 1771 Réflexions sur la
résolution algébrique, inspirera Abel et Galois et peut être considéré comme le tout premier de la
théorie des groupes.
Soient deux ensembles G et G0 munis d’un même type de structure (topologique, groupe,
anneau, espace vectoriel, algèbre...). Un morphisme (ou homomorphisme) de G ! G0 est une
application f qui respecte cette structure.
Pour ce faire, cette application doit vérifier certaines conditions, notamment une certaine
« linéarité » vis-à-vis des lois des G et G0 (on pourrait également remplacer le terme linéarité par
« capacité à faire sortir de la fonction »).
Un morphisme entre deux espaces topologiques est tout simplement une application continue
(voir chapitre suivant sur la continuité). C’est d’ailleurs ce dernier terme qui est utilisé en topologie,
pas celui de morphisme (mais cela revient bien au même).
Un morphisme de groupe entre .G; / et .G0 ; ?/ satisfait à l’égalité suivante qui est bien une
« condition de linéarité par rapport à la loi » : 8.x; y/ 2 G; f .x y/ D f .x/ ? f .y/. En particulier,
si e et e 0 sont les éléments neutres de G et G0 , alors : f .e/ D e 0 . Une autre conséquence directe est
que : 8x 2 G; f .x 1 / D Œf .x/ 1 .
Ceci est équivalent à la condition suivante (on parle de « préservation des combinaisons
linéaires ») :
La structure d’un espace et celle de son dual sont très liées. Nous allons détailler quelques
points en nous restreignant aux cas qui nous intéressent (cas réel, dimension finie).
Si l’espace E est de dimension finie n, alors l’espace dual E est isomorphe à E et est donc lui
aussi de dimension n. On a alors le théorème de la base duale (que je ne présente pas, car je n’ai
pas parlé de base... mais peut-être pourrons-nous nous passer de ces rappels dans ce document).
Pour x 2 E, on note h'; xi pour '.x/. Cette notation est appelée crochet de dualité.
Notons que si E est réflexif, i.e. si E00 D E alors la topologie faible et la topologie faible-* coïncident.
Les topologies faibles ne sont en général pas métrisables, et ne peuvent donc se définir par la seule
donnée des notions de convergence de suites. Cependant, c’est bien la notion de convergence faible de
suites qui est utile en pratique.
On place ici des remarques faisant appel aux espaces de Lebesgue définis au chapitre 4.
Soit un ouvert de Rn . Pour 1 < p < 1, les topologies faible et faible-* sur Lp ./ coïncident, et
la notion de convergence associée (convergence contre des fonctions test dans Lp ) est :
Z Z
8g 2 Lp 0 ./; fn g ! fg
Pour p D 1, la topologie faible-* correspond à la convergence contre des fonctions test dans L1 ;
pour p D 1, la topologie faible correspond à la convergence contre des fonctions test dans L1 . On
s’interdira en revanche de considérer la convergence faible-* dans L1 , ou la convergence faible dans L1
(la question de savoir si L1 est le dual de L1 touche à de subtiles questions d’axiomatique, et la réponse
est négative si l’on admet l’axiome du choix...)
Quel est l’intérêt d’appauvrir la topologie ? Une des motivations majeures est que moins il y a
d’ouverts, plus il y a de compacts. Il est beaucoup plus facile d’être compact pour la topologie faible que
pour la topologie forte.
Théorème 7 — Théorème du rang. Il est assez visible que (théorème de factorisation) f induit
un isomorphisme de l’espace vectoriel quotient E= ker.f / sur l’image im.f /. Deux espaces
isomorphes ayant même dimension, il s’en suit la relation, valable pour un espace E de dimension
finie, appelée théorème du rang :
On a également :
— l’image réciproque d’un sous-espace vectoriel de F par f est un sous-espace vectoriel de E ;
— l’image directe d’un sous-espace vectoriel de E par f est un sous-espace vectoriel de F.
2.3 Opérateur
Le terme « opérateur » a été utilisé au paragraphe précédent... regardons d’un peu plus près.
D’une manière générale, un opérateur est une application entre deux espaces vectoriels topolo-
giques.
Un opérateur O W E ! F est linéaire si et seulement si :
Continuité et dérivabilité
Dans le manuscrit de 1673 la Méthode inverse des tangentes ou à propos des fonctions, Leibniz
dit : « J’appelle fonctions toutes les portions des lignes droites qu’on fit en menant des droites
indéfinies qui répondent au point fixe et aux points de la courbe ; comme sont abscisse, ordonnée,
corde, tangente, perpendiculaire, sous-tangente, sous-perpendiculaire... et une infinité d’autres d’une
construction plus composée, qu’on ne peut figurer. » Finalement, au terme d’une correspondance
nourrie entre Leibniz et Jean Bernoulli, celui-ci donne en 1718 la définition suivante : « On appelle
fonction d’une grandeur variable une quantité composée, de quelque manière que ce soit, de cette
grandeur variable et des constantes. ». Il propose la notation x.
La continuité est en quelque sorte contenue, sous-jacente à ces définitions car les fonctions
considérées sont « physiques » et ne présentent au plus qu’un nombre fini de discontinuités.
Dans son Introductio in analysin infinitorum de 1748,
Euler définit une fonction d’une quantité variable comme
« une expression analytique composée d’une manière quel-
conque de cette quantité variable et de nombres ou de
quantités constantes ». Le mot « analytique » n’est pas
davantage précisé. En fait, pour Euler, une fonction est
une combinaison quelconque d’opérations prises dans le Euler Bolzano Heine
stock des opérations et des modes de calcul connus de
son temps, et applicables aux nombres : opérations classiques de l’algèbre, exponentielle, logarithme,
passage d’un arc à ses lignes trigonométriques..., certaines de ces opérations pouvant être itérées un
nombre illimité de fois...
Dans ce même ouvrage, Euler dit qu’une fonction est continue si elle est définie par une seule
expression anlytique (finie ou infinie) et mixte ou discontinue si elle possède plusieurs expression
analytiques suivant ses intervalles de définition.
La définition actuelle est celle due à Bernard Bolzano dans sa théorie des fonctions en 1816 : « La
fonction f .x/ varie suivant la loi de continuité pour la valeur x si la différence jf .x C w/ f .x/j
peut-être rendue plus petite que toute valeur donnée. » Il existe une notion de continuité uniforme
qui est plus forte que la simple continuité et fixée par Heinrich Eduard Heine en 1872.
Cette définition est donnée pour la culture, car nous n’avons pas rappelé la notion de voisinage dans
ce document. Cela n’a pas d’importance dans ce contexte puisque nous travaillerons sur des cas
moins généraux.
Évidemment, une application de E dans F est continue si elle est continue en tout point de E.
Ramenons nous à des choses plus connues et plus en lien avec ce document. Dans le cas des espaces
métriques, la continuité se définit comme suit.
Ainsi f est continue en a si et seulement si la limite de f en a existe (elle vaut alors nécessaire-
ment f .a/).
pas vers 0 lorsque .x; y/ tend vers l’origine. Une telle fonction est dite partiellement continue.
Pour une fonction réelle, on peut définir une fonction continue comme une fonction dont on
peut tracer le graphe sans lever le crayon. Si l’on exclut certains fonctions très particulières (comme
les fractales), alors l’idée générale de la continuité est bien traduite par cette phrase (la fonction ne
présente pas de « saut »).
La classe des fonctions continues est notée C0 . Elle inclut par exemple les fonctions continues
par morceaux ainsi que les constantes (dont la fonction nulle).
Attention à ne pas confondre la classe des fonctions continues C0 avec C0 l’ensemble des
fonctions continues qui s’annulent à l’infini (sous-espace de l’espace des fonctions continues).
La continuité höldérienne d’une fonction dépend donc d’un paramètre réel strictement positif a 2
0 I 1, et prend en compte toutes les variations de la valeur de la fonction sur son ensemble de
définition. Si 0 < a 6 1 est fixé, l’ensemble des fonctions réelles a-höldériennes est un espace
vectoriel, conventionnellement noté C 0;a .E; R/.
Théorème 8 Toute application f qui est a-höldérienne est continue. Mieux, elle est uniformé-
ment continue, dans le sens suivant :
Si " > 0, alors pour D ."=C/1=a ; d .x; y/ < ) d .f .x/; f .y// < " (3.5)
Le réel dépend de " mais est indépendant de la variable x parcourant l’espace de définition de
l’application.
Toute fonction lipschitzienne (en tant que fonction höldérienne) est uniformément continue. Toute
fonction réelle lipschitzienne est (absolument continue donc à variation bornée donc) dérivable
presque partout pour la mesure de Lebesgue et sa dérivée est essentiellement bornée.
3.3 Dérivée
Le nombre dérivé en un point d’une fonction à variable et valeurs réelles est le coefficient directeur
de la tangente au graphe de cette fonction en ce point, ou aussi le coefficient directeur de l’approxi-
mation affine de cette fonction en ce point : ce nombre n’est donc défini que si cette tangente, ou
cette approximation, existe. La dérivée d’une fonction f est une fonction qui, à tout nombre pour
lequel f admet un nombre dérivé, associe ce nombre dérivé.
Définition 22 — Dérivée d’une fonction. Soit f une application d’un ouvert du corps K
(resp. un intervalle de R) dans un espace affine normé F (resp. R), alors on peut donner un sens,
f .a C h/ f .a/
f 0 .a/ D lim 2F (3.6)
h¤0;h!0 h
aCh2
phénomènes physiques et techniques conduit à la création au siècle suivant de l’analyse en tant que
branche des mathématiques abordant les notions de dérivation, intégration et équations différentielles.
Les échanges, les conceptions et la compréhension des infinité-
simaux ont animé le monde scientifique pendant bien longtemps, la
discussion était tout autant philosophique que mathématique, ce qui est
somme toute assez normal compte tenu du sujet... Les fondateurs incon-
testés de l’analyse sont Newton et Leibniz. La portée de leurs travaux est
considérable car ils vont permettre non seulement la compréhension des
courbes (puis le calcul des aires), mais aussi celle du mouvement des Newton Leibniz
corps. C’est véritablement une révolution où l’on passe d’une science de
la statique à une science de la dynamique.
Un débat houleux sur la paternité du calcul différentiel entre mathématiciens britanniques et
allemands fit rage, mais Newton et Leibniz s’en tinrent à l’écart. Toutefois, cette polémique entre les
deux camps ne sera pas sans effet. Celle-ci fera que les anglais resteront à l’écart du développement
général des mathématiques au XVIIIe siècle, la tradition newtonienne dominante ayant abouti à une
certaine stagnation scientifique. Au début du XIXe siècle, avec la diffusion des notations symboliques
leibniziennes du calcul infinitésimal, un petit groupe de mathématicien de Cambridge se mettra à
réfléchir sur le rôle et l’importance de la symbolique.
L’opinion générale est aujourd’hui que Leibniz s’est presque certainement inspiré de certaines
des idées de Newton (dont les travaux sont antérieurs, mais non publiés au moment où Leibniz
publie), mais que sa contribution reste suffisamment importante pour que le mérite de l’invention du
calcul différentiel soit accordé aux deux hommes.
L’approche de Newton est proche de la physique : il parle en termes de mouvement physique et
ses notations ne sont employées que très rarement en dehors de la sphère des physiciens.
L’approche de Leibniz en revanche est plus géométrique et conduit à une présentation plus natu-
relle, encore utilisée aujourd’hui, et qui va rapidement être adoptée en Europe. C’est lui également
qui introduit la notation dy=dx encore utilisée.
Avant Newton et Leibniz, Descartes s’était intéressé au pro-
blème des tangentes à une courbe et Fermat avait introduit une
notion assez proche de la dérivée lors d’une recherche d’extremum.
En remontant encore plus dans l’Histoire, on peut dire que la mé-
thode d’exhaustion, telle qu’elle a été développée et utilisée par
Archimède (i.e. avec toute la rigueur nécessaire, en usant du procédé
Archimède Oresme
d’encadrement évitant le recourt aux ") est réellement la première
utilisation des infinitésimaux. Malgré, ou peut-être à cause de, la
grande originalité de ses travaux, Archimède n’a été que peu suivi dans le monde grec et n’a pas eu
de disciple direct. Les savants arabes commenceront d’ailleurs dès le IXe siècle à s’intéresser aux
procédés infinitésimaux mis en œuvre par le génial Alexandrin.
Notons également que le plus grand progrès théorique réalisé au Moyen-Âge est l’étude quanti-
tative de la variabilité par Nicole Oresme. Éveillé par la diffusion des méthodes infinitésimales des
Anciens, et d’Archimède en particulier, nourri par les spéculations scolastiques sur l’infini, le goût
des mathématiciens pour les considérations infinitésimales se développe peu à peu. Ils s’appliquent
d’abord à imiter le modèle archimédien, puis s’émancipent et essaient de trouver des substituts pour
la méthode d’exhaustion, jugée trop lourde.
En 1748, Euler définit les deux nouvelles opérations que sont la différentiation et l’intégration
dans Introductio in analysin infinitorum ; puis dans Calculi differentialis (1755) et Institutiones
calculi integralis (1770) il essaie de mettre au point les règles d’utilisation des infiniment petits et
développe des méthodes d’intégration et de résolution d’équations différentielles.
Enfin, dans la seconde moitié du XIXe siècle, on s’intéressera aux propriétés des fonctions
La dérivabilité est elle-aussi une notion topologique et non métrique (même si on sait l’écrire
en termes métrique comme ci-dessus), elle ne dépend donc pas de la norme choisie (du moment
que celle-ci est compatible avec la topologie choisie...).
3.5 Différentielle
Il s’agit de généraliser la formule (3.6) à des applications de Rn dans Rp . La variable étant mainte-
nant un vecteur de Rn , il n’est plus question de diviser par h... il faut donc modifier la définition
pour supprimer le dénominateur. Et cela est très facile : il suffit d’approcher l’accroissement de la
fonction par une application linéaire.
Cela signifie qu’au voisinage de a, l’application f se comporte à peu près comme l’application
affine x 7! f .a/ C Df .a/.x a/, pour laquelle on pourra utiliser les outils de l’algèbre linéaire.
Géométriquement, on retrouve l’idée qu’une courbe est, au voisinage d’un point, à peu près
confondues avec sa tangente 1 (droite tangente, plan tangent...). C’est pourquoi l’application
linéaire Df .a/ est également appelée application linéaire tangente à f en a.
Évidemment, dans le cas où E D R et F D R, tout ce qui a été dit coïncide avec ce que nous
avions vu avec la dérivée classique.
@f f .a1 ; : : : ; ai 1 ; ai C h; ai C1 ; : : : ; an / f .a1 ; : : : ; an /
.a/ D lim (3.9)
@xi h!0 h
Remarques. Pour en revenir au cadre général qu’est celui de la différentielle, on parle ici des dérivées
partielles pour désigner les dérivées dans la direction des vecteurs de base, ou plus exactement des dérivées
selon les vecteurs de base.
@f
Nous utiliserons, comme tout le monde, la notation @x , bien qu’elle présente quelques défauts : tout
i
d’abord elle est typographiquement lourde, ensuite elle a les apparences trompeuses d’un quotient, et
enfin la mention de la variable xi peut être source de confusions dans les calculs de dérivées de fonctions
composées. On devrait donc toujours préférer la notation @i f ou fi0 .
La formulation générale de la différentielle donnée au paragraphe 3.5 permet de retrouver tous les
cas :
— le cas « classique » d’une fonction d’une variable réelle à valeur réelle correspond au cas E D R
et F D R ;
— le cas d’une fonction numérique d’une fonction vectorielle correspond au cas E D Rn et F D R
et vient d’être présenté. Dans ce cas, la différentielle Df .a/ est la forme linéaire sur Rn de
composantes T @1 f.a/; :::; @n f.a/ dans la base canonique.
Si Rn est muni d’un produit scalaire, alors il existe un unique vecteur, appelé gradient de f en a et
noté grad f .a/ 2 Rn tel que : Df .a/h D grad f .a/:h. Pour une définition plus opérationnelle du
gradient, on se réfèrera au paragraphe 3.8. Celui-ci a pour composantes les dérivées partielles par
rapport à la base considérée ;
— le cas d’une fonction vectorielle d’une variable réelle correspond au cas E D R et F D Rp . Il est à
nouveau possible de diviser par h. La différentiabilité 2 de f équivaut alors à la dérivabilité de ses
composantes fi .
— le cas d’une fonction vectorielle d’une variable vectorielle correspond au cas le plus général E D Rn
et F D Rp . La différentielle est l’application linéaire définie dans les bases canoniques de Rn
et Rp par la matrice jacobienne. La ligne i de la matrice jacobienne correspond à la différentielle
de la composante fi de f . La colonne j de la matrice jacobienne correspond à la dérivée de f
dans la direction du j e vecteur de la base. Nous y reviendrons de manière plus pragmatique au
paragraphe 12.5.2.
@f @f
Attention : Même si toutes les dérivées partielles @x1
.a/; : : : ; @x n
.a/ existent en un point a, la
fonction peut ne pas être continue en ce point. Si l’on considère à nouveau la fonction f définie
sur R2 par :
8 xy
<
2 2
si .x; y/ ¤ .0; 0/
f .x; y/ D x C y (3.10)
0 si .x; y/ D .0; 0/
:
@f @f
on voit qu’elle vérifie @x
.0; 0/ D @y
.0; 0/ D 0 mais elle n’a pas de limite en .0; 0/ comme nous
l’avons vu avant.
2. On ne confondra pas « différencier » et « différentier »
@2 f @2 f
D (3.11)
@xi @xj @xj @xi
@f @f X @f
df D dx1 C C dxn D dxi (3.12)
@x1 @xn @xi
i
3.7 Retour sur les classes Ck pour une fonction de plusieurs variables
Ce paragraphe met l’accent sur le cas multi-dimensionnel des classes Ck (ce qui est le cas avant,
mais peut-être moins visiblement). Nous en profitons également pour introduire des notations et
notions dont nous nous servirons plus loin.
Compléments sur Ck . Ck ./ est l’espace vectoriel des fonctions f W ! R telles que 8˛, j˛j 6 k,
x 7! @˛ f .x/ existe et appartient à C0 ./.
Pour f et g dans Ck ./, on définit :
sup j@˛ f .x/ g.x/j
X
1
X 1 j˛j6k Ki
d.f; g/ D (3.15)
2i 1 C sup j@˛ f .x/ g.x/j
X
iD1
Ki
j˛j6k
qui est une distance sur Ck ./ qui en fait un espace complet.
L’espace C1 ./ peut se définir par :
C1 ./ D Ck ./
[
(3.16)
k2N
qui est également complet pour la même distance.
Cette définition permet de clairement voir les inclusions successives des espaces Ck .
qui est une norme sur Ckb ./ et en fait un espace de Banach (i.e. normé complet pour cette norme).
Définition 28 On définit, pour k 2 N [ f1g, l’ensemble Ckc ./ par : f 2 Ckc ./ si f 2 Ck ./
et si f est de support compact inclus dans .
La dérivée d’un champ scalaire est un champ vectoriel appelé gradient (voir plus bas).
Nabla, noté r, est un pseudo-vecteur servant à noter un opérateur différentiel :
@
0 1
@
0 1 0
@
1
B @x C B @ C
B 1@r@ C
B C
B @ C
B1 @ C
B C
rDB C DB DB (3.19)
B C B C
B @y C B r @ C
C C
B @' C
@ @ A @ @ A @ 1 @ A
@z cartésiennes @z polaires r sin @' sphériques
grad f rf (3.20)
— La divergence :
div A r A (3.21)
rot A r ^ A (3.22)
f D r 2 f (3.23)
2
A D r A (3.24)
0 @A
z @Ay 1
@y @z
rotationnel : rot A D @ @A
B x @Az C
— @z @x A
@Ay @Ax
@x @y
@2 f 2 2
— Laplacien scalaire : f D @x 2 C @@yf2 C @@zf2
0 1
Ax
— Laplacien vectoriel : A D Ay A
@
Az
De plus :
On appelle normale extérieure une normale qui pointe vers l’extérieur du domaine en tout point.
@u
.x/ D ru.x/ n.x/ (3.31)
@n
Il s’agit d’un produit scalaire car ru est un vecteur, tout comme n.x/.
AD grad U (3.32)
Le même résultat se vérifie dans l’espace entier à condition que, vers l’infini, le champ décroisse
« assez » rapidement.
Le potentiel scalaire n’est pas unique, il est défini à une constante près.
Un champ de vecteurs A régulier et de divergence nulle dérive d’un potentiel vectoriel B :
AD rot B (3.33)
Le même résultat se vérifie dans l’espace entier à condition que, vers l’infini, le champ décroisse
« assez » rapidement.
Le potentiel vecteur n’est pas unique, il est défini à un « gradient » près (i.e. B0 D B C grad f
convient également).
où est la frontière de , ^ est le produit vectoriel et n.x/ est la normale dirigé vers l’extérieur.
Une autre identité remarquable met en relation l’intégrale de surface du rotationnel d’un champ
vectoriel et l’intégrale curviligne (ou circulation) du même champ sur la frontière. Elle découle du
théorème de Green qui, pour une surface S (généralement non fermée) de frontière C, implique :
“ I
rot A ds D A dl (3.36)
S C
Si S est fermée, C est vide (ou réduit à un point) et le membre de droite est nul.
L’orientation de la surface et celle de la courbe frontière sont liées puisque le changement d’une
orientation modifie le signe de l’intégrale correspondante. En fait, la relation est satisfaite lorsque
Il existe une autre façon de noter ce théorème. On se place sur un domaine compact lisse du
plan , de bord @, en notant la forme différentielle !. Alors la dérivée extérieure de ! s’écrit :
@Q @P
d! D dx ^ dy (3.38)
@x @y
(Le cercle sur l’intégrale précise que le bord décrit une courbe fermée).
@u
Z Z Z
.u/v D ru rv C v (3.41)
@n
f r 2 g C rf rg D
f rg n.x/ (3.44)
Ces formules sont exploitées pour obtenir des formulations faibles associées à des problèmes aux
dérivées partielles.
Espaces de Lebesgue
Comme promis au paragraphe 1.3, nous fournissons maintenant quelques compléments sur l’intégra-
tion. Nous avons mentionné la mise en défaut de l’intégrale de Riemann dans le cas de l’indicatrice
des rationnels. Ce contre-exemple a été historiquement fourni par [Link] toutefois
d’un peu plus près d’où tout cela provient.
Histoire
En termes modernes, on dira que pour qu’une fonction bornée soit Riemann-intégrable, il faut et il
suffit que l’ensemble des points de discontinuité de f soit de mesure nulle.
Le cadre classique le plus simple pour définir une intégrale est celui des fonctions en escalier
sur un intervalle Œa; b, ou sa complétion pour la topologie de la convergence uniforme, l’espace
des fonctions réglées (admettant une limite finie à droite et à gauche). La théorie de Riemann
permet d’atteindre une plus grande généralité. Cependant, Riemann lui-même a conscience que
l’intégrabilité « au sens de Riemann » impose encore des conditions relativement fortes : il démontre
qu’une fonction f W Œa; b ! R est intégrable si et seulement si, pour tout ˛ > 0 donné, on peut
choisir une décomposition de Œa; b en sous-intervalles suffisamment fins pour que la somme des
longueurs des sous-intervalles sur lesquels l’oscillation de la fonction dépasse ˛ soit arbitrairement
petite.
Les conditions énoncées au paragraphe précédent
peuvent sembler faibles, mais, dès la fin de son mémoire
de 1829 Sur la convergence des séries trigonométriques
qui servent à représenter une fonction arbitraire entre des
limites données, Dirichlet donne l’exemple, d’une nature
toute nouvelle, discontinue en tous ces points : la fonc-
tion f .x/ qui vaut une constante c si x est un rationnel Riemann Dirichlet Lebesgue
et qui vaut une autre constante d si x est irrationnel (on
appelle une telle fonction, « fonction de Dirichlet », bien que celui-ci l’ai toujours exhibée comme un
monstre, existant, mais non représentatif de ce que devrait être une fonction. Lorsque c D 1 et d D 0,
on parle alors plutôt de l’indicatrice de Q). On voit alors directement que la Riemann-intégrabilité
est mise à mal. Le travail consistera alors à définir précisément ce que l’on peut négliger, i.e. les
parties de mesure nulle.
Dans la théorie de Lebesgue, on va élargir la classe des fonctions intégrables (par exemple,
toute fonction bornée « raisonnable », disons qui peut être décrite par un énoncé mathématique, est
Lebesgue-intégrable). Pour cela, Lebesgue va commencer par définir un ensemble mesurable : c’est
un ensemble dont la mesure extérieure (i.e. la borne inférieure de la mesure des ouverts le contenant)
est égale à sa mesure intérieure (i.e. la borne supérieure de la mesure des fermés qu’il contient) :
est-il besoin de rappeler la différence entre borne inférieure et minimum et entre borne supérieure
et maximum ? La borne supérieure (ou le supremum) d’une partie d’un ensemble partiellement
ordonné est le plus petit de ses majorants. Une telle borne n’existe pas toujours, mais si elle existe
alors elle est unique. Elle n’appartient pas nécessairement à la partie considérée. Dualement, la borne
inférieure (ou l’infimum) d’une partie est le plus grand de ses minorants...
L’intégrale de Lebesgue peut être définie de manière géométrique : pour une fonction positive f
définie sur Œa I b, elle est égale à la mesure de dimension deux de l’ensemble f.x; y/ 2 R2 =a 6
y 6 f .x/; a 6 x 6 bg quand cette mesure existe. D’une manière analytique, cette définition de
l’intégrale de Lebesgue devient :
I SURVOL MATHÉMATIQUE 46
Définition 34 — Intégrale de Lebesgue. soit f une fonction réelle, bornée, définie sur Œa I b.
Supposons m 6 f .x/ 6 M pour x 2 Œa I b. Pour tous ; tels que 6 , on définit :
V; D fx 2 Œa I b= 6 f .x/ 6 g. Si pour tous , , V; est mesurable, alors f est mesurable.
En considérant la partition m D 1 < 2 < ::: < n < nC1 D M et Vi D fx=i 6 f .x/ 6
iC1 g, on définit les sommes
n
X n
X
i [Link] / et i C1 [Link] / (4.3)
iD1 iD1
où m.V/ est la mesure de l’ensemble V. L’intégrale de Lebesgue est la limite commune de ces
deux sommes, dont on peut montrer qu’elle existe lorsque f est mesurable.
Un espace Lp est un espace vectoriel (des classes) de fonctions dont la puissance d’exposant p
est intégrable au sens de Lebesgue, où p est un nombre réel strictement positif. Le passage à la
limite de l’exposant aboutit à la construction des espaces L1 des fonctions bornées.
Théorème 17 Chaque espace Lp est un espace de Banach lorsque son exposant p est > 1.
Théorème 18 Lorsque 0 < p < 1, l’intégrale définit une quasi-norme qui en fait un espace
complet.
Définition 37 — Exposants conjugués. Soient p et q 2 Œ1; C1. On dit que p et q sont des
exposants conjugués si :
1 1
C D1 (4.5)
p q
Théorème 19 Soient p et q des exposants conjugués, alors il existe une dualité entre les espaces
d’exposants p et q.
4.2 Construction de Lp
On considère un ouvert de Rn . Les fonctions f seront considérées de dans R ou C.
p
On appelle Lp D LM .; C/ pour p < 1 l’espace des fonctions f mesurables, de ! C,
telles que .jf jp / < 1 où est une mesure sur .
Comme une fonction s’annulant presque partout est d’intégrale nulle, on peut définir l’es-
pace Lp .X; A; / comme le quotient de l’espace des fonctions p intégrables Lp .X; A; / par le
sous-espace vectoriel des fonctions presque nulles.
Ce quotient identifie donc les fonctions qui sont presque partout égales, autrement dit qui ne
diffèrent que sur un ensemble de mesure nulle.
Théorème 20 Si l’espace X est fini et est muni d’une mesure finie .X/ < 1, alors tous les
espaces Lp (resp. Lp ), pour 1 6 p 6 1 sont les mêmes.
4.3 Espace L0
L’espace L0 est l’ensemble des fonctions mesurables. L0 est l’espace obtenu en quotientant L0 par
les fonctions nulles.
Soit ' une fonction mesurable strictement positive et -intégrable, alors :
kf gk
Z
.f; g/ D 'd (4.10)
1 C kf gk
Rappel : la notion de convergence (de suite) est une propriété topologique et non métrique
(l’écriture avec les " n’est que la traduction métrique de l’écriture topologique avec les boules).
4.5 Espace L2
Par définition, si est un ouvert donné de Rn , L2 ./ est l’espace des fonctions (réelles ou
complexes) qui sont de carré intégrable au sens de l’intégrale de Lebesgue.
La formule des exposants conjugués conduit dans le cas p D 2, à q D 2, i.e. que L2 s’identifie
à son dual.
Si l’on reprend la remarque précédente sur la dualité, cela devient dans L2 (et dans tout espace
de Hilbert H) : en associant à tout v 2 H l’application 'v .u/ D hu; vi, on peut identifier l’espace
de Hilbert H a son dual H0 .
Comme mentionné, l’intérêt d’avoir un espace de Hilbert est de disposer d’un produit scalaire,
donc de pouvoir décomposer un vecteur. Ce que l’on retrouve dans le théorème 3, que nous
réécrivons ci-dessous, et vrai dans L2 en tant qu’espace de Hilbert :
On rappelle que :
1 1 1
C D (4.13)
p q r
Lq Lp (4.15)
Terminons ce cours chapitre en revenant un peu aux fonctions continues et Ck , dont nous avons
déjà parlé au paragraphe 3.7. Nous commençons d’ailleurs par rappeler certaines définitions déjà
données.
Définition 41 Soit un espace topologique. On définit Cb ./ comme l’espace des fonctions
continues bornées de dans R ; C0 ./ comme l’espace des fonctions continues sur , tendant
vers 0 à l’infini ; Cc ./ comme l’espace des fonctions continues à support compact dans . Ces
trois espaces, munis de la norme sup :
kf k1 D sup jf j (4.16)
sont des espaces vectoriels topologiques ; les espaces Cb ./ et C0 ./ sont des espaces de
Banach, alors que Cc ./ ne l’est pas en général.
On note parfois cet espace Ck ./ pour insister sur le fait que l’on prend le supremum
jusqu’au bord. Une notation plus appropriée serait sans doute Ckb ./.
Définition 43 Soit un ouvert de Rn . L’intersection des espaces Cc ./ et Ck ./, pour tout k,
est appelée espace des fonctions C1 à support compact ; on la note C1
c ./ ou D./. Nous en
reparlerons au chapitre suivant...
Puisque Cc ./ n’est déjà pas un espace de Banach, il semble inutile de tenter de normer D./.
On dispose des résultats de densité suivants :
— L’espace C0c ./ des fonctions continues à support compact est dense dans Lp ./ pour 1 6
p < 1 (mais n’est pas dense dans L1 ./)
Dans le cas p D 1, cela se traduit par : 8f 2 L1 ./; 8" > 0; 9g 2 C0c ./; tel que kf
gkL1 ./ 6 "
— L’espace C1 p
c ./ des fonctions infiniment dérivables à support compact est dense dans L ./
pour 1 6 p < 1 (mais n’est pas dense dans L ./)1
— L’espace C1 1
c ./ est dense dans le sous espace de L ./ des fonctions bornées qui tendent
vers 0 en l’infini.
— L’ensemble des fonctions continues à support compact de Lp .Rn / est dense dans Lp .Rn /,
pour p ¤ 1.
Les arguments de densité ont une importance pratique dans certaines démonstrations : si l’on
doit démontrer qu’une propriété est vérifiée sur Lp ./, alors on peut commencer par montrer que
cette propriété est vraie sur C1
c ./, avant de passer à la limite en utilisant l’argument de densité :
1 p
« Cc ./ est dense dans L ./ ».
L’espace C1 c ./ est également noté D./ et appelé espace des fonctions tests. Les distributions
sont définies comme les éléments de D0 ./, dual topologique de D./, muni d’une topologie
I SURVOL MATHÉMATIQUE 4.6 Compléments et retour sur les fonctions continues et différentiables 52
adéquate. Nous allons les aborder dès a présent, au chapitre 5.
Espaces de Sobolev
Résumé — Les espaces de Sobolev sont des espaces fonctionnels. Plus précisément, un
espace de Sobolev est un espace vectoriel de fonctions muni de la norme obtenue par la
combinaison de la norme Lp de la fonction elle-même ainsi que de ses dérivées jusqu’à
un certain ordre. Les dérivées sont comprises dans un sens faible, au sens des distributions
afin de rendre l’espace complet.
Les espaces de Sobolev sont donc des espaces de Banach. Intuitivement, un espace
de Sobolev est un espace de Banach ou un espace de Hilbert de fonctions pouvant être
dérivées suffisamment de fois (pour donner sens par exemple à une équation aux dérivées
partielles) et muni d’une norme qui mesure à la fois la taille et la régularité de la fonction.
Les espaces de Sobolev sont un outil très important et très adapté à l’étude des
équations aux dérivées partielles. En effet, les solutions d’équations aux dérivées partielles,
appartiennent plus naturellement à un espace de Sobolev qu’à un espace de fonctions
continues dont les dérivées sont comprises dans un sens classique (mais rien n’empêche
d’avoir de la chance).
Le XXe siècle avait commencé par la thèse de Lebesgue intégrale, longueur, aire, qui constitue
Histoire
vraiment le début de la théorie de la mesure. La théorie de Lebesgue mène à l’étude des espaces Lp ,
qui permettront, sur les traces de Hilbert, Riesz et Banach, l’étude des opérateurs différentiels.
Cauchy avait publié nombre d’applications
de sa théorie dans des recueils d’exercices, no-
tamment concernant l’évaluation d’intégrales
réelles, qu’il n’hésita pas à généraliser en ce
qu’on appelle aujourd’hui la « valeur principale
de Cauchy », un peu moins d’un siècle avant
que Hadamard en ait besoin dans sa résolution Lebesgue Hadamard Schwartz Sobolev
des équations aux dérivées partielles dans un
problème d’hydrodynamique par les « parties finies » et que Laurent Schwartz n’en vienne aux
distributions.
L’étude des conditions de régularité des solutions des équations aux dérivées partielles permet à
Sergueï Sobolev et ses continuateurs de définir ses espaces de fonctions et les théorèmes de trace en
fonction des propriétés géométriques du domaine.
Quelques mots s’imposent au sujet de la théorie des distributions. Parmi tous les travaux ayant
valus à leur auteur la médaille Fields, la théorie des distributions de Laurent Schwartz est l’une des
très rares à être abordable par les étudiants dès leur premier cycle universitaire. C’est une partie
de l’explication du réel engouement pour les distributions, l’autre partie étant sa puissance, sa
commodité d’usage et surtout sa très grande beauté.
Si cette théorie fournit aux analystes un cadre général et un formalisme agréable pour l’étude
des espaces fonctionnels et des équations aux dérivées partielles dans lesquels ils aiment à se placer,
elle n’est pourtant pas indispensable. D’une part, les spécialistes d’équations aux dérivées partielles
parviennent toujours à trouver des formulations bien adaptées à leurs problèmes, en s’inspirant des
idées sous-jacentes à la théorie des distributions, mais sans y avoir recours explicitement. D’autre
part, on peut étudier la plupart des espaces fonctionnels intéressants sans la notion de distribution :
5.1 Distributions
Une distribution (également appelée fonction généralisée) est un objet qui généralise la notion de
fonction et de mesure. La théorie des distributions étend la notion de dérivée à toutes les fonctions
localement intégrables et au-delà.
Histoire
Définition 45 — Distribution. Une distribution est une forme linéaire continue sur D./.
L’ensemble des distributions est donc le dual topologique de D./ et est par conséquent
noté D0 ./.
Remarque. L’espace des distributions est extrêmement grand, et contient tous les espaces fonctionnels que
nous avons mentionnés jusqu’à présent et leur dual, y compris les espaces à poids et les espaces locaux
(non présentés dans ce document).
Complément topologique. L’espace D n’est pas métrisable. Définir sa topologie n’est pas équivalent à
définir la convergence des suites. Cependant, dans la pratique c’est presque la seule notion de convergence
que l’on utilise.
Soit un ouvert de Rn . On muni D0 ./ de la topologie faible-* induite par D./. Alors :
0
– D ./ est un espace topologique complet (et un espace de Montel) ;
– Le dual topologique de D0 ./ s’identifie à D./, qui est donc réflexif ;
– La topologie de D0 ./ n’est pas métrisable.
Notons que la convergence dans le dual de n’importe quel espace de Sobolev implique la convergence
au sens des distributions.
Toujours aussi classiquement, si T est une distribution et ' une fonction test de D./ alors
on note T.'/ D hT; 'i. (où h; i désigne comme d’habitude le crochet de dualité). Dans D0 .R/,
l’application qui à ' associe '.0/ est une distribution et c’est la distribution de Dirac :
hı; 'i D '.0/ (5.2)
Une propriété fondamentales est que toute fonction localement intégrable f représente aussi une
distribution Tf définie par la forme intégrale suivante :
Z
hTf ; 'i D f .x/'.x/ dx; 8' 2 D./ (5.3)
R
Cette définition étend la notion classique de dérivée : chaque distribution devient indéfiniment
dérivable et l’on peut montrer les propriétés usuelles des dérivées.
Cette notion peut aussi se définir comme la dérivée du produit de dualité hS; 'i. On note 'h W
x ! '.x h/ la translatée de ' par la valeur h 2 R, alors, en utilisant la linéarité et la continuité
par rapport au deuxième terme :
dS hS; 'h i hS; 'i 'h ' d'
h ; 'i D lim D hS; lim iD hS; i (5.6)
dx h!0 h h!0 h dx
Exemples de calculs avec les distributions. Considérons la fonction f .x/ D jxj sur R. Alors f 0 est la
fonction « signe », f 00 D 2ı0 , f 000 vaut 2 fois l’application « évaluation de la dérivée en 0 ». Attention
donc à ne pas confondre ces dérivées au sens des distributions avec les dérivées presque partout, qui sont
nulles à partir du rang 2.
Soit f .x/ D x log jxj x dans R, alors on peut écrire, au sens des distributions, f 0 .x/ D log jxj,
00
f .x/ D v:p:.1=x/ (où v:p: est la valeur principale).
Soit .un / la suite de fonctions définie par un .x/ D .1=n/si [Link]/. Alors la famille .u0n / converge au
sens des distributions vers 0.
où ˛ est un multi-indice tel que 0 6 j˛j 6 m , D˛ u est une dérivée partielle de u au sens faible
(i.e. au sens des distributions) et Lp un espace de Lebesgue.
est une norme équivalente à la précédente. On utilisera donc indifféremment l’une de ces normes,
et on notera la norme employée k kWm;p ou plutôt k km;p .
Si p < 1, alors Wm;p ./ est identique à la fermeture de l’ensemble fu 2 Cm ./I kukm;p <
1g par rapport à la norme kkm;p où Cm ./ désigne l’espace de Hölder des fonctions de classe m
sur .
Remarque. Le théorème de Meyers-Serrin dit « H D W » (qui est le titre de leur article de 1964 [57])
montre l’équivalence de deux définitions des espaces de Sobolev Wm;p ./ D Hm;p ./, i.e. dans notre
cas, donne une définition équivalente, par complétion de l’espace vectoriel normé :
0 11=p
1 ˛ p
X
fu 2 C ./I kukHm;p < 1g avec kukHm;p WD @ kD ukLp A (5.10)
j˛j6m
où D˛ u est une dérivée partielle de u au sens classique (u 2 C1 ./).
Définition 48 — Espace Hm ./. Dans le cas p D 2, on note Hm ./ l’espace Wm;2 ./, défini
par la relation (5.7).
Les espaces de Sobolev Hm ont un intérêt particulier car il s’agit alors d’espaces de Hilbert.
Leur norme est induite par le produit intérieur suivant :
X
D˛ u; D˛ v
.u; v/m D (5.11)
06j˛j6m
où :
Z
.u; v/ D u.x/v.x/ dx (5.12)
est le produit intérieur dans L2 ./ ; produit scalaire dans le cas réel, hermitien dans le cas complexe.
Théorème 29 Si est Lipschitzien, alors l’ensemble C1 ./ des fonctions régulières jusqu’au
bord de est dense dans Hm ./.
Complément. De plus, dans le cas où la transformation de Fourier peut être définie dans L2 ./, l’es-
pace Hk ./ peut être défini de façon naturelle à partir de la transformée de Fourier (voir cours sur le
traitement du signal). Par exemple si D Rn , si b u est la transformée de Fourier de u :
Z
Hm Rn D fu 2 L2 Rn I juO ./ ˛ j2 d < 1
(5.13)
Rn
pour j˛j 6 mg, ou ce qui est équivalent :
Z m
Hm Rn D u 2 L2 Rn I juO ./ j2 1 C 2
d < 1 (5.14)
Rn
et que :
Z m
.u; v/m D O v./
u./ O 1 C 2 d (5.15)
Rn
est un produit hermitien équivalent à celui défini plus haut. Encore, si D 0 I 1Œ, on vérifie que :
C1
( )
m 2 2 4 2m 2
X
H .0 I 1Œ/ D u 2 L .0 I 1Œ/I .1 C n C n C C n /jb un j < 1 (5.16)
nD 1
Définition 49 — Espace Hm m m
0 ./. On définit H0 ./ comme l’adhérence dans H ./ de D./,
ensemble des fonctions de classe C1 et à support compact dans .
m;p
Plus exactement et de manière générale, W0 est l’adhérence de D./ dans l’espace Wm;p . La
définition ci-dessus correspond donc au cas p D 2.
Remarque. Il faut bien comprendre que cela signifie quelque chose de simple : c’est l’ensemble constitué
des fonctions qui sont à la fois D./ et Hm ./. Comme on s’intéresse à la frontière (puisque l’on prend
l’adhérence), et que est un ouvert, la valeur que prennent les fonctions de cet ensemble sur la frontière
est celle que prennent les fonctions de classe C1 et à support compact dans , soit 0 !
muni de la norme :
0 1 12
X
kf kH m D inf @ kf˛ kL22 A (5.20)
j˛j6m
l’infimum étant pris sur toutes les décompositions possibles de f sous la forme intervenant dans
la définition de H m .
Ainsi défini, l’espace H m ./ est un espace de Hilbert isomorphe au dual topologique de Hm
0 ./,
et le crochet de dualité s’écrit :
X Z
˛
hf; uiH m ;H0 D
m . 1/ f˛ @˛ udx; 8f 2 H m ./; 8u 2 Hm 0 ./ (5.21)
j˛6m
C1 ./ est donc l’espace des fonctions C1 sur , prolongeables par continuité sur @ et
dont le gradient est lui-aussi prolongeable par continuité. Il n’y a donc pas de problème pour
définir la trace de telles fonctions.
— On montre que, si est un ouvert borné de frontière @ « assez régulière », alors C1 ./ est
dense dans H1 ./.
— L’application linéaire continue, qui à toute fonction u de C1 ./ associe sa trace sur @, se
prolonge alors en une application linéaire continue de H1 ./ dans L2 .@/, notée
0 , qu’on
appelle application trace. On dit que
0 .u/ est la trace de u sur @.
Pour une fonction u de H1 ./ qui soit en même temps continue sur , on a évidemment
0 .u/ D
uj@ . C’est pourquoi on note souvent par abus simplement uj@ plutôt que
0 .u/.
De manière analogue, on peut définir
1 , l’application trace qui permet de prolonger la définition
usuelle de la dérivée normale sur @. Pour u 2 H2 ./, on a @i u 2 H1 ./, 8i D 1; : : : ; n et on
peut donc définir
0 .@i u/. La frontière @ étant « assez régulière » (par exemple, idéalement, de
classe C1 ), on peut définir la normale n D .n1 ; : : : ; nn /T en tout point de @. On pose alors :
n
X
1 .u/ D
0 .@i u/ni (5.23)
i D1
Cette application continue
1 de H2 ./ dans L2 .@/) permet donc bien de prolonger la définition
usuelle de la dérivée normale. Dans le cas où u est une fonction de H2 ./ qui soit en même temps
dans C1 ./, la dérivée normale au sens usuel de u existe, et
1 .u/ lui est évidemment égale. C’est
pourquoi on note souvent, par abus, @n u plutôt que
1 .u/.
Cas p D 2 et D Rn
Dans ce cas, l’espace de Sobolev Hs .Rn /, s > 0, peut être défini grâce à la transformée de Fourier :
Z
s n 2 n 2 s 2
H .R / D u 2 L .R /W .1 C jj / ju./j
O d < C1 : (5.24)
Rn
Cas p ¤ 2 et D 0 I 1Œ
On définit un opérateur Ds de dérivation d’ordre fractionnaire s par :
1
X
s
D uD .i n/sb
u.n/ei nt (5.27)
nD 1
n’irons guerre au delà de m D 2, et même le plus souvent nous nous contenterons de m D 1. Nous
donnons ici quelques compléments dans ce dernier cas.
En pratique, il est très important de savoir si les fonctions régulières sont denses dans l’espace
de Sobolev H1 ./. Cela justifie en partie la notion d’espace de Sobolev qui apparaît ainsi très
simplement comme l’ensemble des fonctions régulières complétées par les limites des suites de
fonctions régulières dans la norme de l’énergie. Cela permet de démontrer facilement de nombreuses
propriétés en les établissant d’abord sur les fonctions régulières puis en utilisant un argument de
densité.
Par définition de H10 ./, et en prenant en compte une remarque précédente :
On voit que sur cet espace, la condition de Dirichlet est satisfaite automatiquement sur tout le
pourtour D @. La frontière est généralement partitionnée en deux sous-frontière D et N sur
lesquelles on satisfait les conditions de Dirichlet et de Neumann respectivement : D D [ N
Grâce à H 1 ./, on pourrait définir une nouvelle notion de dérivation pour les fonctions de L2 ./,
plus faible encore que la dérivée faible. Devant l’afflux de notions de dérivations, rassurons le
lecteur en disant qu’elles sont toutes des avatars de la dérivation au sens des distributions (c’est
l’intérêt de la théorie des distribution que d’avoir unifié ces divers types de dérivation).
L’exemple le plus simple à retenir, et le plus utile pour la suite, est que tout élément f
de H 1 ./ s’écrit, au sens des distributions, sous la forme :
f D u C div G (5.33)
On se contente ici de donner les résultats concernant l’espace H1 ./ même si des résultats
similaires peuvent êtres démontrés pour les espaces Hm ./ ou les espaces W1;p ./, avec p ¤ 2
(voir théorème 32). Soit un ouvert borné Lipschitzien de Rn .
1
— si n D 1, on a une injection continue de H1 ./ dans l’espace de Hölder C0; 2 ./ ;
— si n D 2, on a une injection continue de H1 ./ dans l’espace Lp ./ pour tout p < 1 (et
donc pas dans L1 ./.
— si n > 3, on a une injection continue de H1 ./ dans l’espace Lp ./ avec p D n2n2 .
De plus les injection non critiques sont compactes.
Comme on va s’intéresser par la suite à la discrétisation de problème aux dérivées partielles, le
cadre de domaines bornés nous suffira.
Définition 51 — Espace [Link]; /. L’espace [Link]; / est un espace intermédiaire entre L2 ./
et H1 ./, qui est défini par :
En particulier pour p D 2 :
Cette inégalité est assez puissante puisqu’elle contient, dans son membre de gauche, toutes les
dérivées partielles de v, alors que son membre de droite ne fait intervenir que certaines combinaisons
linéaires des dérivées partielles.
L’inégalité inverse de l’inégalité de Korn étant évidente, on en déduit que les deux membres
définissent des normes équivalentes.
Histoire
En mécanique, l’inégalité de Korn dit que l’énergie élastique (qui est proportionnelle à la norme
du tenseur des déformations dans L2 ./) contrôle (i.e. est supérieure à) la norme du déplacement
dans H1 ./ à l’addition près de la norme du déplacement dans L2 ./. Ce dernier point permet
de prendre en compte les « mouvements de corps rigides », i.e. les déplacements u non nuls mais
d’énergie élastique nulle.
n
X
Le produit scalaire défini par hx; yi D xi yi a pour norme induite la norme k k2 .
i D1
Soit E un espace vectoriel et .xn /n une suite de E. .xn /n est une suite de Cauchy ssi 8" >
0; 9N; 8p > N; 8q > N W kxp xq k > ".
Toute suite convergente est de Cauchy. La réciproque est fausse.
Un espace vectoriel est complet si et seulement si toute suite de Cauchy y est convergente.
Un espace de Banach est un espace normé complet.
Un espace de Hilbert est un espace préhilbertien complet.
Un espace euclidien est un espace de Hilbert de dimension finie.
Espaces fonctionnels
Un espace fonctionnel est un espace vectoriel dont les éléments sont des fonctions.
Ces formes sont des normes sur Lp ./ et en font des espaces de Banach (i.e. complets).
Dans le cas particulier p D 2, on ontientq l’espace L2 ./ des fonctions de carré sommable pp.
R
sur . On remarque que le norme kukL2 D 2
u est induite par le produit scalaire .u; v/L2 D
2
R
uv, ce qui fait de l’espace L ./ un espace de Hilbert.
Dérivée généralisée
Les éléments des espaces Lp ne sont pas nécessairement des fonctions très régulières. Leurs
dérivées partielles ne sont donc pas forcément définies partout. C’est pourquoi on va étendre la
notion de dérivation. Le véritable outil à introduire pour cela est la notion de distribution. Une idée
simplifiée en est la notion de dérivée généralisée (certes plus limitée que les distributions, mais
permettant de sentir les aspects nécessaires pour aboutir aux formulations variationnelles).
Fonctions tests
On note D./ l’espace des fonctions de vers R, de classe C1 , et à support compact inclus
dans . D./ est parfois appelé espace des fonctions-tests.
Théorème : D./ D L2 ./
Dérivée généralisée
Soit u 2 C1 ./ et v 2 D./. Par intégration par parties (ou Green) on a l’égalité :
Z Z Z Z
@i u' D u@i ' C u' n D u@i '
@
(car ' est à support compact donc nulle sur @).
Le terme u@i ' a un sens dès que u 2 L2 ./, donc @i u' aussi, sans que u ait besoin
R R
Trace
Dans les problèmes physiques que nous rencontrerons (voir partie 2), nous devrons pouvoir imposer
des conditions aux limites. Ceci est vrai, que l’on s’intéresse aux formulations forte ou faible.
Pour cela, il faut que la valeur d’une fonction sur la frontière soit définie. C’est justement ce
que l’on appelle sa trace. La trace est le prolongement d’une fonction sur le bord de l’ouvert .
De manière analogue, il est possible de prolonger la définition de la dérivée normale sur le
contour de , ce qui permet de prendre en compte des conditions aux limites de type Neumann par
exemple.
Inégalité de Poincaré : Si est borné dans au moins une direction, alors il existe une
constante C./ telle que 8u 2 H10 ./ ; kuk0 6 C./juj1 . On en déduit que j:j1 est une norme
sur H10 ./, équivalente à la norme k k1 .
Corollaire : Le résultat précédent s’étend au cas où l’on a une condition de Dirichlet nulle
seulement sur une partie de @, si est connexe.
On suppose que est un ouvert borné connexe, de frontière C1 par morceaux. Soit V D
fv 2 H1 ./I v D 0 sur 0 g où 0 est une partie de @ de mesure non nulle. Alors il existe une
constante C./ telle que 8u 2 V ; kuk0;V 6 C./juj1;V , où k k0;V et j:j1;V sont les norme
et semi-norme induites sur V. On en déduit que j:j1;V est une norme sur V équivalente à la
norme k k1;V .
73 PROBLÈME CONTINU II
Chapitre 6
6.1 Introduction
Une équation différentielle est une relation entre une ou plusieurs fonctions inconnues et leurs
dérivées. L’ordre d’une équation différentielle correspond au degré maximal de dérivation auquel
l’une des fonctions inconnues a été soumise.
Pour les méthodes explicites de résolution des équations différentielles, on ira voir au chapitre C.
Histoire
Les problèmes posés ou menant à des équations différentielles sont aussi vieux que l’analyse elle-
même (XVIIe-XVIIIe, voir notes précédentes sur les fonctions et sur la continuité et la dérivabilité).
Avant même qu’on ait complètement élucidé la question des infiniment petits l’on se préoccupe déjà
de résoudre des problèmes de tangente, qui mènent invariablement à une équation différentielle. Dès
les débuts de la mécanique classique, dite newtonienne, on est confronté à l’intégration de systèmes
d’équations différentielles du second ordre (voir exemples ci-dessous).
On s’habitue progressivement à ce que l’inconnue d’une équation puisse être une fonction,
même si la fonction a encore à l’époque un statut flou.
On résout les équations différentielles
par la méthode des séries entières sans s’en-
combrer des notions de convergence, bien
que l’on entrevoit parfois certaines diffi-
cultés dues justement à ces problèmes de
convergence... et on en arrive presque aux
séries asymptotiques qui ne seront concep- Riccati Clairaut D’Alembert Euler
tualisées qu’au XIXe siècle.
Un autre problème reste d’écrire, à l’aide de fonctions simples, les solutions des équations
différentielles. La liste des équations différentielles que l’on sait résoudre de cette manière est plutôt
maigre : Leibniz et sa méthode de décomposition en éléments simples, les Bernouilli et les équations
différentielles linéaires du premier ordre, puis Riccati...
La découverte par Clairaut en 1734 de l’existence d’une solution singulière à l’équation y
xy 0 C f .y 0 / D 0 (à la famille de droites y D Cx C f .C/ qui est l’expression générale des courbes
intégrales, il faut adjoindre l’enveloppe de cette famille pour avoir toutes les solutions analytiques de
l’équation) relance la dynamique. D’Alembert, en 1748, trouve un second cas d’intégrale singulière.
Si l’on considère une fonction f , sa dérivée f 0 exprime sa variation : positive f est croissante
(et plus sa value est grande, plus la croissance est rapide), négative f est décroissante...
À partir de là on peut considérer une population de personnes. Le nombre de total de personnes à
un instant t est donné par f .t /. Plus la population est nombreuse, plus elle se reproduit : en d’autres
termes, la vitesse de croissance de la population est proportionnelle à la taille de la population. On
peut donc écrire :
f 0 .t / D k f .t / (6.1)
C’est une équation différentielle très simple, mais qui modélise le problème dit de dynamique de
population.
Un autre problème simple, toujours en dynamique des populations, est celui de deux populations
interdépendantes : les proies g.t / et le prédateurs m.t /. Les proies se reproduisent et sont mangées
par les prédateurs, alors que les prédateurs meurent sauf s’ils peuvent manger des proies.
Cela conduit au système de Lotka-Volterra :
(
g 0 .t / D A g.t / B g.t / m.t /
(6.2)
m0 .t / D C m.t / C D g.t / m.t /
où il est nécessaire de résoudre les équations simultanément (on dit que les équations sont couplées).
Pour revenir à des exemples plus connus des lecteurs, on n’oubliera pas que la relation fonda-
mentale de la dynamique de Newton :
xR D c xP !2x (6.4)
Concernant les équations aux dérivées partielles, elles sont initialement résolues par des méthodes
ad-hoc, le plus souvent géométriques. On n’écrit pas l’équation aux dérivées partielles linéaires. On
ne s’étonnera donc pas que la première équation aux dérivées partielles n’apparaisse qu’assez tard.
La théorie de l’intégration des équa-
tions aux dérivées partielles du premier
ordre ne date que de 1734. L’idée d’Eu-
ler est de ramener ladite intégration à celle
des équations différentielles ordinaires. Eu-
ler montre qu’une famille de fonctions dé-
pendant de deux paramètres vérifie une Cauchy Lipschitz Dirichlet Neumann
équation aux dérivées partielles du premier
ordre en éliminant ces paramètres entre les dérivées partielles. Inversement, une telle équation admet
une solution dépendant d’une fonction arbitraire. Il parvient ainsi à intégrer plusieurs équations de
ce type.
Lagrange, dans un mémoire de 1785, résume les connaissances de l’époque sur ces questions. Il
ne sait intégrer que onze types d’équations aux dérivées partielles du premier ordre.
La solution de l’intégration de ces équations allait venir d’un mathématicien peu connu, mort en
1784, Paul Charpit. Son mémoire, présenté en 1784 à l’académie des sciences n’a jamais été publié
et est resté longtemps une énigme. Une copie de ce mémoire a été trouvée en 1928.
@2 u @2 u @2 u @u @u
˛ 2
C ˇ C
2
Cı C C
u D f (6.6)
@x @x@y @y @x @y
où ˛, ˇ,
, ı, ,
sont des scalaires. Cette équation est dite :
— elliptique si ˇ 2 4˛
> 0 ;
— parabolique si ˇ 2 4˛
D 0 ;
— hyperbolique si ˇ 2 4˛
< 0.
On peut dire que :
— les problèmes elliptiques vont concerner les problèmes stationnaires tels que ceux de la
mécanique, la thermique, l’électrostatique, les membranes élastiques, l’écoulement potentiel ;
— les problèmes paraboliques vont concerner les modèles instationnaires tels que la diffusion
thermique (équation de la chaleur), chimique, neutronique, les fluide visqueux incompres-
sibles (on parlera d’irréversibilité, de décroissance, de principe du maximum, de propagation
à vitesse infinie) ;
— les problèmes hyperboliques vont concerner les modèles instationnaires tels que la propa-
gation des ondes, l’électromagnétisme, l’élastodynamique (on parlera de réversibilité, de
conservation de l’énergie, de propagation à vitesse finie).
Le déplacement des atomes, ions ou molécules dans un milieu, que celui-ci soit solide (cristallin ou
amorphe), liquide ou gazeux, est appelé de manière générale « migration ».
Au sens large la diffusion désigne des transferts obéissant aux lois de
Fick, i.e. dont la résultante macroscopique vérifie l’équation de diffusion.
La turbulence entraine ainsi une forte diffusion dans les fluides. La
diffusion moléculaire est la migration sous l’effet de l’agitation thermique,
à l’exception des autres phénomènes. Elle intervient par exemple dans des
procédés d’amélioration des caractéristiques mécaniques (traitements de
surface comme la nitruration ou cémentation), la résistance à la corrosion Brown Fick
et les procédés d’assemblage par brasage.
En 1827, le botaniste Robert Brown observe le mouvement erratique de petites particules de
pollen immergées dans de l’eau. Il ne s’agit pas d’un phénomène de diffusion, puisque ce qui bouge
est une particule macroscopique, mais cette « marche aléatoire », que l’on appellera « mouvement
brownien », servira de modèle pour la diffusion : Le mouvement brownien, ou processus de Wiener,
est le mouvement aléatoire d’une « grosse » particule immergée dans un fluide et qui n’est soumise à
aucune autre interaction que des chocs avec les « petites » molécules du fluide environnant. Il en
résulte un mouvement très irrégulier de la grosse particule. Ce mouvement, qui permet de décrire
avec succès le comportement thermodynamique des gaz (théorie cinétique des gaz), est aussi très
utilisé dans des modèles de mathématiques financières.
En 1855, Adolph Fick propose des lois phénoménologiques, empiriques, inspirées de la loi
de Fourier pour la chaleur (établie en 1822). C’est Albert Einstein qui démontrera les lois de Fick
en 1905 avec ses travaux sur la loi stochastique. La première loi de Fick énonce que « le flux de
diffusion est proportionnel au gradient de concentration ».
' D 0 (6.7)
Elle est elliptique. Elle apparaît dans de nombreux problèmes physiques : astronomie, électrosta-
tique, mécanique des fluides, propagation de la chaleur, diffusion, mouvement brownien, mécanique
quantique. Les fonctions solutions de l’équation de Laplace sont appelées les fonctions harmoniques.
Toute fonction holomorphe est harmonique.
L’équation de Poisson est l’équation :
' D f (6.8)
où f est une fonction donnée. Elle est elliptique. Par exemple, en électrostatique, on exprime le
potentiel électrique V associé à une distribution connue de charges (dans le vide) par la relation :
V D (6.9)
0
En gravitation universelle, le potentiel gravitationnel ˆ est relié à la masse volumique par la
relation :
ˆ D 4 G (6.10)
En dépit des railleries de Hugo et Stendhal Joseph Fourier laissera son nom à la postérité à plus d’un
titre :
en tant qu’égyptologue, une spécialité qu’il bouleverse à la suite de sa
participation à l’expédition napoléonienne de 1798 (en 1810, il créera
l’Université Royale de Grenoble, dont il deviendra le recteur, et y re-
marquera Jean-François Champollion. Ils sont enterrés au cimetière
du Père-Lachaise à côté l’un de l’autre) ; en tant qu’homme politique
puisqu’il est préfet d’Isère sous Napoléon Bonaparte et sous la Restau-
ration ; en tant qu’administrateur, comme réorganisateur des statistiques Fourier Champollion
françaises ; et enfin en tant que scientifique, la consécration de sa car-
rière étant son élection au poste de secrétaire perpétuel de l’Académie des Sciences. Pour tous les
mathématiciens, il reste le fondateur de l’analyse de Fourier, qu’il crée au début du XIXe siècle pour
étudier mathématiquement la répartition de chaleur dans un corps conducteur.
Quelques repères chronologiques : En 1807, Fourier écrit l’équation de la chaleur. En 1808,
son mémoire est accepté par l’Académie des Sciences en dépit de sérieuses critiques concernant la
rigueur des démonstrations (Fourier a en particulier une controverse avec Poisson). En 1822 a lieu la
publication du Traité analytique de la chaleur, version remaniée et augmentée de son mémoire, qui
s’imposera comme l’un des ouvrages scientifiques majeurs du dix-neuvième siècle.
uP u D f (6.17)
Elle est [Link] que le problème soit bien posé, il faut spécifier :
— une condition initiale :
— de Neumann :
@T.x; t /
8x 2 @; D n.x/ rT.x; t / D 0 (6.20)
@n
où n.x/ est le vecteur normal unitaire au point x.
L’équation de la chaleur, introduite initialement pour décrire la conduction thermique, permet de
décrire le phénomène de diffusion. La diffusion est un phénomène de transport irréversible qui se
traduit à terme par une homogénéisation de la grandeur considérée (par exemple la température
dans un domaine, la concentration en produits chimiques dans une solution...). D’un point de vue
phénoménologique, et au premier ordre, ce phénomène est régi par une loi de Fick (par exemple
sorption d’eau dans les matériaux composites, diffusion d’actifs au travers de la peau). Dans
l’équation précédente, T représente alors la répartition de la grandeur considérée (eau dans un
composite, concentration d’un constituant chimique...) et le terme source dans est souvent nul
(i.e. f D 0).
Histoire
Des 23 problèmes de Hilbert énoncés en 1900, la moitié est aujourd’hui complètement résolue, l’un
a été démontré indécidable, cinq ne le sont pas du tout, les autres étant partiellement traités.
Des sept problèmes de l’institut Clay
(ou problèmes du prix du millénaire, do-
tés chacun d’un million de dollars améri-
cains) énoncés en 2000, seule la conjecture
de Poincaré a été démontrée en 2003 par
Grigori Perelman (qui a refusé le prix, ainsi
que la médaille Field. D’ailleurs il n’a pas Navier Stokes Perelman Ricci
publié son travail dans un journal, mais l’a
rendu publique et gratuit sur internet...). Notons au passage que Grigori Perelman est analyste, i.e.
un spécialiste des équations aux dérivées partielles, un sujet que beaucoup de topologues avaient
l’habitude de regarder de loin. Ils ont changé d’avis depuis !
u C rp D f (6.23)
On rappelle que si de plus le fluide est incompressible, alors on a également div u D 0 dans .
Contrairement à l’équation de Navier-Stokes, l’équation de Stokes est linéaire. Les écoulements
solutions de cette équation possèdent par conséquent des propriétés bien particulières :
unicité : pour des conditions aux limites données (valeur de la vitesse au niveau des parois et/ou à
l’infini), il existe un et un seul écoulement vérifiant l’équation de Stokes ;
additivité : les solutions de l’équation de Stokes vérifient le principe de superposition : si u1 et u2
sont solutions, alors toute combinaison linéaire 1 u1 C 2 u2 le sera aussi (ceci n’est pas
incompatible avec la propriété d’unicité : seul l’écoulement vérifiant les bonnes conditions
aux limites sera observé) ;
réversibilité : si un champ de vitesse u est solution de l’équation, alors u l’est aussi, à condition
de changer le signe des gradients de pression, ainsi que des vitesses aux parois et à l’infini ;
cette propriété est contraire à l’intuition, fondée sur notre expérience des écoulements
macroscopiques : la réversibilité des écoulements à bas nombre de Reynolds a ainsi poussé
les êtres vivants de très petite taille à développer des moyens de propulsion originaux.
paradoxe de Stokes : il faut prendre garde au fait que les solutions mathématiques de l’équation
de Stokes, dans un cas donné ou dans certaines régions du domaine de solution, peuvent
être physiquement fausses. Ceci est dû au « paradoxe de Stokes » à savoir que les conditions
physiques permettant de ramener l’équation de Navier-Stokes à l’équation de Stokes ne sont
pas nécessairement réalisées dans tout le domaine de solution, à priori. On aboutit alors à
des solutions présentant des comportements potentiellement aberrants dans certaines limites.
C’est le cas par exemple « à l’infini » où souvent le terme inertiel finit par l’emporter sur le
terme visqueux, sans qu’on puisse le préjuger à priori.
uP C rp D f (6.24)
Bien que ces équations correspondent à une simplification de l’équation de Navier-Stokes, elles
n’en sont pas plus stables... bien au contraire.
Si l’on considère le cas scalaire, i.e. celui d’une fonction inconnue u définie sur un domaine
ouvert R2 , une équations aux dérivées partielles elliptique linéaire se met sous la forme :
a.x; y/uxx .x; y/Cb.x; y/uyy .x; y/Cc.x; y/ux .x; y/Cd.x; y/uy .x; y/Ceu.x; y/ D f .x; y/
(6.26)
n
X n
X
Lu D aij .x/uxi xj .x/ C bi .x/uxi .x/ C c.x/u.x/ (6.28)
i D1;j D1 i D1
avec les valeurs propres de la matrice A D aij qui sont toutes non nulles et de même signe. On
rappelle que le matrice A peut toujours être choisie symétrique du fait de la symétrie uxi xj D uxj xi
et donc que les valeurs sont réelles.
D’ailleurs, en dimension finie, le fait que toutes les valeurs propres sont strictement positives
est équivalent au fait que la matrice A est définie positive.
II PROBLÈME CONTINU 6.6 Équations de la mécanique des milieux continus des solides 86
Flux conservatif
Les équations aux dérivées partielles elliptiques conduisent à la notion physique de flux conservatif
donné par un gradient. Considérons le cas scalaire. On associe à la fonction u.x/ un flux q.x/ qui
est nul vers l’extérieur :
Z
q.x/ nx ds D 0; 8 volume ! (6.29)
ı!
R
En utilisant la formule de divergence-flux, il vient : ! divx q.x/dx D 0, et en supposant la
fonction q.x/ suffisamment régulière on en déduit divx q.x/ D 0, 8x 2 .
En supposant que ce flux est une fonction linéaire du gradient ru, et qu’il est orienté dans
la direction opposée (physiquement, les flux se font souvent de façon opposée au gradient d’une
grandeur), on écrit q.x/ D a.x/ru et on obtient une loi de conservation à l’équilibre du type :
div.a.x/ru.x// D 0 (6.30)
Pour a.x/ 1, on retrouve la plus simple des équations elliptiques, l’équation de Laplace :
Le raisonnement précédent peut être mené en considérant les lois de conservations dans les milieux
continus. Les ingrédients principaux qui vont nous conduire à une grande famille d’équations
elliptiques sont les suivants :
— Conservation d’une grandeur ;
— Milieu immobile et régime stationnaire ;
— Loi de comportement linéaire entre le flux de cette grandeur et le gradient.
Les deux premières notions sont souvent associées à l’équilibre d’un système.
@
.u/ C div.uv/ D 'u div jU (6.33)
@t
Cette équation de conservation local peut, en utilisant la conservation de la masse, se mettre sous la
forme :
Du
D 'u div jU (6.34)
Dt
Du @u @u
D C ru.x; t / v.x; t / D (6.35)
Dt @t @t
@u
D 'u div jU (6.36)
@t
@u
Si de plus, on est en régime stationnaire, i.e. D 0, alors :
@t
En reportant cette loi dans l’équation de conservation, il vient div A.x; t /ru.x; t / D 'u, soit :
Cette équation est de type elliptique dès que l’opérateur A.x/ possède de bonnes propriétés de
positivité.
Équilibre général
j D A.x/ru (6.40)
où :
— la fonction scalaire inconnue u W 2 Rn ! R est appelée un potentiel ;
— la fonction vectorielle inconnue j W 2 Rn ! Rn est appelée un flux ;
— la fonction scalaire connue f W 2 Rn ! R correspond au termes sources du potentiel ;
— les fonctions scalaires connues A.x/ et c.x/ sont les données du problème.
Dans le cas vectoriel :
— u devient une fonction vectorielle 2 Rn ! Rn ;
— j un tenseur d’ordre 2 ;
— A.x/ un tenseur d’ordre 4 ;
— et c.x/ une fonction vectorielle.
II PROBLÈME CONTINU 6.6 Équations de la mécanique des milieux continus des solides 88
Retour sur les exemples précédents
Si u ! T est la température, j ! q le flux de chaleur, f la source de chaleur, A.x/ ! .x/
la conductivité du matériaux et c.x/ 0, on retrouve l’équation de la chaleur.
Si u ! V est le potentiel électrostatique, ru ! E le champ électrostatique, j ! D le
déplacement électrique ou flux de densité de courant, f ! la densité de charge, A.x/ ! .x/
le tenseur diélectrique et c.x/ 0, on retrouve la loi de Faraday.
Si u ! V est le potentiel électrique, ru ! E le champ électrique, j ! J la densité de
courant, f 0, A.x/ ! .x/ le tenseur de conductivité et c.x/ 0, on retrouve la loi d’Ohm.
Si u ! C est la concentration, j le flux molaire, f ! le taux d’absorption ou le taux
de réaction, A.x/ ! D.x/ la constante de diffusion et c.x/ la coefficient de réaction, alors on
retrouve le cas de la diffusion moléculaire.
Si u est le déplacement orthogonal, ru ! " la déformation, j ! la contrainte dans le
plan, f la force volumique extérieure, A.x/ ! H.x/ le tenseur des rigidités et c.x/ D 0, on se
retrouve dans le cas de l’élasticité plane (n D 2) linéaire.
limite élastique
composite
polymère
ε
F IGURE 6.1: Quelques lois types de comportement pour différents matériaux.
Dans le cas de matériaux isotropes (i.e. possédant les mêmes propriétés matérielles quelque
soit l’orientation), les relations constitutives s’écrivent à l’aide des coefficients de Lamé et :
II PROBLÈME CONTINU 6.6 Équations de la mécanique des milieux continus des solides 90
Si l’on considère un matériau constitué d’alvéoles (de
type hexagonal), mais dont les pointes sont rentrantes et
non sortantes (figure 6.3), alors le coefficient de Poisson
macroscopique de ce matériau est négatif : quand on le
comprime, il ne gonfle pas mais rétrécit.
On parle alors d’un matériau auxétique.
Matériau conventionnel Matériau auxétique
On voit bien que dans un tel cas, la différence de
comportement n’est due qu’à l’effet structurel, qu’à l’as- F IGURE 6.3: Matériau auxétique
semblage. Le matériau constituant les alvéoles lui, peut
(et est) tout à fait conventionnel : on trouve ce genre de
structures en carton, en aluminium et en acier, qui sont tous des matériaux homogènes isotropes.
Par exemple, la rigidité globale d’un assem-
Fibres
blage de couches de composites renforcés de fibres,
comme montré à la figure 6.4 peut conduire, selon
les orientations des plis, à des valeurs du coeffi-
cient de Poisson très variables, d’autant plus que
celles-ci sont fonctions du repère dans lequel elles
Plis
sont données (repère global, repère des fibres...)
Couche
Un stratifié unidirectionnel T300 fibre de car-
bone/résine époxyde 934 a un coefficient de Pois-
F IGURE 6.4: Matériau composite
son de seulement 12 D 0; 04.
Un stratifié unidirectionnel T300 fibre de car-
bone/résine N5208 a un coefficient de Poisson de 12 D 0; 3 dans le sens des fibres, mais de
12 D 0; 77 à 45˚.
En utilisant du verre (E D 69 GPa, D 0; 25) et de la résine époxyde (E D 3; 5 GPa, D 0; 4),
on peut réaliser un stratifié unidirectionnel verre-époxyde ayant pour propriétés El l D 45 GPa,
E t t D 12 GPa, Glt D 4; 5 GPa, lt D 0; 3 pour une fraction volumique de verre de 60%, et un
stratifié [0/90] ayant pour propriétés Ex D Ey D 20 GPa, Gxy D 2:85 GPa et xy D 0:13 pour
une fraction volumique de verre de 50%.
Notons également que si ¤ 0 et 3 C 2 ¤ 0, alors on peut inverser la relation contrainte/dé-
placement :
1
"ij D ij tr. /ıij (6.48)
2 3 C 2
Nous reviendrons plus tard sur ces aspects non linéaires au chapitre 19.
Une solution u du problème précédent est alors également solution du problème suivant
appelé formulation faible (On appelle v une fonction test) :
Z Z
Trouver une fonction u; définie sur ; vérifiant A.u/v D f v 8v définie sur :
(7.2)
Justification. Tout réside dans le fait de « multiplier par une fonction test et d’intégrer ». Pourquoi peut-on
faire ça ?
Partons de l’équation forte :
A.u/ D f; 8x 2 (7.3)
Nous pouvons multiplier les deux membres de l’équation par une fonction v à condition que cette fonction v
ne soit pas identiquement nulle. Dans ce cas, la solution de
A.u/v D f v; 8x 2 (7.4)
est la même 8v suffisamment sympathique. Comme cette seconde équation est vraie en tout point de ,
alors elle est encore vraie sous forme intégrale, i.e. sous forme faible, à condition que v soit suffisamment
régulière pour pouvoir être intégrée.
On entrevoit alors clairement la suite des opérations
R : si en plus d’être régulière, la fonction v est
différentiable, alors on va pouvoir réaliser sur le terme A.u/v une intégration par parties pour diminuer
les conditions de dérivabilité portant sur u, mais en augmentant celles portant sur v.
Si l’ordre de dérivation de u est paire (disons 2k), alors on pourra faire en sorte par un nombre
suffisant de manipulations que u et v aient à être différentiables à l’ordre k... et on pourra essayer de faire
en sorte que u et v appartiennent au même espace fonctionnel...
Remarque. Petit aparté sur d’autres manières de présenter « la chose » : en mécanique, on parle également
du « principe du travail virtuel » (ou des puissances virtuelles), de « minimisation de l’énergie potentielle »,
ou encore de « méthode des résidus pondérés », Toutes ces méthodes sont « équivalentes » à celle exposée
ici, mais souvent leur démonstration se fait uniquement « par les mains », ce qui peut être frustrant.
Principe des puissances virtuelles Le Principe des travaux virtuels est un principe fondamental
en mécanique, qui postule un équilibre de puissance dans un mouvement virtuel. Il s’agit d’une
formulation duale du principe fondamental de la dynamique. On raisonne comme suit : si un solide
est à l’équilibre (statique du solide), la somme des efforts est nulle. Donc si l’on fait faire un
déplacement fictif (virtuel) à l’objet, la somme des puissances des forces et moments est nulle. Ce
principe constitue la base d’une démarche de modélisation pour les milieux continus (théorie du
premier gradient, théorie du second gradient). On parle parfois du principe des travaux virtuels
qui est sensiblement identique.
Histoire
L’origine de ce principe revient à Jean Bernoulli, qui énonce en 1725 le principe des vitesses
virtuelles, qui consiste à considérer la perturbation de l’équilibre d’un système mécanique par un
mouvement infinitésimal respectant les conditions de liaison du système, un mouvement virtuel, et
d’en déduire une égalité de puissance.
Minimisation de l’énergie potentielle L’énergie potentielle d’un système physique est l’énergie
liée à une interaction, qui a le « potentiel » de se transformer en énergie cinétique.
Cette énergie est une fonction de ce système, dépendant des coordonnées d’espace, et éventuel-
lement du temps, ayant la dimension d’une énergie et qui est associée à une force dite conservative
dont l’expression s’en déduit par dérivation (ce qui veut dire au passage que dans les cas où
toutes les forces ne sont pas conservatives, par exemple s’il y a du frottement, alors la méthode ne
fonctionne pas et doit être « adaptée »). La différence entre les énergies potentielles associées à
deux points de l’espace est égale à l’opposé du travail de la force concernée pour aller d’un point
à l’autre, et ce quel que soit le chemin utilisé.
Méthode des résidus pondérés Si l’on considère l’équation A.u/ D f; 8x 2 sur , on ap-
pelle résidu des équations d’équilibre la quantité R.x/ D A.u/ f . Aux bords, on a les conditions
B.u/ D h; 8x 2 D @, et on appelle résidu des équations de bords la quantité R.x/ D B.u/ h.
En écrivant le problème dans et sur son contour de manière intégrale, on retombe sur une
formulation faible.
La méthode des résidus pondérés consiste à calculer les composantes de la solution approchée
par la projection sur des couples de fonctions de projection, appelées fonctions de pondération.
Selon les fonctions de pondération on retrouve les méthodes de collocation par sous-domaine
(fonctions constante sur chaque sous-domaine), de collocation par point (fonctions non nulles en
un seul point), des fonctions splines (en ajoutant des conditions des raccordements des fonctions et
des dérivées), de Galerkine (mêmes fonctions pour les fonctions de forme et de pondération).
Selon la nature du problème (donc selon la nature de l’opérateur A), la formulation faible (7.2)
peut être transformée (par exemple par intégration par parties), ce qui permet de faire apparaître une
forme symétrique ayant la nature d’un produit scalaire (ou d’un produit hermitien dans le cas de
fonctions complexes). Les deux fonctions u et v appartiennent alors à un même espace fonctionnel.
La formulation variationnelle d’un problème régi par des équations aux dérivées partielles
correspond à une formulation faible de ces équations qui s’exprime en termes d’algèbre linéaire
dans le cadre d’un espace de Hilbert.
Il n’y a donc fondamentalement pas de différence entre formulation faible et formulation
variationnelle, sauf que cette dernière appellation est généralement réservée à des cas « qui vont
bien »... et de manière très pragmatique aux cas où la formulation faible peut être transformée en
un problème de minimisation : on dit alors qu’il existe une fonctionnelle … dont la variation ı…
correspond au problème faible, d’où le terme de formulation variationnelle (d’un point de vue
physique, la fonctionnelle représente l’énergie du système, et le lieu où sa variation est nulle
ı… D 0 l’état recherché).
Par ailleurs, d’autres contraintes sur la frontière de peuvent (doivent) être imposées à u (et
à v). Ce sont les conditions aux limites. Les conditions aux limites types ont été présentées au
paragraphe 6.2.
Définition 56 — Coercitivité d’une forme linéaire. Une forme bilinéaire a est dite coercitive
si elle vérifie :
Autre démonstration. Une autre voie consiste à étudier directement le problème de minimisation de
1
la fonctionnelle définie par J D 2 a.v; v/ f .v/ (ce qui est une façon de montrer le théorème de
Riesz) en considérant une suite minimisante et en montrant qu’elle est de Cauchy par l’identité du
parallélogramme.
Ce théorème est l’un des fondements de la méthode des éléments finis dits « classiques » ou
« en déplacement ». Dans le cas d’éléments finis plus « exotiques » tels que les éléments mixtes,
hybrides, mixtes-hybrides et des multiplicateurs de Lagrange, on a besoin d’un peu plus... c’est ce
que nous allons voir aux paragraphes suivants.
La forme variationnelle du théorème de Lax-Milgram a déjà été abordée dans la seconde
« preuve » proposée du théorème de de Lax-Milgram. Il s’agit, dans les cas où cela est possible,
d’aborder le problème sous l’angle de la minimisation d’une fonctionnelle.
Pour les problèmes issus de la mécanique des solides déformables, la quantité 12 a.v; v/ repré-
sente l’énergie de déformation du solide associé à un « déplacement virtuel » v, et f .v/ l’énergie
potentielle des forces extérieures appliquées au solide pour le même déplacement virtuel v.
Histoire
Nous venons de voir le théorème de Lax-Milgram, fondement de la méthode des éléments finis,
publié dans leur article Parabolic equations de 1954 (Contributions to the theory of partial differential
equations, Annals of Mathematics Studies, no. 33, Princeton University Press, Princeton, N. J., 1954, pp.
167–190).
Plusieurs améliorations ont permis
d’étendre encore la portée de celle-ci. En
1971, dans son article Error-bounds for
finite element method (Numerische Mathematik
16 : 322–333), Babuška a généralisé le
théorème de Lax-Milgram aux problèmes
non nécessairement coercitifs (i.e. la forme Lax Babuška Brezzi Ladyjenskaïa
bilinéaire n’est plus supposée coercitive),
élargissant le champ d’application de la méthode des éléments finis.
Puis, en 1974, dans On the existence, uniqueness and approximation of saddle-point problems
arising from lagrangian multipliers, Brezzi a traité le cas où le problème n’est plus décrit à l’aide
d’une seule variable (u jusqu’à présent), mais à l’aide de plusieurs variables, ouvrant ainsi la porte
aux formulations mixtes. Ce faisant, il a eu besoin d’imposer des conditions, connues sous le nom
de conditions inf-sup, ou conditions BBL pour Babuška-Brezzi-Ladyjenskaïa (voir par exemple : O.
[Link]ženskaya, On integral estimates, convergence, approximate methods, and solution in functionals for
elliptic operators. (Russian) Vestnik Leningrad. Univ. 13 1958 no. 7, 60–69).
Enfin, Brezzi a remplacé sa condition de V-ellipticité par des conditions de Babuška dans son
théorème généralisé.
Nous allons écrire la démarche de la démonstration car elle est typique de la résolution de ce
genre de problème.
Démonstration. On va se servir de H D V ˚ V? avec V D ker B où évidemment B est l’opérateur définit
de H dans M0 tel que : hBu; iM D b.u; /, 8 2 M. Pour le produit scalaire induit, V et V? sont des
Hilbert.
La solution cherchée est u D u0 C u? .
Les inconnues seront donc u0 , u? et .
On a les équations pour v0 , v? et . Pour v? l’équation est :
< a.u? ; v? / C a.u0 ; v? / C b.v? ; / D hf; v? i 8v? 2 V?
8
Or dans la seconde équation, b.v0 ; / D 0 car v0 2 ker b et dans l’équation 3, b.u0 ; / D 0 car u0 2 V.
On est ramené à 3 sous-problèmes :
— sous-problème 1 : chercher u? 2 V? tel que b.u? ; / D g./, 8 2 M, ce qui est effectivement
le cas d’après le théorème de Babuška ;
— sous-problème 2 : chercher u0 2 V tel que a.u0 ; v0 / D hf; v0 i a.u? ; v0 /, 8v0 2 V. Il suffit de
vérifier que le second membre définit un opérateur tel que l’on a bien les hypothèses du théorème
de Lax-Milgram ;
— sous-problème 3 : chercher 2 M tel que b.v? ; / D hf; v? i a.u; v? /, 8v? 2 V? . Il faut un
peu plus travailler pour vérifier que l’on est dans le cas d’application du théorème de Babuška, mais
ça se fait (une quinzaine de lignes en détail).
On peut remplacer les conditions de V-ellipticité sur a.; / par des conditions de Babuška.
a.u; v/
9˛ > 0 tel que : sup > ˛kukH2 ; 8u 2 V D ker B (7.21)
v2Vnf0g kvk
Pour vérifier la condition LBB, il est souvent plus simple de vérifier le critère suivant : la condition
inf-sup sur b.; / est équivalente à :
1
8 2 M; 9v 2 H (uniquement dans V? ) tel que : b.v; / D kk2 et kvk 6 kk (7.23)
ˇ
et que l’on pose les nouvelles variables U D .u; / et V D .v; /, alors on revient à un problème
de type Lax-Milgram : trouver U 2 H M tel que pour tout V 2 H M, A.V; U/ D L.V/.
Résumé — Nous allons maintenant reprendre les problèmes types exposés sous forme
forte dans le chapitre 6, pour leur appliquer les méthodes du chapitre précédent. Nous
obtiendrons alors les formulations faibles et variationnelles des mêmes problèmes.
Laplacien de Dirichlet
On considère l’équation de Poisson (Laplace si f D 0) avec les conditions aux limites de Dirichlet,
i.e. le problème suivant :
(
u D f dans
(8.1)
uD0 sur D @
Multiplier par une fonction test v et intégrer par parties (ou utiliser la formule de Green) conduit à :
@u
Z Z Z
ru rv vD fv (8.2)
@n
Comme, rappelons le, le but est de faire en sorte que u et v jouent un rôle symétrique, on veut donc
qu’ils aient même régularité et mêmes conditions aux limites. La régularité nécessaire est u et v
sont H1 ./, et la prise en compte des conditions aux limites conduit finalement à chercher u et v
dans H10 ./.
Le problème est donc :
Z Z
Trouver u 2 H10 ./ tel que 8v 2 H10 ./; ru rv D fv (8.3)
dont le second membre n’a de sens que si f 2 L2 ./, ce qui sera supposé. Notons que le problème
est également bien défini dans le cas où f 2 H 1 ./, car v 2 H10 ./.
Avec cette formulation, le théorème de Lax-Milgram permet de conclure à l’existence et à
l’unicité de la solution.
Comme nous l’avons vu au chapitre précédent, l’intérêt des problèmes symétriques, c’est qu’ils
peuvent s’interpréter en terme de minimisation.
La solution du problème initial est également la fonction u qui réalise le minimum de :
1
Z Z
2
min jruj fu (8.5)
u2H10 ./ 2
Laplacien de Neumann
On considère l’équation de Laplace-Poisson avec les conditions aux limites de Neumann, i.e. le
problème suivant :
8
< u D f dans
(8.8)
: @u D g sur D @
@n
On obtient immédiatement la formulation variationnelle suivante :
Z Z Z
Trouver u 2 H1 ./ tel que 8v 2 H1 ./; ru rv D fv C gv (8.9)
sous la forme :
Z Z Z
Trouver u 2 V tel que 8v 2 V; ru rv D fv C gv (8.11)
En utilisant l’inégalité de Poincaré-Wirtinger, on montre que le problème est bien coercitif et
donc qu’il admet une unique solution.
Remarquons que si la condition de compatibilité n’est pas vérifiée, alors la solution de la
dernière formulation résout bien le problème initial avec un f ou un g modifié d’une constante
additive, de manière à vérifier la condition de compatibilité.
Il y aurait lieu à quelques discussions sur l’interprétation des résultats, mais celle-ci serait
beaucoup plus abstraite que dans le cas du Laplacien de Dirichlet. Le lecteur désireux de rentrer
(à raison) dans ce genre de détails trouvera aisément un cours de M2 sur le sujet. Cela dépasse le
cadre que nous nous sommes fixé dans ce document.
On opère comme précédemment. Toutefois, nous faisons jouer à la variable de temps t et la variable
d’espace x des rôles différents : on considère la fonction u.t; x/ comme une fonction du temps,
à valeur dans un espace fonctionnel en x, et on prend comme fonctions tests des fonctions qui
dépendent seulement de la variable x. On obtient alors facilement la formulation variationnelle
suivante :
Z Z Z
1 d
Trouver u tel que 8v 2 H0 ./; uv C ru rv D fv (8.17)
dt
Pour donner un sens à cette formulation variationnelle, une certaine régularité de u est requise :
Pour résoudre le problème, il faut utiliser une méthode numérique. On montre l’existence puis
l’unicité par les estimateurs à priori. Nous détaillerons cela plus loin.
et
b W H10 ./n L2 ./ ! R
(8.23)
Z
.u; q/ 7! b.u; q/ D .div v/q
On introduit également l’espace L20 ./ des fonction L2 ./ à moyenne nulle :
Z
L20 ./ D q 2 L2 ./; qD0 (8.24)
La formulation variationnelle mixte (vitesse/pression) est alors : Trouver .u; p/ 2 H10 ./ L20 ./
tels que :
Notons que dans la seconde équation, la fonction test q pourrait de manière équivalente être
cherchée dans L2 ./ (pas besoin de la moyenne nulle).
Par contre, il est important de chercher p 2 L20 ./. C’est cette condition de pression nulle qui
assure l’unicité de la pression (comme dans le cas précédent du Laplacien de Neumann).
Il est possible de formuler ce problème de manière˚non mixte, i.e. sous une forme compatible
avec le théorème de Lax-Milgram. Considérons V D u 2 H10 ./n = div u D 0 dans . Alors u
ˆ
ˆ b.u.t /; q/ D 0 8q 2 Y (8.31)
:
u.0/ D u0
où c W X X X ! R est la forme trilinéaire définie par :
n Z
@zi
Z X
c.w; z; v/ D Œ.w r/z v D wj vi (8.32)
@x j
i;j D1
Cette forme trilinéaire a certaines propriétés, par exemple lorsque l’on permute des termes où
lorsque l’on considère que certaines variables sont à divergence nulle...
Pour des raisons de stabilité, la forme c est remplacée par le forme cQ définie par :
1
cQ D .c.w; z; v/ c.w; v; z// (8.33)
2
Q
qui est antisymétrique et possède par conséquent la propriété c.w; z; z/ D 0, 8w; x; z 2 X.
Nous ne rentrons pas plus avant dans la méthode, car il faudrait alors parler dès à présent du
schéma numérique associé, ce qui sera fait plus tard.
< @u C .u r/u C rp D f
8
dans RC
@t (8.34)
div u D 0 dans RC
:
On considérera qu’il s’agit d’un cas simplifié de Navier-Stokes, et on le traitera comme tel.
uDu
8
ˆ sur D déplacement imposé
ˆ
ˆ
< n.x/ D g sur N condition sur le bord
N
(8.36)
ˆ
ˆ
ˆ P D 0/ D uP 0
u.t condition initiale en vitesse
:
u.t D 0/ D u0 condition initiale en déplacement
Nous sommes restés très prudents sur l’appartenance des différentes grandeurs à des espaces, car
cela dépend du choix des variables.
Notons que pour prouver la coercitivité de la forme bilinéaire, il faut recourir à l’inégalité de Korn.
Si en plus le matériau est isotrope, cette expression se simplifie en :
Z Z Z
div.u/ div.v/ C 2".u/ ".v/ D fv C gN v; 8v 2 HD1 ./ (8.42)
N
II PROBLÈME CONTINU 8.3 Équations de la mécanique des milieux continus des solides 112
D’une manière générale, toute formulation faible, quelque soit le nombre de champs inconnus,
s’écrit en mécanique sous la forme de minimisation d’une fonctionnelle, dont on sait trouver la
signification physique.
et ne comporte aucune condition (puisqu’elles sont toutes formellement présente dans la formula-
tion).
Il est possible d’obtenir une formulation n’ayant que le champ de contraintes comme inconnue.
Il s’agit de la fonctionnelle duale :
1
Z Z Z
J./ D S .v/S0 C u .v/ n (8.45)
2 D
1
Z Z Z Z
J.u; / D S .div C f /u C u n . n gN /v (8.47)
2 D N
En notant :
Z
a.; / D S (8.48)
Z Z Z
b.; u/ D div u C u n C u n (8.49)
D N
Z Z
L.u/ D fuC gN v (8.50)
N
Sa variation :
Z Z Z
ıJ.u; "/ D S C u .v/ n C v .u/ n gN v (8.54)
N
p C k 2 p D 0 (8.56)
Nous aurons l’occasion de reparler de cela au chapitre 18 qui portera spécifiquement sur l’acous-
tique.
Q D b.w/;
a.p; w/ Q 8w 2 HD0 ./ (8.58)
avec Q la conjugaison complexe, a.p; q/ W HD1 HD1 ! C la forme sesquilinéaire définie par :
Z Z
a.p; q/ D @xi p@xi qQ k 2 p qQ C ickAn qQ (8.59)
R
Et comme on peut vérifier que l’on est dans un cas qui va bien, on définit l’énergie potentielle
totale du domaine par :
1
Wp D a.p; p/
Q Q
b.p/ (8.61)
2
et la résolution du problème devient :
ou de manière équivalente :
Trouver p 2 HD1 ./ tel que ıWp D 0; 8ıp 2 H10 ./ (8.63)
En fonction des conditions aux limites et des conditions initiales du problème, la varia-
tion de pression solution du problème peut être réelle, i.e. être de la forme p 0 .x; y; z; t / D
<.p.x; y; z/e i !t /. Une telle solution est appelée monochromatique. On pourra alors résoudre
le problème en restant dans R au lieu de passer dans C. Cela se produit lorsque l’on résout
l’équation d’Helmholtz (8.56) assortie uniquement de conditions aux limites de Dirichlet.
Résumé — Dans cet exemple, portant sur formulation d’un problème de Neumann, nous
allons essayer de montrer comment la réflexion mathématique se fait et évolue « au fil
de l’eau » pour transformer un problème initial donné et obtenir les bonnes conditions
d’existence et d’unicité de la solution sur les « espaces qui vont bien » (et qui eux, feront
ensuite l’objet d’une discrétisation numérique).
Supposons que nos « investigations » sur un problème nous aient conduit à sa formulation
sous forme d’un problème de Neumann (avec f 2 L2 ./, ce qui n’est pas une condition si forte
finalement, et c’est ce que l’on aimerait avoir) :
8
< u D f dans
(9.1)
: @u D 0 sur D @
@n
Nous allons donc nous poser la question de l’existence et de l’unicité de la solution de ce problème.
qui conduit à :
Z
f D0 (9.3)
i.e. f 2 L20 ./. Donc f 2 L20 ./ est bien une condition nécessaire d’existence.
De même, d’après ce qui précède, on voit que si u est solution, alors u C c est solution, ce qui
prouve qu’il n’y a pas unicité de la solution.
@u
Z Z Z
u v D ru rv v (9.4)
@n
Comme nous voulons la formulation variationnelle dans H1 ./ (i.e. u 2 H1 .//, il vient :
Z Z
ru rv D f v; 8v 2 H1 ./ (9.6)
ce qui implique, pour que ce soit vrai pour tout d que f D 0, et on retrouve bien que f 2 L20 ./
R
sur H1 ./, donc V est un sous-espace fermé de H1 , donc est un Hilbert. Et donc maintenant
que l’on sait que V est un Hilbert, on aimerait bien appliquer Lax-Milgram (voir 7.3), c’est
pourquoi il faut s’intéresser à la forme bilinéaire définissant le problème. C’est ce qui suit.
— Si l’on considère la forme bilinéaire :
Z
a.u0 ; v0 / D ru0 rv0 (9.9)
Cela conduit à :
jv0 j21 1
2 2
> ; 8v 2 V (9.12)
jv0 j1 C jv0 j0 1 C c 2 ./
et la loi de comportement :
Z Z Z
0D .p ru/ q D pqC u div q D hu; q ni (9.15)
On n’a pas d’information sur u. Par contre, on a un terme p n D 0. Et il serait bien de disposer de
son « symétrique », i.e. d’un terme q:n D 0 qui serait imposé par le choix de l’espace dans lequel
vivrait u. Ce que l’on a en tête c’est bien sûr de pouvoir utiliser le théorème de Brezzi.
En correspondance avec les notations du paragraphe 7.5, on choisit alors les espaces : H D
fq 2 [Link] /I q nj D 0g et M D L2 ./. Le problème est donné, dans la formulation de Brezzi,
par les deux formes bilinéaires
Z Z
a.p; q/ D p q et b.v; q/ D v div q (9.16)
Quant aux formes linéaires définissant le problème, la première est nulle et la seconde vaut :
Z
fv (9.17)
Maintenant que nous nous sommes rapprochés du cadre du théorème de Brezzi, il va nous
falloir vérifier si toutes les conditions sont vérifiées. Toutefois, d’après ce qui a été vu dans les
paragraphes précédents, il semble naturel de se demander s’il faudra bien conserver L2 ./, ou s’il
faudra passer à L20 ./, et ce qu’il se passe quand u D u0 C c.
Pour l’instant, on a :
Z Z
v div p D f v; 8v 2 L2 ./ (9.18)
Z Z
pqC u div q D 0; 8q 2 H (9.20)
R
Or q 2 H implique div q D 0, ce qui implique aussi que si u est solution, alors u C c l’est aussi.
1
Z
8v 2 M; 9q 2 HI v div q D jvj20 et [Link]/ 6 jvj0 (9.27)
ˇ
@'
q 2 H, div q D v et q D r', q n D 0 donc @n
D 0, et :
8
< ' D v
(9.28)
: @' D 0
@n
Comme au paragraphe sur la formulation variationnelle, 9Š' 2 H1 ./ \ L20 ./ car v 2 L20 ./ :
Z Z
2 2 2
jqj0 D jr'j D j'j1 D v' 6 jvj0 j'j0 6 c./jVj0 j'j1 (9.29)
d’où :
1
ˇDp (9.31)
1 C c 2 ./
En introduction à cette partie, il nous semblait important d’en exposer sa structure, car elle peut
sembler un peu décousue à la simple lecture de la table des matières.
Après les pré-requis exposés dans les deux premières parties, la partie III va s’ouvrir « tout
naturellement » sur une présentation générale de la méthode des éléments finis au chapitre 10.
Immédiatement après, compte tenu du public visé, le chapitre 11 essayera de mettre en avant
les spécificités et surtout la complexité de la mécanique comme champ d’application de la méthode
des éléments finis.
Cette mise en garde, au regard de l’expérience, nous semble importante : on a tendance souvent
à considérer que la mécanique est quelque chose de très bien maîtrisé, et ce n’est pas le cas. Bien
des sujets restent pointus, voire ouverts. Il convient donc de rester prudent, surtout pour ceux ayant
une expérience de calcul importante qui les porte parfois à sous estimer les difficultés.
Le chapitre 12, qui fait suite logiquement en terme de présentation de la méthode des éléments
finis, au chapitre 10, est souvent le mieux maîtrisé par le public d’ingénieurs, au moins concernant
les éléments de Lagrange. Nous l’avons complété de remarques sur les modèles à plusieurs champs
(qui peut faire pour partie écho au chapitre 11). Encore une fois, pour le public visé, c’est surtout le
paragraphe 12.4 sur la validation pratique des éléments finis qui aura sans doute le plus de valeur
ajoutée. Le contenu de ce paragraphe fait souvent partie des choses oubliées.
Comme nous en serons sur des choses un peu oubliées, nous en profiterons au chapitre 13 pour
continuer dans le même sillon et rappeler quelques méthodes d’amélioration de la performance du
calcul. Nous y présentons des choses qui sont utilisées plus ou moins fréquemment par notre public.
Seuls les paragraphes 13.6 et 13.7 (sur les méthodes de réanalyse et aux dérivées d’ordre élevé)
sont généralement moins bien connus. Nous avons voulu les introduire ici plutôt qu’au chapitre 22
car ils sont vraiment en lien avec les préoccupations directes de notre public.
Le chapitre 14 permettra une petite pause en exposant brièvement les stratégies de maillage, et
plus particulièrement la triangulation de Delaunay.
Pour continuer sur les choses pouvant avoir un réel intérêt pratique pour le public d’ingénieurs
mécaniciens (ou acousticiens), nous aborderons au chapitre suivant 15 les méthodes d’homogé-
néisation. Si les méthodes les plus « physiques » sont bien connues, l’approche mathématique
(généralisante) est bien souvent quasi totalement inconnue.
À ce niveau du document, il nous semble que nous aurons parcouru bon nombre des applications
typiques de la méthode des éléments finis, surtout dans le cadre de la mécanique... mais uniquement
sous l’angle statique !
Il sera donc temps d’aborder « le temps », i.e. les problèmes non stationnaires, au chapitre 16.
La propagation des ondes, quant à elle, ne sera traitée qu’au chapitre 17. C’est d’ailleurs dans ce
chapitre que seront abordés également les modes propres, dont nous n’aurons pas parlé jusqu’alors
(à quelques exception près lors de remarques diverses et variées... mais rien de sérieux). Le
chapitre 18 continuera la spécification en se focalisant sur l’acoustique.
Dès lors, on pourra considérer qu’une présentation assez complète de la méthode des éléments
finis a été faite. Toutefois, ce qui était encore de l’ordre de la recherche il y a une décennie fait
aujourd’hui partie des phénomènes que notre public doit prendre en compte de plus en plus souvent.
Ces phénomènes, un peu plus complexes, un peu plus exotiques, sont souvent liés à ce que l’on
Résumé — Une fois le travail précédent accompli, i.e. une fois que l’on dispose d’une
formulation faible, « il n’y a plus qu’à » calculer la solution ! La méthode des éléments
finis est l’un des outils numérique développé pour cela.
La méthode des éléments finis se propose de mettre en place, sur la base de formula-
tions faibles, un algorithme discret (discrétisation) permettant de rechercher une solution
approchée d’un problème aux dérivées partielles sur un domaine compact avec conditions
aux bords et/ou dans l’intérieur du compact.
Il s’agit donc de répondre aux questions d’existence et d’unicité de la solution, de
stabilité et de convergence des méthodes numériques, ainsi que d’apprécier l’erreur entre
la solution exacte et la solution approchée (indicateurs et estimateurs d’erreur, a priori et
a posteriori).
10.1 Introduction
Les chapitres des parties précédentes ont eu pour but de nous permettre de décrire un problème
physique à partir d’équations aux dérivées partielles, mais à aucun moment nous n’avons encore
parlé de la manière de résoudre ces équations (que ce soit la formulation forte ou la faible). Ce sera
le but de ce chapitre, dans lequel la mise en œuvre de la méthode des éléments finis va être exposée.
Nous avons vu dans la partie précédente, comment passer d’une formulation forte à une
formulation faible. Dans de nombreux cas, nous avons montré qu’il y a équivalence (existence et
unicité de la solution) entre l’étude des deux formulations. Il faut toutefois bien garder en tête que
dans certains cas, l’équivalence entre les deux formulations n’est pas évidente.
Une condition qui n’a pas été vraiment évoquée jusque là est le cas où le bord du domaine
n’est pas suffisamment régulier, par exemple s’il possède des points singuliers (par exemple point
d’inflexion ou de rebroussement, ou dans le cas des fissures en MMC).
Cela peut également se produire si l’on ne peut pas assurer que la fonction f (celle agissant
dans tout le domaine) n’est pas suffisamment dérivable. Comme en général on suppose dans les
problèmes physiques que la solution est C1 , on n’est pas confronté à ce genre de problème.
L’idée de la méthode des éléments finis est de décomposer (on dit discrétiser) le domaine
en un certain nombre de sous-domaines (les éléments). Les éléments recouvriront l’intégralité
du domaine (de la frontière pour les éléments de frontière qui est une autre méthode) et sans
chevauchement entre eux (les éléments peuvent se chevaucher dans la méthode des volumes finis).
De plus, on va chercher la fonction solution u comme étant interpolée par des « bouts » de solutions
définis sur chaque élément.
Le problème étant interpolé sur les éléments, on se doute que le nombre d’éléments va jouer
sur la qualité de cette approximation de la solution. On se doute également que, comme on résout
un problème comportant des dérivées, c’est plutôt dans les endroits où la solution va varier vite
qu’il sera nécessaire de « resserrer » le maillage.
C’est l’ingénieur américain Ray William Clough qui, semble-t-il, a utilisé le terme de méthode des
éléments finis le premier dans un article de 1960 intitulé The Finite Element Method in Plane Stress
Analysis. Le mot rigidité (Stiffness) apparaissait dans le titre de son article Stiffness and Deflection
Analysis of Complex Structures datant de 1956 (et coécrit avec M. Turner, H. C. Martin et L. J.
Topp).
Si on veut replacer très brièvement cela dans un contexte
plus global, on peut dire que l’analyse des structures est née
vers 1850.
La RdM, recourant au calcul manuel, était développée
par Maxwell, Castigliano, Mohr.
Le concept d’éléments finis est né vers 1940, avec des
figures comme Newmark, Hrennikoff (1941), Mc Henry, Clough Courant Zienkiewicz
Courant (1942)...
Son réel essor ne commence toutefois que dans les années 60 avec le développement du calcul
numérique sur ordinateur.
La méthode des éléments finis (MEF) prend ses origines dans le besoin de résoudre des pro-
blèmes complexes d’élasticité et d’analyse de structures en ingénierie civile et aéronautique. Son
développement remonte aux travaux d’Alexander Hrennikoff (1941) et de Richard Courant (1942).
Bien qu’utilisant des approches différentes, ces deux pionniers partagent la même caractéristique
essentielle à savoir la discrétisation par maillage du domaine continu en sous-domaines discrets,
que l’on appelle éléments. C’est Olgierd Zienkiewicz de l’Imperial College qui synthétisa ces deux
méthodes en ce que l’on peut appeler la méthode des éléments finis et qui fit la première formalisation
mathématique de la méthode.
Dans ses travaux, Hrennikoff discrétise le domaine
en utilisant une analogie avec les treillis, tandis que l’ap-
proche de Courant divise le domaine en sous-régions finies
triangulaires pour résoudre les équations aux dérivées par-
tielles elliptiques du second ordre issues du problème de
la torsion d’un cylindre. On peut dire que la contribution
Rayleigh Ritz Galerkine
de Courant était une évolution s’appuyant sur un vaste cor-
pus de résultats antérieurs pour les équations aux dérivées
partielles développés par Rayleigh, Ritz et Galerkin.
Le développement de la méthode des éléments finis
a véritablement commencé au milieu de années 1950 pour
l’analyse structurale et aéronautique, et prit de l’ampleur
à l’Université de Stuttgart grâce au travail de John Argyris
et à Berkeley grâce au travail de Ray W. Clough [21].
Ray Clough est également l’un des pionniers du génie
para sismique et s’est vu décerné en 2008, à la World Argyris Strang Fix
Conference of Earthquake Engineering en Chine, le titre
de « légende du génie para sismique » (“Legend of Earthquake Engineering”).
À la fin des années 50, les concepts clés de matrice de rigidité et d’assemblage d’éléments
existaient quasiment sous la forme actuelle. La NASA publia une demande de propositions pour le
développement du logiciel d’éléments finis NASTRAN en 1965.
La base mathématique rigoureuse de la méthode des éléments finis a été consolidée en 1973
avec la publication de Strang et Fix de An Analysis of The Finite Element Method. Elle a depuis
été intégrée comme une branche des mathématiques appliquées à la modélisation numérique des
systèmes physiques dans une large variété de disciplines. Pour une discussion plus approfondie des
apports et contributions relatives des différents pionniers de cette méthode, on pourra se référer à
[33].
Bien que la méthode des éléments finis soit théoriquement généralisable à toutes les dimensions
d’espace et à tous les ordres de dérivation, dans la pratique, l’augmentation de ces paramètres
ρ(K)
En fait, on expose généralement les méthodes à partir de la méthode dite conforme selon le
tableau ci-dessus, alors qu’en pratique on travaille souvent dans le second cas, i.e. approximation
conforme en espace (Vh V) mais non conforme concernant les formes, car au moins les
intégrations sont réalisées numériquement.
Toujours est-il que dans tous les cas, l’espace V est remplacé par un espace Vh , de dimension
finie Nh , dont .e1 ; : : : ; eNh / en est une base. L’approximation uh de u dans cette base s’écrit :
Nh
X
uh D qi ei (10.2)
i D1
Le problème de la méthode des éléments finis devient donc (on écrit ici en utilisant ah , mais on
pourrait utiliser a) :
Nh
X
Trouver q1 ; : : : ; qNh tels que qi ah .ei ; vh / D fh .vh / ; 8vh 2 Vh (10.3)
i D1
III ÉLÉMENTS FINIS 10.3 Principe de la méthode : résolution d’un système matriciel 132
que l’on note (dans la tradition mécanicienne) :
Kq D F ou encore K q D F (10.6)
Telle que formulée, la matrice K semble pleine à priori. L’astuce consiste à choisir des fonctions
de base ei dont le support sera petit, i.e. chaque fonction ei sera nulle partout sauf sur quelques
mailles. Ainsi les termes ah .ei ; ej / seront le plus souvent nuls et la matrice sera creuse. De plus,
on ordonnera les ei de sorte que K soit à structure bande, avec une largeur de bande la plus faible
possible. Il existe un moyen avantageux de stocker une matrice creuse qui s’appelle le stockage en
ligne de ciel ou skyline.
Les difficultés majeures en pratique sont de trouver les ei et de les manipuler pour les calculs
d’intégrales nécessaires à la construction de K. Indiquons d’ores et déjà que la plupart de ces
difficultés seront levées grâce à trois idées principales qui seront détaillées au chapitre sur la
formulation pratique des éléments finis :
Principe d’unisolvance : On s’attachera à trouver des degrés de liberté (ou ddl) tels que la donnée
de ces degrés de liberté détermine de façon univoque toute fonction de Vh . Il pourra s’agir
par exemple des valeurs de la fonction en quelques points. Déterminer une fonction reviendra
alors à déterminer ses valeurs sur ces degrés de liberté.
Définition des ei : On définira les fonctions de base par ei D 1 sur le i ème degré de liberté,
et ei D 0 sur les autres degrés de liberté. La manipulation des ei sera alors simplifiée, et
les ei auront par ailleurs un support réduit à quelques mailles.
Notion de « famille affine d’éléments » : Le maillage sera tel que toutes les mailles soient iden-
tiques à une transformation affine près. De ce fait, tous les calculs d’intégrales pourront
se ramener à des calculs sur une seule maille de référence, par un simple changement de
variable.
Notons que la matrice K est appelée matrice de rigidité par analogie avec la mécanique des solides.
Si la forme bilinéaire a est coercive, alors K est symétrique, définie positive donc inversible.
On obtient donc l’existence et l’unicité de q D K 1 F.
De nombreuses méthodes permettent de résoudre le système matriciel (d’inverser la matrice de
rigidité) : Gauß, mais également toutes sortes de décomposition de la matrice de rigidité comme
LU, LDLT, LLT (Cholesky). Lorsque K est symétrique définie-positive, la méthode de Cholesky
est la meilleure. Elle consiste à décomposer la matrice K en le produit LT L, où la matrice L est une
matrice triangulaire inférieure.
Lemme 1 — Lemme de Céa. La forme bilinéaire a.; / étant continue, de constante de majora-
tion M, et coercitive, de constante de minoration ˛, il est aisé d’obtenir la majoration de l’erreur
appelée lemme de Céa :
M M
ku uh k 6 ku vh k; 8vh 2 Vh c’est-à-dire ku uh k 6 d.u; Vh / (10.7)
˛ ˛
où d est la distance induite par la norme k k.
Cette majoration donnée par le lemme de Céa ramène l’étude de l’erreur d’approximation
ku uh k à celle de l’erreur d’interpolation d.u; Vh /.
Lemme 2 — Lemme 1 de Strang. S’il existe une constante ˛ > 0 indépendante de h telle
que :
Lemme 3 — Lemme 2 de Strang. Si la norme k kh est définie sur Vh C V, et s’il existe deux
constantes indépendantes de h, ˛ > 0 et M > 0 telles que :
III ÉLÉMENTS FINIS 10.6 Convergence de la méthode des éléments finis en approximation conforme interne 136
— ku h ukm 6 C0 hkC1 m juj
kC1 : assemblage des majorations locales.
M
ku uh k 6 ku h uk (10.22)
˛
et le calcul est ramené à un calcul sur chaque élément, pour toutes les semi-normes k kl;K pour
l D 0,. . . , m.
Ce résultat n’est rien d’autre qu’un résultat de changement de variable dans une intégrale.
On a même :
8v 2 Hl .K/; jvj
O l;K 6 C.l; n/kB 1 l
k2 j det Bj1=2 jvjl;KO (10.25)
hK 1 hO
kBk 6 et kB k6 (10.26)
O K
O /;
9C.K; O 8vO 2 HkC1 .K/;
O jvO O vj
O j;KO 6 Cjvj
O kC1;KO (10.27)
ku uh km 6 ChkC1 m
jujkC1 (10.30)
Quelques remarques :
— On rappelle que la formule précédente a été obtenue pour un domaine polygonal. Si ce
n’est pas le cas, elle n’est plus valable. Les éléments linéaires conduisent à une erreur en
O.h/. Utiliser des éléments plus raffinés (de degré 2 au lieu d’éléments linéaires) permet,
même si la géométrie n’est également ici pas décrite de manière exacte, d’obtenir une erreur
asymptotique en O.h2 /. Quelque soit le degré de l’approximation, si le domaine n’est pas
représenté de manière exact, alors cela revient à une modification des conditions aux limites.
— Les calculs, et plus précisément les intégrations, ont été réalisées sans erreur. Si celles-
ci sont faites de manière numérique, il convient d’introduire encore un terme correctif
supplémentaire, appelé erreur de consistance due au remplacement des formes par leur
approximation (a.; / par ah .; /, f ./ par fh ./...). Toutefois, cette erreur supplémentaire
peut être estimée en [Link] / selon la précision du schéma d’intégration utilisé.
— Le résultat de majoration d’erreur est souvent utilisé dans le cas m D 1. Comme l’espace des
polynômes Pk .k/O H1 .K/,O alors si O est bien défini sur HkC1 .K/, O on a :
III ÉLÉMENTS FINIS 10.6 Convergence de la méthode des éléments finis en approximation conforme interne 138
Chapitre 11
La mécanique est sans doute aussi vieille que l’homme. Aussi bien pour des aspects pratiques (faire
des outils pour chasser...), que pour des aspects plus philosophiques et spirituels visant notamment à
expliquer les mouvements des astres...
Archimède, outre ces travaux en mathématiques, pour-
rait sans conteste être considéré comme le saint patron de
la mécanique. Il est tout au moins indubitablement le père
de la mécanique statique.
Il s’intéressa aussi bien aux aspects « théoriques »
portant sur le principe du levier et la recherche de centre
de gravité dans De l’équilibre des figures planes, sur le Aristote Archimède Galilée
principe d’Archimède pour les corps plongés dans un li-
quide dans Des corps flottants ; qu’aux aspects « pratiques » au travers de nombreuses inventions :
machines de traction (où il démontre qu’à l’aide de poulies, de palans et de leviers, l’homme peut
soulever bien plus que son poids), machines de guerre (principe de la meurtrière, catapultes, bras
mécaniques utilisés dans le combat naval), l’odomètre (appareil à mesurer les distances), la vis sans
fin et la vis d’Archimède, le principe de la roue dentée...
Le siège de Syracuse, les miroirs d’Archimède et sa mort ne font qu’ajouter à sa légende.
Bien qu’Aristote posa le premier (avant Archimède) les bases d’une véritable théorie mécanique
(alors encore très imparfaite), les fondements de la mécanique, en tant que science et au sens moderne
du terme, sont posées par Galilée en 1632 dans les Dialogues et en 1638 dans les Discours.
La mécanique n’est alors pas dissociée des arts mécaniques, i.e. des techniques de construction
des machines. La distinction entre la mécanique en tant science et la mécanique en tant que technique
ne se fera qu’au XIXe siècle.
u D Nu q (11.1)
De la même manière, si les déformations " sont choisies comme champ inconnu, elles se-
ront approchées à partir des déformations nodales (sous forme de vecteur avec la convention de
l’ingénieur)
par l’intermédiaire de fonctions de formes rangées dans la matrice N" .
Cela vaut également pour les contraintes qui, si elles sont choisies comme champ inconnu,
seront approchées à partir des contraintes nodales (sous forme de vecteur avec la convention de
l’ingénieur) par l’intermédiaire de fonctions de formes rangées dans la matrice N .
D N (11.3)
Notons que dans la pratique, rien n’empêche de prendre les mêmes fonctions de forme pour les
différents champs. De manière analogue, le vecteur des multiplicateurs de Lagrange sera interpolé
de la façon suivante :
D N L (11.4)
On aura donc :
" D N"
approximation en déformations
D Lu relation déformations déplacements
D LNu q approximation en déplacements (11.5)
D S loi de Hooke inverse
D SN approximation en contraintes
D N approximation en contraintes
D H" loi de Hooke généralisée
D HN"
approximation en déformations (11.6)
D HLu loi de Hooke en élasticité linéaire
D HLNu q approximation en déplacements
où les conditions subsidiaires sur les déplacements sont prises en compte directement par l’espace
dans lequel les déplacements sont recherchés.
La fonctionnelle d’Hellinger-Reissner est sans doute la plus connue des fonctionnelles mixtes.
Elle utilise les champs de déplacements et de contraintes comme variables indépendantes. Son
expression est la suivante, U désignant les déplacement imposés :
Z Z Z
1
…HR D 2 S ij;j u C f u d T T u d UT d (11.10)
u
bien que l’on puisse la trouver sous une autre forme, obtenue par intégration par parties de celle-ci.
Elle n’a à satisfaire à aucune condition subsidiaire. La stationnarité de cette fonctionnelle conduit
aux équations d’équilibre, à la loi de comportement, et aux conditions aux limites en forces et
déplacements.
Cette fonctionnelle conduit à un élément ayant les champs de déplacements et de contraintes
comme inconnues nodales. Il s’en suit que toutes les composantes de ces champs sont continues.
Elle conduit à la résolution d’un système du type mixte (Brezzi) :
A B O
D (11.11)
BT O q F
L’un des moyens pour obtenir une fonctionnelle hybride est d’introduire une condition sur une
partie du contour par l’intermédiaire de multiplicateurs de Lagrange comme présenté un peu plus
loin.
Tout comme la fonctionnelle d’Hellinger-Reissner est associée à la méthode mixte, celle de
Pian et Tong est indissociable de l’adjectif hybride, même si, en toute rigueur, elle est une méthode
mixte (deux champs) hybride (différents domaines d’interpolation) :
Z Z Z
1
…PT D 2 S d C Tu d Tu d (11.13)
Sa variation est :
Z Z Z
ı…PT D S d C ıTu C ıuT d ıuT d (11.14)
Toutefois, le champ de contraintes n’étant défini que dans , il est possible d’effectuer une
condensation statique du champ de contraintes, i.e. de transformer le système (11.11) en écrivant :
1
DA Bq (11.17)
Sur le plan du déroulement du calcul, on commence par résoudre une « forme classique »
mais en utilisant la matrice de rigidité équivalente Keq , puis le calcul des contraintes dans chaque
élément se fait par la relation (11.17).
Les multiplicateurs de Lagrange peuvent être utilisés pour introduire des conditions supplémen-
taires directement à la fonctionnelle. Ces conditions sont généralement imposées sur X , tout ou
partie de D @. Pour ce faire, il suffit d’ajouter à la fonctionnelle … du problème le terme :
Z
T
˙ condition (11.19)
X
La figure 11.1 propose trois modélisations d’un même problème. Il s’agit d’une poutre encastrée
à une extrémité et soumise à une force décentrée à l’autre. Vaut-il mieux modéliser l’intégralité du
volume de la poutre, la traiter comme une plaque, ou peut-on se contenter d’un modèle de poutre ?
L’étude de la contrainte axiale xx est donnée à la figure 11.2. Si les cartographies présentent
bien la même répartition, seul le modèle 3D permet de mettre en évidence la concentration de
contrainte due à la force ponctuelle.
Une comparaison plus fine des trois modélisations concernant cette même contrainte normale
III ÉLÉMENTS FINIS 11.2 Plusieurs modélisations d’un même problème 144
F
F
F
(a) (b)
F IGURE 11.2: Contrainte axiale sur la peau supérieure : (a) modèle 3D, (b) modèle 2D.
le long de trois lignes est reportée à la figure 11.3. Les différences de résultats proviennent des
L3
L1
F
L2
1.80
Plaque 1.80
Plaque 1.80
Plaque
Poutre Poutre Poutre
1.60 Volumique 1.60 Volumique 1.60 Volumique
1.40 1.40 1.40
0.40
L1 0.40
L2 0.40
L3
0.20 0.20 0.20
ABS ABS ABS
0.00 0.00 0.00
0.00 0.20 0.40 0.60 0.80 1.00 0.00 0.20 0.40 0.60 0.80 1.00 0.00 0.20 0.40 0.60 0.80 1.00
hypothèses cinématiques et des composantes accessibles dans les différentes théories utilisées. On y
voit que loin des extrémités, les trois solutions concordent parfaitement, conformément au principe
de Saint-Venant. 1 Les différences proviennent des effets de bord, i.e. des variations locales des
contraintes et déformations au voisinage des conditions aux limites.
Les hypothèses cinématiques sont :
— la théorie des poutres (Figure 11.4a) suppose que chaque section droite suit un mouvement
de solide rigide. Les sections ne peuvent donc pas se déformer, elles peuvent uniquement se
translater et tourner dans l’espace (sans forcément rester perpendiculaires à la ligne moyenne,
car dans ce cas on considère une théorie avec cisaillement transverse) ;
— la théorie des plaques (Figure 11.4b), moins restrictive, suppose que chaque segment per-
pendiculaire au plan moyen de la plaque suit un mouvement de solide rigide (là encore, sans
1. Le principe de Saint-Venant précise que le comportement en un point quelconque de la poutre, pourvu que ce
point soit suffisamment éloigné des zones d’applications des forces et des liaisons, est indépendant de la façon dont
sont appliquées les forces et de la façon dont sont physiquement réalisées les liaisons ; le comportement dépend alors
uniquement du torseur des forces internes en ce point.
F IGURE 11.4: Mouvements d’une section droite : (a) théorie des poutres ; (b) théorie des plaques
ou coques ; (c) théorie volumique
forcément rester perpendiculaire au plan moyen). Elle permet donc de modéliser certaines
formes de déformations des sections, planes (flexion et cisaillement dans le plan de la section,
traction dans le sens de la largeur) ou hors plan (certains types de gauchissement). Cepen-
dant, les segments ne pouvant pas changer de longueur, elle ne permet pas de modéliser
l’écrasement de l’épaisseur ;
— enfin, la théorie 3D (Figure 11.4c) ne comporte aucune de ces restrictions et peut modéliser
n’importe quelle forme de gauchissement ou d’écrasement, à condition que le maillage
employé soit suffisamment fin.
Il est essentiel de noter que toutes les théories ne permettent pas d’accéder à toutes les com-
posantes du champ des contraintes : de manière générale, seuls les efforts qui travaillent dans les
déplacements permis par la théorie sont accessibles (les théories sont ainsi faites afin de pouvoir
respecter le premier principe de la thermodynamique dont l’écriture nécessite de calculer le travail
de tous les efforts permis par la théorie). Ainsi :
— la théorie des poutres ne permet pas de calculer les contraintes dans le plan transversal (yy ,
zz et yz ) du fait de l’indéformabilité des sections droites ;
— la théorie des plaques ne permet pas de calculer la contrainte normale au plan de la plaque
zz , du fait de l’indéformabilité des segments perpendiculaires à ce plan ;
— enfin, la théorie 3D permet de calculer les six composantes du tenseur des contraintes.
À ces limitations théoriques peuvent s’ajouter des limitations techniques propres à chaque
logiciel, susceptibles de rendre d’autres grandeurs physiques inaccessibles. Avant d’effectuer une
modélisation par éléments finis, il est donc indispensable de s’assurer que le logiciel utilisé et son
cadre théorique permettent bien d’accéder au résultat voulu (en plus d’être pertinents vis-à-vis de la
géométrie du produit, de son environnement et du comportement attendu).
De plus, si certains logiciels peuvent avoir des « limitations techniques », ils peuvent également
disposer de méthodes de post-traitement susceptibles d’améliorer ou de permettre d’accéder à
certaines composantes (sous certaines hypothèses). Cela aussi doit être pris en compte.
III ÉLÉMENTS FINIS 11.2 Plusieurs modélisations d’un même problème 146
11.3 Exemple : retour sur le calcul de poutre du paragraphe 11.1
avec C AST 3M
Nous allons reprendre pour partie l’exemple de la poutre encastrée traitée au paragraphe précédent
et présenter le listing C AST 3M correspondant.
11.3.1 Modélisation 2D
Nous commençons par définir les données du problème : longueur, largeur, épaisseur, nombre
d’éléments selon chacune de ces directions, et force appliquée :
1 DONNEES
2 geometrie
3 long1 =22.0;
4 larg1 =8.;
5 ep1 =4;
6 maillage
7 nlong1 =22;
8 nlarg1 =8;
9 effort
10 Forc1 = -41.;
Nous définissons les points ki (dénommés ainsi pour rappel de la syntaxe A NSYS des keypoints
k,i,...), puis les lignes Li , et la surface S1 .
12 k1 = 0. 0. 0.;
13 k2 = long1 0. 0.;
14 k3 = long1 larg1 0.;
15 k4 = 0. larg1 0.;
16
17 L1 = DROI nlong1 k1 k2 ;
18 L2 = DROI nlarg1 k2 k3 ;
19 L3 = DROI nlong1 k3 k4 ;
20 L4 = DROI nlarg1 k4 k1 ;
21
22 S1 = DALLER L1 L2 L3 L4 ;
Le modèle est un modèle de mécanique élastique isotrope utilisant l’élément de coque COQ4 pour
le maillage S1 :
1 Model1 = MODL S1 MECANIQUE ELASTIQUE ISOTROPE COQ4 ;
Enfin on résout le problème après avoir fourni les données matérielles et les conditions aux
limites :
2 Mater1 = MATER Model1 YOUNG 70000.0 NU 0.33 RHO 2700.0;
3 Car1 = CARAC Model1 EPAI ep1 ;
4 Mater1 = Mater1 ET Car1 ;
5 MR1 = RIGID Model1 Mater1 ;
6 CL1 = BLOQ DEPL L4 ;
7 CL2 = BLOQ ROTA L4 ;
8 FOR1 = FORC (0. 0. Forc1 ) k3 ;
9 MTOT1 = MR1 ET CL1 ET CL2 ;
10 Depl1 = RESO MTOT1 FOR1 ;
11.3.2 Modèle 3D
Nous allons maintenant construire le modèle tridimensionnel. Comme nous voulons travailler à
l’économie, nous repartons du fichier précédent que nous adaptons.
1 DONNEES
2 geometrie
3 long1 =22.0;
4 larg1 =8.;
5 ep1 =2;
6 maillage
7 nlong1 =22;
8 nlarg1 =8;
9 nep1 =3;
10 effort
11 Forc1 = -41.;
Cette fois, nous nous servons de l’élément volumique CUB8 au lieu de l’élément surfacique QUA4.
1 OPTION DIME 3 ELEM CUB8 ;
2
3 k1 = 0. 0. 0.;
4 k2 = long1 0. 0.;
5 k3 = long1 larg1 0.;
6 k4 = 0. larg1 0.;
7
8 L1 = DROI nlong1 k1 k2 ;
9 L2 = DROI nlarg1 k2 k3 ;
10 L3 = DROI nlong1 k3 k4 ;
11 L4 = DROI nlarg1 k4 k1 ;
12
13 S1 = DALLER L1 L2 L3 L4 ;
À partir de la surface S1 , qui est la même que précédemment, nous allons créer le volume V1
par translation selon le vecteur Vect1 .
Nous en profitons également pour créer Face1 sur laquelle porterons les conditions aux limites.
Notons par exemple que la ligne L4 appartient bien au maillage V1 , puisque ce dernier est construit
dessus. Par contre, la surface f ace1 n’appartient pas à V1 , il est donc nécessaire de la lier à V1 en
utilisant la commande ELIM.
1 Vec1 =0. 0. ( -1.0* ep1 ) ;
2 V1 = S1 VOLU nep1 TRAN Vec1 ;
3 Face1 = L4 TRAN nep1 Vec1 ;
4 ELIM 0.0000001 V1 Face1 ;
III ÉLÉMENTS FINIS 11.3 Exemple : retour sur le calcul de poutre du paragraphe 11.1 avec C AST 3M 148
Cette fois, le modèle correspond au maillage V1 et utilise l’élément CUB8.
1 Model1 = MODL V1 MECANIQUE ELASTIQUE ISOTROPE CUB8 ;
Ceci est juste une remarque en passant, pour définir le vocabulaire en somme, nous n’en
reparlerons plus dans la suite du document.
Résumé — L’intégralité de la méthode des éléments finis a été présentée au chapitre 10.
Dans ce chapitre et dans les suivants, nous allons détailler certains aspects. Nous
proposons dans ce chapitre d’exposer un peu plus complètement les notions d’interpo-
lation sur un élément, ainsi que le lien entre approximation locale (sur un élément) et
approximation globale (construction de la base de Vh ).
Nous avons dit vouloir interpoler le problème sur chaque élément. Pour ce faire, il faut prendre
une base sur chaque élément. Plusieurs choix sont possibles, mais en général, les fonctions de base
utilisées pour les éléments finis sont dites interpolantes, c’est-à-dire que les valeurs nodales sont les
valeurs des grandeurs inconnues aux nœuds, et que c’est à partir de ces valeurs que l’on effectue
l’interpolation.
La méthode la plus simple consiste à utiliser les polynômes de Lagrange. Dans cette méthode
les fonctions de base valent 1 à un nœud du maillage et 0 à tous les autres. La fonction de base i est
alors la fonction valant 1 au nœud i et 0 sur les autres nœuds et qui est polynomiale sur chaque
élément. Il y a autant de fonctions de base par élément que de nombre de nœuds. On appelle
élément la donnée d’une géométrie (souvent polygonale en deux dimensions, polyédrique en trois
dimensions) et de fonctions de base associées à cette géométrie.
D’autres solutions peuvent exister pour les fonctions de base. Par exemple, les éléments finis
d’Hermite ont la particularité d’avoir deux fonctions de base associées à chaque nœud. La valeur
de la solution est alors ajustée avec la première fonction alors que la deuxième permet d’ajuster
la valeur de la dérivée. Ce type de fonctions de base peut avoir un intérêt pour la résolution de
certaines équations aux dérivées partielles (telle que l’équation des plaques en Mécanique des
Milieux Continus), même si elle nécessite d’avoir deux fois plus de fonctions pour un maillage
donné.
12.1.1 Unisolvance
Définition 58 — Unisolvance. Soit † D fa1 ; : : : ; aN g un ensemble de N points distincts de Rn .
Soit P un espace vectoriel de dimension finie de fonctions de Rn à valeurs dans R. On dit que †
est P-unisolvant si et seulement si pour tous réels ˛1 ,. . . , ˛N , il existe un unique élément p de P
tel que 8i D 1; : : : ; N, [Link] / D ˛i .
Définition 59 — Éléments finis de Lagrange. Un élément fini de Lagrange est un triplet .K; †; P/
tel que :
— K est un élément géométrique de Rn , compact, connexe, et d’intérieur non vide ;
— † D fa1 ; : : : ; aN g est un ensemble de N points distincts de Rn ;
— P est un espace vectoriel de dimension finie de fonctions réelles définies sur K, et tel
que † soit P-unisolvant (donc dim P D N).
Les fonctions de bases locales de l’élément fini de Lagrange .K; †; P/ sont les N fonctions de P
telles que pi .aj / D ıij pour 1 6 i; j 6 N.
Remarque. .p1 ; : : : ; pN / est une base de P.
Théorème 46 k v est l’unique élément de P qui prend les mêmes valeurs que v sur les points
de †.
On notera Pk l’espace vectoriel des polynômes de degré total inférieur ou égal à k.
— sur R, Pk D vectf1; X; : : : ; Xk g et dim Pk D k C 1 ;
— sur R2 , Pk D vectfXi Yj ; 0 6 i C j 6 kg et dim Pk D .kC1/.kC2/
2 ;
3 i j l .kC1/.kC2/.kC3/
— sur R , Pk D vectfX Y Z ; 0 6 i C j C l 6 kg et dim Pk D 6 .
On notera Qk l’espace vectoriel des polynômes de degré inférieur ou égal à k par rapport à chaque
variable.
— sur R, Qk D Pk ;
— sur R2 , Qk D vectfXi Yj ; 0 6 i; j 6 kg et dim Qk D .k C 1/2 ;
— sur R3 , Qk D vectfXi Yj Zl ; 0 6 i; j; l 6 kg et dim Qk D .k C 1/3 .
Remarque. Les fonctions de base pour l’élément P1 sont définies par pi .aj / D ıij . Ce sont les coordon-
nées barycentriques : pi D i .
a2 a2
a4 a3
a12
a23
a1 a1
a13 a1 a2
a3 a3
a3 a8 a7
a3 a6 a5
a5 a6 a4
a34
a13 a23
a4
a3
a4 a4 a3
a2
a14 a24
a2 a2
a1 a1 a12 a1 a2 a1
Élément Q1
K parallélépipède de sommets fa1 ; : : : ; a8 g de côtés parallèles aux axes.
† fai g16i 68
P Q1
Élément Q1
K prisme droit de sommets fa1 ; : : : ; a6 g
† fai g16i 66
P fp.X; Y; Z/ D .a C bX C cY/ C Z.d C eX C f Y/; a; b; c; d; e; f 2 Rg
Définition 61 — Famille affine d’éléments finis. On appelle famille affine d’éléments finis
une famille d’éléments finis tous affine-équivalents à un même élément fini appelé élément de
référence.
D’un point de vue pratique, le fait de travailler avec une famille affine d’éléments finis permet de
ramener tous les calculs d’intégrales à des calculs sur l’élément de référence. Voir figure 12.3 pour
une illustration d’une transformation affine. Les éléments finis de références sont ceux obtenus,
a3
F
(0,1)
K
a1
K
a2
(0,0) (1,0)
F IGURE 12.3: Transformation affine sur un triangle
Revenons à notre problème décrit par une équation aux dérivées partielles sous forme faible dans
un domaine sur lequel on réalise un maillage Th à partir d’une famille affine de Ne éléments
finis .Ki ; †i ; Pi /i D1;:::;Ne .
Par unisolvance, la solution approchée uh sera entièrement définie sur chaque élément fini
par ses valeurs sur les points de †i , nommés nœuds du maillage. Notons .a1 ; : : : ; aNh / les nœuds
du maillage (Nh < Ne cardi ). Le problème approché revient à déterminer les valeurs de uh aux
points ai : ce sont les degrés de liberté du problème approché. On va construire une base de Vh en
associant à chaque degré de liberté ai un vecteur de la base. On définit ainsi les fonctions de base
globales 'i .i D 1; : : : ; Nh / par :
F IGURE 12.4: Base de Vh : exemple de fonction de base globale 'i sur un maillage avec des
éléments triangulaires P1 .
De plus, sur un élément K dont ai est un nœud, 'i vaut 1 en ai et 0 aux autres nœuds de K.
Donc 'i jK est une fonction de base locale de K. On voit donc que la fonction de base globale 'i est
construite comme réunion des fonctions de base locales sur les éléments du maillage dont ai est un
nœud.
Remarque. Ce qui précède est vrai dans le cas d’un maillage conforme, i.e si l’intersection entre deux
éléments est soit vide, soit réduite à un sommet ou une arête en dimension 2, ou à un sommet, une arête ou
une face en dimension 3.
Définition 62 — Élément fini d’Hermite. Un élément fini d’Hermite ou élément fini général
est un triplet .K; †; P/ tel que :
— K est un élément géométrique de Rn , compact, connexe, et d’intérieur non vide ;
— † D f1 ; : : : ; N g un ensemble de N formes linéaires sur l’espace des fonctions définies
sur K, ou sur un sous-espace plus régulier contenant P ;
— P est un espace vectoriel de dimension finie de fonctions réelles définies sur K, et tel
que † soit P-unisolvant.
Théorème 47 k v est l’unique élément de P qui prend les mêmes valeurs que v sur les points
de †.
On remarque immédiatement que si i .p/ D [Link] /; i D 1; : : : ; N, on retrouve les éléments finis
de Lagrange. Cette généralisation permet d’introduire des opérations de dérivation dans †, et donc
d’améliorer la régularité des fonctions de Vh . Les fonctions de base globales 'i , .i D 1; : : : ; Nh /
sont définies par :
Suivant les éléments utilisés, ces fonctions de base pourront être de classe C1 ou même plus, et
il en sera donc de même pour la solution approchée uh .
a4 a3
a12 a13
a123
a2
a3 a1 a2
a3 a23
a2
F IGURE 12.6: Éléments finis d’Hermite triangulaire cubique, élément d’Argyris et élément rectan-
gulaire Q3
Élément
K n triangle de sommets fa1o; a2 ; a3 g
@p
† [Link] /; @x .ai /; i D 1; 2; 3 [ fp.a0 /g
P P3
Régularité C0 , mais pas C1
Élément
K n triangle de sommets fa1 ; a2 ; a3 g o n o
@p 2 2 @2 p
† [Link] /; @x
.ai /; @p
@y
.ai /; @@xp2 .ai /; @@yp2 .ai /; @x@y .ai /; i D 1; 2; 3 [ @p
.a /; 1
@n ij
6i <j 63
P P5
Régularité C1
Élément Q3
K rectangle nde sommets fa1 ; a2 ; a3 ; a4 g de côtés parallèleso aux axes
@p 2
† [Link] /; @x .ai /; @p
@y
@ p
.ai /; @x@y .ai /; i D 1; : : : ; 4
P P3
Régularité C1
Ω1
Continuité :
Ω 1 σi
- des déplacements Ui n
ΓI n
- des contraintes Ti σi
Ω2
(σi désigne le vecteur des contraintes) Ω2
chacun de leur propre loi de comportement, ont une interface commune I . On supposera que
le champ de déplacement est continu le long de I . La continuité de la composante normale du
déplacement à l’interface traduit le fait qu’il n’y a pas décollement entre les deux domaines ; celle
de la composante tangentielle qu’il n’y a pas glissement entre eux. On peut donc dire qu’en tout
point de I , le déplacement est continu, ce que l’on note u1i D u2i .
L’état d’équilibre des forces doit lui-aussi être vérifié le long de l’interface. Cela implique
que les composantes normales des contraintes doivent être également continues le long de cette
interface, i.e. que la trace du tenseur des contraintes doit être continue le long de l’interface. Par
contre, les autres composantes peuvent (et doivent selon les cas) être discontinues.
Plusieurs stratégies sont envisageables :
Post-traitement : Dans cette méthode, on effectue le calcul en déplacements, de manière classique.
On obtient les contraintes de manière toute aussi classique, mais celles-ci n’ont évidemment
pas les continuités et discontinuités souhaitées.
On utilise une méthode de post-traitement qui va modifier le calcul des contraintes sur les
éléments situés de part et d’autre de l’interface, par exemple en imposant la vérification
des équations d’équilibre (voir par exemple la méthode de Reissner locale, qui applique la
fonctionnelle de Reissner uniquement le long de l’interface).
Élément mixte : La fonctionnelle d’Hellinger-Reissner possède les champs de déplacement et de
contrainte comme inconnues. Elle conduit à un système de type mixte donné par l’équation
(11.11), avec une matrice qui n’est plus définie-positive.
Il s’en suit que ces deux champs sont continus, et notamment que toutes les composantes des
contraintes le sont, ce qui ne satisfait pas les conditions d’équilibre.
Si l’on souhaite utiliser ce type d’élément, il faudra par exemple, faire une condensation
statique des composantes devant être discontinues, ou utiliser une méthode de post-traitement.
Élément hybride : La fonctionnelle de Pian et Tong (mixte hybride) présentée précédemment
Cette condition est facile à écrire et facile à implémenter (mais on perd la définie-positivité
de la matrice de rigidité).
De plus, l’interprétation physique de ces multiplicateurs de Lagrange montre que ceux-ci
sont égaux aux composantes normales des contraintes.
Voilà quelques stratégies. D’autres peuvent exister, surtout si en plus on veut passer sur des modèles
poutre ou plaque. Le but était de montrer que non seulement il est possible de développer de
nombreux éléments finis, mais que même à partir d’éléments existants, il est toujours possible de
construire une méthode numérique permettant d’obtenir les résultats souhaités. Les méthodes de
post-traitement, que nous n’aborderons pas dans ce document, sont très riches et permettent de
faire beaucoup de choses. Elles ont en outre l’avantage d’être parfois plus faciles à implémenter
dans des codes industriels que de nouveaux éléments.
Kuin D 0 (12.6)
Un mode rigide quelconque étant une combinaison linéaire des modes uin , on construit un problème
éléments finis dans lequel on impose mr valeurs du vecteur uin et on détermine numériquement
les n mr composantes (non imposées) de uin . Celles-ci doivent être identiques aux valeurs
théoriques associées au mode rigide i . Si ce n’est pas le cas, c’est qu’il existe des modes rigides.
12.4.3 Patch-tests
Le domaine choisi doit posséder au moins un nœud à l’intérieur du domaine. On définit un champ
de déplacement conduisant à un état de déformation désiré (constant ou non), que l’on introduit
comme condition aux limites dans le modèle élément fini. On vérifie qu’au point intérieur, l’état de
déformation obtenu est bien celui désiré.
y
ξ
-1 0 1 x
(a) Élément de référence (b) Segment réel (3D)
ses seules extrémités. Nous avons donc deux nœuds de coordonnées .x1 ; y1 ; z1 / et .x2 ; y2 ; z2 /.
Un point courant de cet élément sera obtenu par interpolation linéaire entre les deux nœuds,
paramétrée par . On Cherche cette interpolation sous la forme x./ D T N./xn , y./ D T N./yn
et z./ D T N./zn , avec T N./ D T .N1 ./; N2 .// et N1 ./ D 21 .1 /, N2 ./ D 12 .1 C /. De
manière plus « compacte », on peut écrire :
1
Ni ./ D .1 C i /; i D 1; 2; i D ˙1 .et 2 Œ 1; 1/ (12.7)
2
Dans cette interpolation (comme dans toute interpolation), on a une relation entre dx et d. No-
tons dx D ad. Alors :
T 1
a D T x; ; y; ; z; D T .x21 ; y21 ; z21 /
(12.8)
2
III ÉLÉMENTS FINIS 12.5 Exemple : quelques variations sur le thème des éléments unidimensionnels 162
Théorème 48 — Théorème d’inversion locale. Soit f une application de U dans F, où U est
un ouvert d’un espace de Banach réel et F un espace de Banach et soit x un point de U.
Si f est de classe Cp , avec p un entier strictement positif et si la différentielle de f au point
x (définie au paragraphe 3.5) est un isomorphisme bicontinu, alors il existe un voisinage ouvert V
de x et un voisinage ouvert W de f .x/ tels que f se restreigne en une bijection de V dans W
dont la réciproque est de classe Cp .
Comme illustré au paragraphe précédent, le jacobien sert surtout pour effectuer des changement de
variables dans le calculs des intégrales.
Si l’on considère un « petit » domaine, le volume de l’image de ce domaine par la fonction f sera
celui du domaine de départ multiplié par la valeur absolue du jacobien.
alors il vient :
... ...
1 i n
n
Y r
Ni ./ D (12.19)
rD1 r
i
r¤i
x
ξ
-1 0 1 x1
(a) Élément de référence (b) Élément réel (infini)
dans ce paragraphe, mais on pourrait étudier le cas de l’élément linéaire unidimensionnel à n nœuds,
ou des éléments bidimensionnels par exemple. Nous allons toujours avoir une interpolation de la
fonction solution sous la forme :
1 1C
uD u1 C u2 (12.21)
2 2
car concrètement u2 est connu (et fini) : c’est le « lieu » où l’on situe l’infini dans le modèle
éléments finis. Par contre, l’interpolation du point courant sera donnée par :
1C
x D x1 C ˛ (12.22)
1
où ˛ est une constante. On aura également :
u2 u1
u;x D .1 /2 (12.23)
4˛
avec les conditions aux limites u.0/ D 0 et v.0/ D 0. Cette forme variationnelle s’écrit pour
l’élément (ou pour chaque élément si on en avait plusieurs) :
W D T v .Ku f/ (12.25)
Nous avons vu que la géométrie, représentée par un élément, est approximée par :
1 1C
xD x1 C x2 ; 2 Œ 1I 1 (12.26)
2 2
III ÉLÉMENTS FINIS 12.5 Exemple : quelques variations sur le thème des éléments unidimensionnels 164
La représentation de la fonction solution, est :
1 1C
uD u1 C u2 ; 2 Œ 1I 1 (12.27)
2 2
Nous obtenons la matrice de rigidité élémentaire :
2EA 1 1
KD (12.28)
L 1 1
et le vecteur des forces nodales élémentaires :
L 1
f D gA (12.29)
4 1
En prenant en compte les conditions aux limites u1 D 0 et v1 D 0, l’expression W D 0 donne :
gL2
u2 D (12.30)
8E
1T
La déformation est donnée par " D u;x D T Bun avec T B D L . 1; 1/. On obtient :
gL
"D (12.31)
8E
Elle est constante sur l’élément. Quant à la contrainte, elle est obtenue par D E", soit :
gL
D (12.32)
8
Elle est également constante sur l’élément.
1 2 3 4
F IGURE 12.13: Discrétisation de la barre en trois éléments
longueur Le ) possède une matrice de rigidité élémentaire ke et un vecteur des forces nodales
élémentaires fe définis par :
2EA 1 1 Le 1
Ke D et fe D gA (12.33)
Le 1 1 4 1
K D K1 C K2 C K3 et F D f1 C f2 C f3 (12.34)
u1 et u,x1 u2 et u,x2
x
ξ
-1 0 1 x1 x2
(a) Élément de référence (b) Segment réel (2 nœuds, 4ddl)
1 1C
xD x1 C x2 ; 2 Œ 1I 1 (12.35)
2 2
Mais cette fois, l’approximation de la fonction solution est voulue sous la forme :
u D T .N1 ; N2 ; N3 ; N4 / un T
un D T u1 ; u;x1 ; u2 ; u;x2
avec les ddl (12.36)
Les fonctions Ni choisies sont des fonctions cubiques de type Hermite assurant une continuité de u
et u;x aux nœuds 1 et 2. On a :
1 L 2 1 L 2
N1 D .1 /2 .2C/I N2 D . 1/.1 /I N3 D .1C/2 .2 /I N2 D . 1/.1C/
4 8 4 8
(12.37)
Le champ u;x , qui est la déformation ", est donné par " D T Bfun g avec :
3 2 1 2 3 1 2
T
BD T
. 1/; .3 2 1/; .1 2 /; .3 C 2 1/ (12.38)
2L 4 2L 4
La contrainte, dans le cas où le coefficient d’élasticité H est constant est donnée par :
x D HT Bun (12.39)
III ÉLÉMENTS FINIS 12.5 Exemple : quelques variations sur le thème des éléments unidimensionnels 166
u1 u2
σ1 x
ξ
-1 0 1 x1 x2
(a) Élément de référence (b) Élément réel
1 1C
xD x1 C x2 ; 2 Œ 1I 1 (12.40)
2 2
avec :
Z L Z L Z C1
A 1 A N1 L L 1
bD I aD dxI fn D A fx d D Afx (12.44)
0 L 1 0 H 1 N2 2 2 1
On pourrait s’arrêter là avec le système précédent à résoudre. Toutefois, 1 est une variable locale
sans couplage avec les autres éléments. On peut donc l’exprimer en fonction des autres degrés de
liberté de l’élément u1 et u2 . La dernière ligne du système donne : 1 D T bun a1 D 0, 81 ,
d’où :
1T
1 D bun (12.45)
a
avec :
1T
KDb b (12.47)
a
Dans le cas où H est constant, la matrice de rigidité obtenue est identique à celle du modèle en
déplacements du paragraphe 12.5.5.
K D K1 C K2 C K3 et F D f1 C f2 C f3 (12.49)
On traitera le cas où les trois éléments ont la même longueurs (L1 D L2 D L3 D L) et où seul le
nœud 4 est soumis à une force. On obtient alors un système de type Kq D F à résoudre, avec la
condition aux limites u1 D 0, qui s’écrit explicitement dans ce cas :
2 30 1 0 1
1 1 0 0 u1 0
6 1 2 1 0 7 Bu2 C B 0 C
7 B C B
6
40 D C (12.50)
1 2 15 @u3 A @ 0 A
0 0 1 1 u4 F0
Ax D b (12.51)
A22 A21 A111 A12 x2 D b2 A21 A111 b1 que l’on note : Sx2 D c (12.54)
Kq D F (12.57)
Ce système est réduit, i.e. possède moins d’inconnues que le problème initial, mais il est nécessaire
de modifier le chargement extérieur.
Son inconvénient est qu’il est nécessaire d’effectuer un tri explicite au sein des degrés de liberté,
ce qui a pour conséquence de « remplir » le terme K12 , et par suite peut générer un surcoût de calcul
lorsque le nombre de degrés de liberté bloqués est grand.
C’est ainsi que procède le code A BAQUS pour imposer les déplacements.
Séparons maintenant les déplacements imposés des autres degrés de liberté, par une matrice
booléenne D : D est une matrice 1 4, et seul le terme D11 D 1 est non nul : q2 D Dq D u1 , et le
déplacement imposé est ud D 0. Le système s’écrit :
K11 K12 q1 f1
D (12.65)
K21 K22 q2 f2
avec :
2 2 3 2 33 0 0 11 0 0 1 1
2 1 0 1 u2 0
6K11 D 4 1 2 15 K12 D 4 0 57 Bq1 D @u3 AC B f1 D @ 0 A C
CDB
0
6 7B C
0 1 1 0 u 4 F
4 5@ A @ A
K21 D 1 0 0 K22 D 1 q2 D u1 f2 D ud D 0
Par la méthode du complément de Schur, on est encore ramené à la résolution du système (12.64).
On pourra s’amuser à calculer K111 , S... Par la méthode des multiplicateurs de Lagrange, on obtient
le système :
2 30 1 0 1 8
1 1 0 0 1 u1 0 ˆ u1 D 0
< u2 D F 0
ˆ
6 1 2 1 0 0 7 Bu2 C B 0 C
7 B C B C ˆ
ˆ
u3 D 2F0
6
60
6 1 2 1 0 7 Bu3 C D B 0 C
7 B C B C d’où : (12.66)
40 0 1 1 0 5 @u4 A @ F A 0 ˆ
u D 3F0
: 4
ˆ
ˆ
ˆ
1 0 0 0 0 ud D 0 D F0
1 2
3 4 5
de la poutre 1 2 et de la poutre 3 4 5, les points 2 et 3 étant astreints à rester collés ensembles.
En ajoutant la relation u2 D u3 par l’intermédiaire d’un multiplicateur de Lagrange, c’est-à-dire en
ajoutant le terme 1 .u2 u3 ), on obtient le système :
2 30 1 0 1
1 1 0 0 0 0 u1 0
6 1 1 0 0 0 17 7 B u2 C B 0 C
B C B C
6
60 0 1 1 0 17 Bu3 C B 0 C
7 B C B
B CDB 0 (12.70)
6 C
60 0 1 2 1 07 7 B u4 C B F C
C
6
40 0 0 1 1 0 5 @ u5 A @ 0 A
0 1 1 0 0 0 1 ud D 0
auquel il faut également ajouter la condition aux limites u1 =0. Si cette condition aux limites est
imposée par un multiplicateur de Lagrange 2 , alors finalement, il faut résoudre :
8
u1 D 0
2 30 1 0 1
1 1 0 0 0 0 1 u1 0 ˆ
0
ˆ
ˆ u2 D F0
6 1 1 0 0 0 1 0 7 B u2 C B 0 C ˆ
ˆ
ˆ
6 7B C B C
60 0 1 1 0 1 07 7 Bu3 C B 00 C < u3 D F
B C B C ˆ
ˆ
u4 D 2F0 (12.71)
6
60
6 0 1 2 1 0 0 u D
7 B 4C B C
7 B C B F C d’où :
60 0 0 1 1 0 07 7 B u5 C B 0 C u5 D 3F0
B C B C ˆ
ˆ
6 ˆ
ˆ
40 1 1 0 0 0 0 5 @1 A @ 0 A D F0
ˆ
: 1
ˆ
ˆ
ˆ
1 0 0 0 0 0 0 2 0 2 D F0
13.2 Post-traitement
Une manière d’améliorer les résultats est de recourir à des méthodes dites de post-traitement, i.e.
des méthodes qui, à partir des données issues de la résolution du système matriciel correspondant
au problème, fournissent des données complémentaires ou améliorent toute ou partie des données
déjà disponibles.
En fonction des problèmes (donc des données disponibles et des données souhaitées) de
nombreuses méthodes existent. Elles sont généralement dédiées à un problème donné.
173 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
Un exemple déjà abordé dans ce document est celui où, à partir des déplacement nodaux obtenus
par un calcul éléments finis d’une structure mécanique composée de deux matériaux différents dans
une discrétisation classique « en déplacements », on souhaite remonter aux contraintes à l’interface
entre lesdits matériaux.
Une méthode déjà mentionnée est la méthode dite de « Reissner local » qui consiste à intégrer les
équations d’équilibre sous la forme mixte de Reissner, mais uniquement sur des paires d’éléments
ayant des propriété matérielles différentes et possédant une face commune.
On obtient alors un champ de contrainte complètement continu dont on ne retient que les
composantes correspondant à la trace des contraintes, les autres étant calculées comme dans la
méthode classique. On montre que l’on améliore grandement la qualité de l’approximation des
contraintes aux interfaces, mêmes avec un maillage grossier.
Un autre exemple serait, disposant d’un point de pression constant par élément, de calculer
la pression en un point quelconque comme interpolation (linéaire ou non) à partir des points
disponibles.
Un troisième exemple serait à partir des données nodales dans un modèle unidimensionnel
(poutre ou barre selon le cas), de remonter à la répartition des contraintes en un points quelconque
de la structure (donc via les hypothèses faites dans le modèle sur la répartition des contraintes dans
une section + via une interpolation lorsque l’on se trouve dans une section ne passant pas par un
nœud)...
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 174
— partir des fichiers de résultats binaires pour les retravailler : c’est également théoriquement
faisable, mais il faut décortiquer les formats d’écriture desdits fichiers, les relire, réécrire
dans le même format... et c’est donc beaucoup de travail (purement informatique) ;
— programmer la fonction directement sous A NSYS (A NSYS possède un langage assez riche et
puissant) ;
— ou alors, et c’est la voie que nous allons montrer, utiliser A NSYS en « coopération » avec un
petit programme extérieur.
La structure de la macro A NSYS que nous proposons (et qui peut être également incluse dans les
menus d’A NSYS, mais nous n’entrons pas dans ce niveau de détail d’interfaçage) est très simple :
— on récupère les données dont nous avons besoin et on les stocke dans des fichiers externes
(simples) ;
— on lance un programme externe (dont la structure sera exposée plus bas) ;
— on réintègre les résultats dans A NSYS.
On peut alors se servir de toutes les fonctions de visualisation disponibles dans A NSYS avec les
résultats modifiés...
175 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
96 ! / s y s , ’ rm t e m p . dum ’ 111 *vread , CoordZ (1) , temp , sxy , ,
97 ! / s y s , ’ rm t e m p . t m p ’ 112 ( E13 .7)
98 ! 113 *vput , CoordX (1) , node ,1 ,S , x
99 ! l i r e l e s v a l e u r s dans temp . s i j 114 *vput , CoordY (1) , node ,1 ,S , y
100 ! ces f i c h i e r s c o n t i e n n e n t : sxx , syy , 115 *vput , CoordZ (1) , node ,1 ,S , xy
sxy 116 !
101 ! ce s o n t t o u t e s l e s composantes 117 ! Cleanup
continues . 118 CoordX (1) =
102 ! E l l e s n ’ o n t p a s t o u t e s un s e n s : 119 CoordY (1) =
103 ! Sxx d i s c o n t i n u , mais dans l e r e p e r e 120 CoordZ (1) =
104 ! l o c a l , d o n c on l a l a i s s e p o u r qu ’ 121 Mater (1 ,3) =
ANSYS 122 E1 (1 ,9) =
105 ! p u i s s e f a i r e l a ! r o t a t i o n de r e p e r e 123 Nnoeuds =
106 ! si necessaire . 124 Nmater =
107 *vread , CoordX (1) , temp , sxx , , 125 Nelem =
108 ( E13 .7) 126 !
109 *vread , CoordY (1) , temp , syy , , 127 /gop
110 ( E13 .7)
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 176
53 C Average multi computed nodes
54 CALL AVERAGE ( NODMAX , MULTI )
55 C
56 C Last line o f VM2
57 END
La routine READANS permet de lire les données. Celles-ci sont stockées dans les matrices
COORD ( NODMAX ,3), DISP ( NODMAX ,2), XMAT ( MATMAX ,3) et NLT ( NELTMAX ,9).
la routine CHELT détecte les interfaces, i.e. les faces dont les nœuds appartiennent à des éléments
dont les propriétés matérielles sont différentes. Le vecteur NADJ ( NODMAX ,8) contient les numéros
des deux éléments adjacents (positions 1 et 2) ainsi que les numéros des nœuds de l’interface
(positions 3, 4 et 5).
La routine MATER construit, pour chacune des faces de l’interface, la matrice de rigidité (aussi
bien pour le cas isotrope que pour le cas orthotrope).
La routine LOCREISS calcule les contraintes nodales le long de chaque face.
Les contraintes ayant été calculées pour les nœuds de chacune des faces de l’interface, on
moyenne les résultats pour un nœud appartenant à plusieurs faces. C’est ce que fait la routine
AVERAGE . À ce stade, nous disposons donc des contraintes nodales, calculées par la fonctionnelle
mixte d’Hellinger-Reissner, le long de toutes les interfaces présentes dans le modèle. Il ne reste
plus qu’à écrire ces résultats pour qu’ils soient relus par la macros A NSYS.
13.3.2 Poutre en U
R E F
C
A D
177 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
51 l ,4 ,5 , NCORE
1 ! geom : p a r t i e d r o i t e 52 l ,3 ,6 , NCORE
2 long =20 53 larc ,5 ,6 ,13 , Rh1 , NRCUTS
3 H1 =0.2 54 l ,5 ,8 , NSKIN
4 H2 =1.8 55 larc ,8 ,7 ,13 , Rint , NRCUTS
5 H3 =2 56 l ,6 ,7 , NSKIN
6 ! 57 !
7 ! geom c o u d e 58 l ,2 ,9 , NCUTS
8 Rint =10 59 l ,9 ,10 , NSKIN
9 Rh1 = Rint + H1 60 l ,10 ,3 , NCUTS
10 Rh2 = Rint + H2 61 l ,10 ,11 , NCORE
11 Rh3 = Rint + H3 62 l ,11 ,6 , NCUTS
12 ! 63 l ,11 ,12 , NSKIN
13 ! maillage 64 l ,12 ,7 , NCUTS
14 NCUTS =40 65 !
15 NRCUTS =20 66 mat ,1
16 NCORE =16 67 a ,1 ,2 ,3 ,4
17 NSKIN =2 68 a ,5 ,6 ,7 ,8
18 Forc =1 69 a ,2 ,9 ,10 ,3
19 ! 70 a ,6 ,11 ,12 ,7
20 /PREP7 71 mat ,2
21 ! 72 a ,4 ,3 ,6 ,5
22 ET ,1 , PLANE82 73 a ,3 ,10 ,11 ,6
23 ! 74 mat ,1
24 UIMP ,1 , EX , , ,70000 , 75 amesh ,1
25 UIMP ,1 , NUXY , , ,0.3 , 76 amesh ,2
26 UIMP ,1 , EMIS , , ,1 , 77 amesh ,3
27 ! 78 amesh ,4
28 ! 79 mat ,2
29 UIMP ,2 , EX , , ,3400 , 80 amesh ,5
30 UIMP ,2 , NUXY , , ,0.34 , 81 amesh ,6
31 UIMP ,2 , EMIS , , ,1 , 82 !
32 ! 83 DL ,4 , , symm
33 k ,1 ,0 , - Rh3 84 DL ,5 , , symm
34 k ,2 , Rh3 ,0 85 DL ,8 , , symm
35 k ,3 , Rh2 ,0 86 !
36 k ,4 ,0 , - Rh2 87 dk ,8 , uy ,0
37 k ,5 ,0 , - Rh1 88 !
38 k ,6 , Rh1 ,0 89 fk ,9 , fx , Forc
39 k ,7 , Rint ,0 90 !
40 k ,8 ,0 , - Rint 91 FINISH
41 k ,9 , Rh3 , long 92 /SOLU
42 k ,10 , Rh2 , long 93 /STAT , SOLU
43 k ,11 , Rh1 , long 94 SOLVE
44 k ,12 , Rint , long 95 FINISH
45 k ,13 ,0 ,0 96 /POST1
46 ! 97 pldisp ,1
47 larc ,1 ,2 ,13 , Rh3 , NRCUTS 98 *USE , INTERF
48 l ,2 ,3 , NSKIN 99 ...
49 larc ,4 ,3 ,13 , Rh2 , NRCUTS
50 l ,1 ,4 , NSKIN
On utilise la macro en tapant *USE,INTERF, et les contraintes sont directement modifiées dans
A NSYS. On peut alors utiliser les mêmes fonctions qu’habituellement pour visualiser les différentes
composantes des contraintes.
Le listing Cast3M du même problème est le suivant :
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 178
27 k2 = Rh3 0.; 69 PEAUINT = SURF1 ET SURF2 ;
28 k3 = Rh2 0.; 70 ELIMINE PEAUINT ;
29 k4 = 0. MRh2 ; 71 PoutU = PEAUEXT ET AME ET PEAUINT ;
30 k5 = 0. MRh1 ; 72
31 k6 = Rh1 0.; 73 Model1 = MODL peauext MECANIQUE
32 k7 = Rint 0.; ELASTIQUE ISOTROPE QUA4 ;
33 k8 = 0. MRint ; 74 Model2 = MODL ame MECANIQUE ELASTIQUE
34 k9 = Rh3 long ; ISOTROPE QUA4 ;
35 k10 = Rh2 long ; 75 Model3 = MODL peauint MECANIQUE
36 k11 = Rh1 long ; ELASTIQUE ISOTROPE QUA4 ;
37 k12 = Rint long ; 76 ModelTot = Model1 ET Model2 ET Model3 ;
38 k13 = 0. 0.; 77
39 78 Mater1 = MATERIAU Model1 YOUNG 70000.0
40 L1 = CERCLE NRCUTS k1 k13 k2 ; NU 0.3 RHO 2700.0;
41 L2 = DROITE k2 k3 NSKIN ; 79 Mater2 = MATERIAU Model2 YOUNG 3400.0 NU
42 L3 = CERCLE NRCUTS k4 k13 k3 ; 0.34 RHO 1000.0;
43 L4 = DROITE k1 k4 NSKIN ; 80 Mater3 = MATERIAU Model3 YOUNG 70000.0
44 L5 = DROITE k4 k5 NCORE ; NU 0.3 RHO 2700.0;
45 L6 = DROITE k3 k6 NCORE ; 81
46 L7 = CERCLE NRCUTS k5 k13 k6 ; 82 MR1 = RIGIDITE Model1 Mater1 ;
47 L8 = DROITE k5 k8 NSKIN ; 83 MR2 = RIGIDITE Model2 Mater2 ;
48 L9 = CERCLE NRCUTS k8 k13 k7 ; 84 MR3 = RIGIDITE Model3 Mater3 ;
49 L10 = DROITE k6 k7 NSKIN ; 85 Mrigid = MR1 ET MR2 ET MR3 ;
50 86
51 L11 = DROITE k2 k9 NCUTS ; 87 CondL2 = BLOQUER UX ( L4 ET L5 ET L8 ) ;
52 L12 = DROITE k9 k10 NSKIN ; 88 CondL1 = BLOQUER UY k8 ;
53 L13 = DROITE k10 k3 NCUTS ; 89 CondLtot = CondL1 ET CondL2 ;
54 L14 = DROITE k10 k11 NCORE ; 90 FOR1 = FORCE ( Forc 0.) k9 ;
55 L15 = DROITE k11 k6 NCUTS ; 91 Mrigid = Mrigid ET CondLtot ;
56 L16 = DROITE k11 k12 NSKIN ; 92
57 L17 = DROITE k12 k7 NCUTS ; 93 Depl1 = RESOUD Mrigid For1 ;
58 94
59 SURF1 = DALLER L1 L2 L3 L4 ; 95 UX1 = EXCO ’ UX ’ depl1 ;
60 SURF2 = DALLER L2 L11 L12 L13 ; 96 UY1 = EXCO ’ UY ’ depl1 ;
61 PEAUEXT = SURF1 ET SURF2 ; 97 TRACE UX1 PoutU ;
62 ELIMINE PEAUEXT ; 98 TRACE UY1 PoutU ;
63 SURF1 = DALLER L3 L6 L7 L5 ; 99
64 SURF2 = DALLER L6 L13 L14 L15 ; 100 def0 = DEFORMEE poutU Depl1 0.0 BLEU ;
65 AME = SURF1 ET SURF2 ; 101 def1 = DEFORMEE poutU Depl1 ROUG ;
66 ELIMINE AME ; 102 TRACE ( def0 ET def1 ) ;
67 SURF1 = DALLER L8 L7 L10 L9 ; 103
68 SURF2 = DALLER L10 L15 L16 L17 ;
Il reste maintenant à avoir les résultats dans les repères locaux élémentaires, et voir ce qui se
passe aux interfaces...
13.3.3 Résultats
Dans les listings ci-dessus, nous n’avons pas exposé comme atteindre les contraintes, car cela sera
fait en TP. Les résultats directement obtenus sont donnés à la figure 13.2 pour la déformée et les
déplacements selon x et y, et à la figure 13.3 pour les déformations.
Néanmoins nous présentons quelques résultats issus de l’analyse des contraintes aux interfaces
avec ou sans utilisation de la méthode de Reissner local. Nous nous concentrerons sur les résultats
numériques et graphiques permettant d’apprécier la précision de la méthode de Reissner local par
rapport aux résultats d’A NSYS aux points A, B, E et F ainsi que d’étudier l’influence de la valeur
du rayon R.
Une première remarque s’impose : plus la valeur du rayon est faible, plus les résultats sont
mauvais, quelle que soit la méthode employée (et vous devez savoir pourquoi à ce niveau du
document).
Au point A, l’influence des conditions aux limites est encore sensible : en plus de la condition
de symétrie de la structure, le blocage du déplacement vertical induit une légère détérioration des
résultats numériques. L’influence de cette condition aux limites n’est plus visible au point B.
En A et B, la composante discontinue est xx , les composantes yy et xy sont continues.
Par contre, aux points E et F, c’est yy la composante discontinue et xx et xy les composantes
179 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
Déplacement Déplacement
UX UY
0.59 −7.26E−04
0.56 −5.75E−03
0.53 −1.08E−02
0.50 −1.58E−02
0.47 −2.08E−02
0.45 −2.58E−02
0.42 −3.09E−02
0.39 −3.59E−02
0.36 −4.09E−02
0.33 −4.59E−02
0.30 −5.09E−02
0.28 −5.60E−02
0.25 −6.10E−02
0.22 −6.60E−02
0.19 −7.10E−02
0.16 −7.60E−02
0.13 −8.11E−02
0.11 −8.61E−02
7.76E02 −9.11E−02
4.93E02 −9.61E−02
2.10E02 −0.10
continues.
D’après les résultats présentés ci-après, il est visible que la méthode de Reissner local permet
d’améliorer les résultats par rapport à A NSYS.
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 180
xx xx
Méthode R peau âme yy xy
(mm) (MPa) (MPa) (MPa) (MPa)
Ref 10 72,085 4,0121 1,5883 -0,0023339
Reiss 10 71,957 3,9670 1,5987 -0,0155720
A NSYS 10 71,957 3,9670 1,5116 0,0086231
Ref 8 66,500 3,8400 1,8556 -0,0022273
Reiss 8 66,665 3,7863 1,8888 -0,0118667
A NSYS 8 66,665 3,7863 1,7943 0,0111910
Ref 5 57,289 3,6598 2,6317 -0,0018770
Reiss 5 57,727 3,5771 2,6979 -0,0071024
A NSYS 5 57,727 3,5771 2,6066 0,0128210
Ref 3 49,551 3,6827 3,9304 -0,0014368
Reiss 3 49,454 3,5552 4,0110 -0,0046055
A NSYS 3 49,454 3,5552 3,9621 0,0106410
Ref 2 42,699 3,8289 5,4227 -0,0009727
Reiss 2 42,362 3,6593 5,4785 -0,0037443
A NSYS 2 42,362 3,6593 5,5264 0,0059470
Ref 1 27,171 4,3142 9,2850 -0,0004173
Reiss 1 25,876 4,0707 9,1122 -0,0021565
A NSYS 1 25,876 4,0707 9,6345 0,0005231
8 Ref
Reissner local C
7 A NSYS
6
yy (MPa)
C
5
4 C
3
C
2 C
C
1
1 2 3 4 5 6 7 8 9 10
R
181 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
0:015
0:010
0:005 Ref
Reissner local C
xy (MPa)
0:000 A NSYS
C
C C
0:005
C
0:010
C
0:015 C
0:020
1 2 3 4 5 6 7 8 9 10
R
xx xx
Méthode R peau âme yy xy
(mm) (MPa) (MPa) (MPa) (MPa)
Ref 10 -72,248 -3,0668 1,4080 -0,0023570
Reiss 10 -72,542 -3,1193 1,2573 -0,0138800
A NSYS 10 -72,542 -3,1193 1,5464 -0,0198700
Ref 8 -67,723 -2,7905 1,5799 -0,0026011
Reiss 8 -68,001 -2,8457 1,4570 -0,0107700
A NSYS 8 -68,001 -2,8457 1,7111 -0,0197300
Ref 5 -61,079 -2,3199 2,0320 -0,0024099
Reiss 5 -61,335 -2,3813 1,9460 -0,0067001
A NSYS 5 -61,335 -2,3813 2,1535 -0,0186210
Ref 3 -56,823 -1,9080 2,6611 -0,0022154
Reiss 3 -57,072 -1,9756 2,5888 -0,0045050
A NSYS 3 -57,072 -1,9756 2,7738 -0,0176540
Ref 2 -54,772 -1,6200 3,2396 -0,0021234
Reiss 2 -55,026 -1,6902 3,1612 -0,0036124
A NSYS 2 -55,026 -1,6902 3,3438 -0,0172760
Ref 1 -52,491 -1,1687 4,2871 -0,0020797
Reiss 1 -52,777 -1,2391 4,1625 -0,0031088
A NSYS 1 -52,777 -1,2391 4,3710 -0,0174200
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 182
4:5
C
4:0 Ref
Reissner local C
3:5 A NSYS
C
yy (MPa)
3:0
2:5 C
2:0 C
1:5 C
C
1:0
1 2 3 4 5 6 7 8 9 10
R
0:002
C C
0:004 C Ref
Reissner local C
0:006 A NSYS
C
0:008
xy (MPa)
0:010
C
0:012
0:014 C
0:016
0:018
0:020
1 2 3 4 5 6 7 8 9 10
R
183 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
yy yy
Méthode R peau âme xx xy
(mm) (MPa) (MPa) (MPa) (MPa)
ref 10 46,0340 2,410200 0,57548 0,77518
Reiss 10 46,0265 2,397450 0,48986 0,71342
A NSYS 10 46,0265 2,397450 0,62529 0,80572
ref 8 45,8625 2,443000 0,68465 0,83563
Reiss 8 45,8505 2,421450 0,62992 0,76781
A NSYS 8 45,8505 2,421450 0,71910 0,90760
ref 5 45,4490 2,539450 1,00620 1,00390
Reiss 5 45,5025 2,499550 1,00810 0,91559
A NSYS 5 45,5025 2,499550 0,94099 1,19800
ref 3 44,5760 2,694400 1,56820 1,25790
Reiss 3 45,0015 2,644550 1,62260 1,10690
A NSYS 3 45,0015 2,644550 1,30760 1,64980
ref 2 43,1510 2,856600 2,23920 1,50420
Reiss 2 44,1375 2,805000 2,32600 1,25870
A NSYS 2 44,1375 2,805000 1,77080 2,12250
ref 1 37,6330 3,184100 4,02110 1,90380
Reiss 1 40,1735 3,143850 4,04390 1,41220
A NSYS 1 40,1735 3,143850 3,04360 3,15350
4:5
4:0C Ref
Reissner local C
3:5 A NSYS
3:0
yy (MPa)
2:5
C
2:0
C
1:5
1:0 C
C
0:5 C
0:0
1 2 3 4 5 6 7 8 9 10
R
III ÉLÉMENTS FINIS 13.3 Exemple d’implémentation d’un post-traitement dans A NSYS 184
3:5
3:0 Ref
Reissner local C
A NSYS
2:5
xy (MPa)
2:0
1:5C
C
C
1:0 C
C
C
0:5
1 2 3 4 5 6 7 8 9 10
R
yy yy
Méthode R peau âme xx xy
(mm) (MPa) (MPa) (MPa) (MPa)
ref 10 -47,438 -2,1425 0,54431 0,30073
Reiss 10 -47,510 -2,1533 0,39416 0.38936
A NSYS 10 -47,510 -2,1533 0,65814 0,30000
ref 8 -47,4555 -2,11415 0,62480 0,24743
Reiss 8 -47,5605 -2,12825 0,50447 0.33138
A NSYS 8 -47,5605 -2,12825 0,72921 0,21663
ref 5 -47,415 -2,03645 0,83167 0,10558
Reiss 5 -47,5335 -2,05940 0,76054 0,18360
A NSYS 5 -47,5335 -2,05940 0,87666 0,00251
ref 3 -47,223 -1,92305 1,11940 -0,10042
Reiss 3 -47,1885 -1,95170 1,07760 -0,01314
A NSYS 3 -47,1885 -1,95170 1,05090 -0,28747
ref 2 -46,940 -1,81445 1,38480 -0,29901
Reiss 2 -46,6760 -1,84130 1,35230 -0,18934
A NSYS 2 -46,6760 -1,84130 1,21420 -0,55238
ref 1 -46,1935 -1,61030 1,85860 -0,67358
Reiss 1 -45,4620 -1,62185 1,81270 -0,50497
A NSYS 1 -45,4620 -1,62185 1,52280 -1,03400
185 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
2
1:8C Ref
Reissner local C
1:6 A NSYS
1:4 C
yy (MPa)
1:2
C
1
0:8 C
0:6
C
0:4 C
0:2
1 2 3 4 5 6 7 8 9 10
R
0:400 C
C
0:200 C
0:000 C
Ref
Reissner local C
0:200 C
xy (MPa)
A NSYS
0:400
C
0:600
0:800
1:000
1:200
1 2 3 4 5 6 7 8 9 10
R
La simulation est donc décomposée en différents niveaux, chacun représentant une échelle
différente. Pour coordonner ces niveaux entre eux, on utilise généralement une approche dite
descendante : on commence par simuler le comportement global de l’avion, puis les résultats sont
utilisés pour déterminer les conditions aux limites appliquées sur les niveaux inférieurs. Toutefois,
il est parfois également nécessaire de combiner ces approches descendantes à des approches
ascendantes : les résultats des simulations fines sont utilisés pour construire des modèles de
comportements plus grossiers.
De très nombreuses méthodologies multi-échelles existent et reposent toutes sur le fait de
simuler chaque phénomène à l’échelle la plus pertinente. Pour cela, il est nécessaire :
1. de distinguer différentes échelles dans la modélisation et dans la simulation ;
2. de modéliser les relations existant entre ces différentes échelles.
Considérons un problème comportant deux échelles. Le point 1 évoqué ci-dessus se traduit par
le fait de disposer de deux modèles distincts :
— Un modèle macroscopique
Il représente le produit et son environnement extérieur et est constitué d’une géométrie et d’un
modèle de comportement relativement grossiers (puisque n’ayant pas vocation à représenter
les phénomènes microscopiques) ;
— Un modèle microscopique
Il possède une description géométrique et un maillage suffisamment fins ainsi que soit
modélisé le comportement détaillé. Il ne comprend qu’une ou quelques « cellules types » de
petites dimensions.
Il est souvent possible, pour les matériaux notamment, de disposer d’hypothèses simplifica-
trices (telles que l’élasticité linéaire) permettant de se restreindre à une seule cellule type.
Pour le cas des matériaux anisotropes constitués de plusieurs autres matériaux tels que les
matériaux composites et fibreux, on se reportera au chapitre 15.
Toutefois, il peut arriver que l’on ne dispose pas de telles hypothèses et il soit alors nécessaire
de multiplier le nombre de cellules (jusqu’à parfois couvrir tout le modèle macro !)
Le point 2 consiste donc à relier les deux échelles de modélisation.
En effet, les modèles macroscopique et microscopique ne sont pas indépendants puisqu’ils
modélisent la même physique à des échelles différentes. Il est donc nécessaire qu’ils soient cohérents
l’un vis-à-vis de l’autre, en tout point et à chaque instant de la simulation. Pour cette raison, les
modélisations multi-échelles comportent des couplages, i.e. des modèles d’interactions, entre les
187 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
échelles.
Dans le cadre de la mécanique des solides déformables, on procède généralement comme suit :
— Le modèle de comportement macroscopique : qui modélise des phénomènes se produisant
« au cœur du matériau », doit correspondre à la relation contraintes/déformations observée
sur la cellule microscopique ;
— Le modèle de comportement microscopique : qui traduisent la façon dont la cellule est
sollicitée par son environnement extérieur, doivent correspondre à l’état de contraintes ou de
déformations macroscopique.
Le modèle macroscopique ayant une résolution beaucoup plus grossière, la cohérence des
deux modèles ne peut donc pas se traduire par une correspondance exacte, point par point, des
conditions aux limites ou du comportement. Pour cette raison, la plupart des approches multi-
échelles postulent que les quantités macroscopiques doivent correspondre à des « moyennes » des
quantités microscopiques correspondantes, la définition mathématique exacte de cette moyenne
variant fortement d’une approche à l’autre.
On peut donc dire que :
— Le modèle de comportement macroscopique d’un élément de volume est choisi égal au
comportement moyen de la cellule microscopique correspondante, calculé en utilisant une
technique d’homogénéisation ;
— Les conditions aux limites microscopiques sont appliquées en moyenne, car les champs de
contraintes ou de déplacements macroscopiques sont beaucoup trop grossiers.
Le problème est donc d’échanger des données pertinentes entre les différentes échelles, compte
tenu des différentes résolutions des modèles.
Dans le sens microscopique vers macroscopique, on utilisera les techniques d’homogénéisation
dont il sera l’objet au chapitre 15.
Dans le sens macroscopique vers microscopique, il s’agit donc de spécifier des conditions aux
limites à appliquer sur le bord des cellules microscopiques, à partir d’une solution macroscopique.
Les méthodes les plus simples se contentent d’imposer directement le champ de déplacements
(ou de contraintes) macroscopique comme condition aux limites. Le modèle macroscopique ayant
par définition une résolution beaucoup plus grossière que le modèle microscopique, il est incapable
de capturer l’allure microscopique des déplacements ou des contraintes. Ces conditions aux limites
sont donc souvent trop imprécises par rapport aux finalités du modèle microscopique. Cela peut
dégrader fortement la qualité des résultats et donc restreindre le domaine de validité de ces
approches.
Pour palier ce problème, on écrit les conditions aux limites microscopiques de manière plus
subtile en les décomposant en la somme d’un terme moyen et d’un terme de moyenne nulle. Ainsi,
le champ microscopique um sera écrit um D uM C v m où uM est le champ macroscopique (i.e.
la moyenne du champ microscopique), et v m est le reste (à moyenne nulle). Notre problème
devient donc de déterminer le reste v m , et c’est là que les méthodes divergent le plus. On distingue
néanmoins deux grosses familles :
— Condition de périodicité
On postule que l’on connaît l’allure du reste. On peut alors enchaîner un calcul macroscopique
avec un solveur spécifique puis un calcul microscopique avec un solveur classique ;
— Couplage des cellules microscopiques
Lorsqu’il n’est pas si simple de « séparer » les échelles, alors il faut résoudre en même temps
les deux problèmes microscopique et macroscopique avec échange de données.
Ces deux approches vont être un peu plus détaillées maintenant.
189 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
de propager une information fine sur l’ensemble de la pièce et, ainsi, de se passer de l’hypothèse
de séparation des échelles : il n’est plus nécessaire de modéliser séparément les phénomènes
microscopiques et macroscopiques.
Concrètement, les méthodes de décomposition de domaine sont des solveurs qui ont générale-
ment un fonctionnement multi-échelles. Elles partent d’un modèle microscopique du produit décom-
posé en sous-structures, et consistent à coupler les sous-structures en échangeant des contraintes et
des déplacements sur les interfaces. Des versions plus ou moins simples existent. Notons bien qu’il
n’y a pas de modèle macroscopique dans une telle approche.
En résumé :
Dans la décomposition de domaine, la simulation est abordée à l’échelle la plus fine : la
sous-structuration et le problème grossier ne sont utilisés que pour améliorer l’efficacité de la
résolution.
Une simulation par décomposition de domaine conduit toujours à un champ de déplacement
microscopique continu sur toute la structure, et un champ de contraintes microscopique équilibré
(au sens des éléments finis) sur toute la structure.
Actuellement, la simulation multi-échelles est encore relativement peu répandue dans le monde
de l’ingénierie. Elle représente en effet un changement considérable par rapport aux pratiques
usuelles de simulation ; de plus, la plupart des logiciels de calcul multi-échelles sont des outils
développés par des chercheurs, qui n’ont pas encore l’ergonomie et la robustesse des solveurs
utilisés dans l’industrie. Les industriels attendent donc l’apparition d’outils mieux adaptés à leurs
problématiques avant d’envisager une utilisation plus fréquente de ces méthodes.
Cependant, à plus long terme, la simulation multi-échelles suscite un intérêt considérable dans
l’industrie : en rendant accessibles à la simulation des phénomènes qui ne peuvent actuellement être
étudiés qu’expérimentalement, elle constitue un pas important en direction du « virtual testing ».
Pour cette raison, elle fait toujours l’objet de nombreux projets de recherche. Ceux-ci concernent
aussi bien la modélisation, avec notamment la mise au point de « matériaux virtuels » (qui ne sont
rien d’autre que des modèles multi-échelles, notamment pour les matériaux composites), que les
solveurs qui sont en constante évolution.
13.5 Super-éléments
Histoire
Originellement introduit dans l’aéronautique dans les années 60 (d’où le choix de la figure 13.12), le
concept de sous-structuration répondait à trois motivations :
— faciliter la division du travail : des sous-structures avec des fonctions différentes (fuselage,
ailes, train d’atterrissage...) pouvaient être traitées par des groupes d’experts différents.
Chaque groupe pouvait à loisir améliorer, raffiner... sa partie tant que l’interface avec les
autres parties restait inchangée.
— profiter de la répétition : en remarquant qu’une même structure peut contenir plusieurs sous-
structures identiques, il est possible de diminuer le temps d’étude (par exemple symétrie des
ailes...)
— contourner les limitations des ordinateurs : les ordinateurs de l’époque atteignaient vite
leurs limites (par exemple en terme de taille mémoire). Diviser une structure complexe, que
l’on était incapable de calculer en une seule fois, permettait de sauvegarder des résultats
sous-structure par sous-structure puis d’effectuer « l’assemblage » des résultats.
Si les deux premiers points sont toujours d’actualité, le troisième l’est moins, surtout depuis le
recourt aux algorithmes parallèles.
C’est d’ailleurs le développement de procédures pour le calcul parallèle qui a conduit les
mathématiciens appliqués au concept de sous-domaines, alors qu’ils devaient grouper des éléments
pour des raisons de calcul.
Dans ce paragraphe, nous allons parler du concept de super-élément, qui n’est qu’une application
de ce qui a été présenté au paragraphe précédent.
191 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
13.5.2 Remonter aux degrés de liberté internes
Une fois le système condensé résolu, on obtient les valeurs aux nœuds internes en réutilisant la
formule :
1
qi D Kii .fi Kib qb / (13.6)
1T T
MD qKq qF (13.10)
2
Notons que l’on peut écrire les p conditions cinématiques sous la forme : Cq D ı avec C une
matrice p n, et q et ı des vecteurs.
Nous avons également déjà vu que ces conditions peuvent être prises en compte par l’inter-
médiaire de p multiplicateurs de Lagrange dans la fonctionnelle précédente qui devient alors :
1T
M D T T
qKq qF Cq ı (13.11)
2
Le système à résoudre est symétrique, régulier, mais non défini-positif.
Pour pallier la lenteur des algorithmes disponibles pour la résolution numérique d’un tel système
(méthode de Gauß ou de décomposition avec pivotage), il est préférable de régulariser le système.
Pour cela, on s’appuie sur la connaissance de la matrice K : celle-ci est singulière d’ordre r,
où r correspond aux mouvements de corps rigide, ou modes rigides. Cela veut également dire
qu’elle possède .n r/ valeurs propres strictement positives (et r valeurs propres nulles).
Il est évident qu’il faut au moins disposer de p > n conditions aux limites cinématiques pour
que le système puisse admettre une solution. C’est évidemment l’hypothèse que nous ferons (sinon
le problème est mal posé).
Nous considérons la matrice R de dimension n r des r vecteurs propres de K correspondant
à la valeur propre nulle. Quitte à construire ces vecteurs (qui sont orthogonaux), nous les choisirons
normés. On a alors :
(
KR D 0
T
(13.12)
RR D Ir
On introduit la matrice :
K˛ D K C ˛RT R (13.13)
qui, pour ˛ > 0 est régulière car symétrique et définie-positive. Les colonnes de R sont vecteurs
propres de K˛ pour la valeur propre ˛.
Dans le système Kq D F, on remplace K par K˛ RT R, et en introduisant le fait que K˛ RT R D
T
˛R R, on obtient finalement :
RT R q D F
K˛ In (13.14)
K˛ v D F (13.15)
193 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
On peut remarquer que K˛ 1 , que l’on note S˛ , possède les mêmes valeurs propres que K˛ (et
donc les mêmes que K et RT R). On est donc naturellement amené à décomposer S˛ comme :
1 T
S˛ D S C R R (13.17)
˛
avec SR D 0 et T RS D 0.
C’est cette matrice S que l’on appelle quasi-inverse de K, car :
SK D KS D I RT R (13.18)
Elle est de dimension n n, symétrique, semi-définie positive (ses valeurs propres sont de même
signe), admet les mêmes valeurs propres que K et en particulier ceux de la valeur propre nulle
dont R est une base, mais elle ne possède pas de caractère bande.
Voici brossé, en quelques lignes, les idées principales de la méthode. Nous n’irons pas plus loin
dans sa présentation.
On rappelle que l’intérêt de la méthode est qu’une fois une première étape consistant à intégrer
les données relatives à la structure (géométrie, discrétisation, matériaux) réalisée, on peut alors
effectuer autant de fois que nécessaire la seconde étape qui porte sur la prise en compte des
chargements et conditions aux limites.
Notons qu’il est possible de modifier un peu la forme du système secondaire afin de pouvoir
prendre en compte, lors de la seconde étape, des conditions plus complexes, comme des conditions
mixtes par exemple...
Nous considérons la fonction coût J.; u / choisie pour décrire le problème d’optimisation. Cette
fonction dépend donc naturellement du domaine et de la solution du problème sur ce domaine u .
Nous allons donner un sens à la dérivée de la fonction
R coût par rapport aux variations du domaine.
Si l’on considère le cas simple J.; u/ D u, alors il vient :
Z Z
dJ
.; uv ; V/ D u V n C u0Iv (13.21)
d
avec :
uCt V ı .I C tV/.x/ u .x/
uP IV D lim ; 8x 2 (13.24)
t !0 t
On définit alors la dérivée de la fonction coût par rapport aux variations du domaine par :
@J
Z
.; u; V/ D uV n (13.25)
@
La difficulté tient à ce que l’ensemble des domaines ne constitue pas un espace vectoriel. En
effet, les perturbations ne s’ajoutent pas, ou au moins ne sont pas associatives : . C V/ C W ¤
. C W/ C V.
Toutefois, si l’on utilise un paramétrage du domaine, il est alors possible d’utiliser les outils
classiques du calcul différentiel dans les espaces normés. On introduit alors le paramètre F par :
1
J.F; u/ D J.F./; u ı F /; F 2 W1;1 .I Rn / (13.26)
et l’on a alors .F C V/ C W D .F C W/ C V.
195 13. Calcul efficient : qualité des résultats et efforts de calcul ÉLÉMENTS FINIS III
La dérivée partielle de J par rapport au paramètre F est définie sans ambiguïté. Dans notre
exemple, cette dérivée, prise au point F D I est :
@J
Z
.I; u/ V D u div V (13.27)
@F
@J @J
.; u; V/ ¤ .I; u/ V (13.28)
@ @F
ce à quoi nous répondrons que c’est le prix à payer pour se ramener à un espace vectoriel. Cette
démarche peut alors être généralisée aux ordres supérieurs.
Le maillage
Le maillage n’a pas seulement un intérêt en calcul scientifique. Le maillage est un « support » de
représentation tridimensionnel utilisé par exemple dans les jeux vidéo ainsi que dans les animations
3D. Dans ce dernier cas, on s’intéresse à la qualité du maillage en lien avec la qualité de rendu
de l’image générée lorsque l’on applique une texture sur le maillage. On s’intéresse également à
comment « faire bouger » le maillage pour qu’un personnage n’apparaissent pas distordu pendant
une animation...
L’opération de maillage peut se faire à partir de plusieurs données :
— soit à partir de données à utiliser pour la discrétisation : sommets, arêtes...
— soit à partir de données de type CAO décrivant uniquement les entités géométriques.
Nous nous contenterons de présenter quelques outils de maillage relatifs au premier cas, i.e. lorsque
nous disposons de points, arêtes... mais le second cas n’est pas vraiment plus compliqué.
Les méthodes de construction de maillage sont essentiellement :
— Maillage par triangles/tétraèdres :
— Méthodes utilisant le critère de Delaunay : les bases de la méthode seront exposées au
paragraphe 14.1 ;
— Méthodes par avancement de fronts : le principe en sera détaillé au paragraphe 14.2 ;
— Méthodes par décomposition spatiale ;
— Maillage par quadrangles/hexaèdres : quelques remarques seront faites au paragraphe 14.4
— Maillage par avancement de fronts ;
— Maillage par décomposition en domaines.
Parfois, il ne s’agit pas de mailler un domaine, mais de le remailler. Les méthodes d’adaptation
de maillage les plus connues sont :
— Raffinement ;
— Transformation des éléments ;
— Déplacements de nœuds ;
— Simplification de maillage.
L’utilisation de différentes méthodes est illustré à la figure 14.1 issue de [2].
i.e. une arête ai est définie par ses deux sommets .s1i ; s2i / qui sont des points de S, et un sous-
domaine est défini par une arête frontière ai et un sens de parcourt (positif ou négatif).
Définition 65 Les diagrammes de Voronoï sont les polygones convexes Vi , i D 1; :::; Np formés
par l’ensemble des points de R2 plus proches de x i que des autres points x j , soit :
n o
Vi D x 2 R2 =kx x i k 6 kx x j k; 8j 2 f1; :::; Np g (14.4)
Chaque point x i (point générateur) étant considéré comme une île d’où partent des bateaux, la
région de Voronoï du point x i est la région où un bateau issu de l’île x i arrive avant tout bateau issu
d’une autre île. Notons qu’il est possible de jouer avec la métrique pour définir des diagrammes
de Voronoï sur d’autres géométries que la géométrie euclidienne, et qu’il est possible également
d’étendre la définition au cas où « les bateaux ne vont pas tous à la même vitesse ». Dans ce dernier
cas, on parle de diagramme de Voronoï pondéré, que nous n’utiliserons pas dans ce document.
Les diagrammes de Voronoï sont des polygones obtenus comme intersections finies de demi-
espaces et sont donc convexes. De plus, les sommets v k de ces polygones sont à égale distance
k
des points fx ij ; j D 1; :::; nk g de S, où le nombre nk est généralement égal ou supérieur à 3.
À chacun de ces sommets v k , nous pouvons associer le polygone convexe construit avec les
k
points fx ij ; j D 1; :::; nk g en tournant dans le sens trigonométrique. Ce maillage est généralement
formé de triangles, sauf si il y a des points cocycliques.
Le domaine d’un maillage de Delaunay d’un ensemble de points S est l’intérieur du convexi-
fié C.S/ de l’ensemble de points S.
Il existe une propriété qui permet de savoir si un maillage est un maillage de Delaunay. C’est la
propriété de la boule ouverte.
La propriété de la boule ouverte peut être appliquée également au cas d’un quadrilatère convexe,
puisque celui-ci peut être découpé en deux triangles adjacents. Il suffit alors de vérifier la propriété
pour toutes les paires de triangles adjacents formant un quadrilatère convexe.
Si le critère de la boule ouverte n’est pas vérifié pour un quadrilatère convexe, alors on fera un
échange de diagonal pour résoudre le problème.
14.1.3 Remarques
À ce niveau, nous sommes en mesure de générer un maillage... toutefois, celui-ci n’est pas forcément
exempt de problèmes :
— si les seuls points disponibles sont sur la frontière, il va falloir générer des points internes au
maillage. En effet, le maillage généré existe, mais peut présenter des distorsions inadmissibles
en termes de calcul. Le critère le plus naturel pour distribuer les points internes est d’imposer
en tout point x de R2 , le pas de maillage h.x/. En pratique, on ne dispose pas toujours de
cette information et il faut la construire à partir des points disponibles, i.e. des points de la
frontière. D’autre part, dans de nombreuses applications, on préfère donner le nombre de
F IGURE 14.5: a) Points et arêtes, b) maillage de Delaunay ne respectant pas la frontière, c) maillage
respectant la frontière obtenu en échangeant des diagonales
III ÉLÉMENTS FINIS 14.4 Remarques sur le maillage quadrangulaire et hexaédrique 204
Chapitre 15
Homogénéisation
Les théories des milieux effectifs visent à estimer les propriétés effectives (i.e. macroscopiques) d’un
milieu en fonction des propriétés locales de chacun des constituants ainsi que d’un certain nombre
d’informations sur la microstructure.
Les premières théories remontent au XIXe
siècle et sont dues à Mossotti, Maxwell, Pois-
son ou encore Lorentz. Le but de ces modèles
est de fournir soit des bornes pour le com-
portement effectif, soit des approximations
du comportement effectif. Les bornes sont
optimales lorsqu’il existe une microstructure Mossotti Maxwell Poisson Lorentz
particulière qui réalise exactement le modèle
physique. On évalue les « bonnes » propriétés de ces théories en les confrontant à d’autres résultats
théoriques, calculs analytiques et numériques...
Ces théories sont utilisées pour les problèmes de conductivité (milieux diélectriques), en méca-
nique, magnétique, thermique... lorsque l’on a des phases de conductivité, d’élasticité, des coeffi-
cients thermiques... variables. Ces problèmes sont en général très difficiles à résoudre (non-linéaires
et anisotropes) alors qu’en même temps, du point de vue des applications pratiques, il n’est pas
forcément nécessaire de tenir compte de l’ensemble des degrés de liberté de ces systèmes.
L’existence d’un comportement effectif n’est nullement assurée. On montre que, sous certaines
hypothèses (en particulier l’existence d’un volume élémentaire représentatif ), on peut effectivement
remplacer un matériau hétérogène par un milieu homogène équivalent.
Enfin, d’un point de vue purement numérique, ces méthodes peuvent être utilisées pour « simpli-
fier » un système à résoudre, que cela ait un sens physique ou non.
Notons que les techniques d’homogénéisation ne sont pas que des artefacts, des trucs et astuces.
Certaines des grandeurs physiques que nous utilisons tous les jours ne sont que des moyennes.
Le meilleur exemple en est la pression. Bien qu’en un point donné d’un gaz on voit passer des
particules allant en tous sens, on constate également, à une échelle plus macroscopique, un schéma
d’ensemble qui permet de définir par exemple la pression qu’exerce ledit gaz sur une paroi...
pourtant rien de « cohérent » ne se dégage à l’échelle microscopique.
On rappelle également que si est un ouvert connexe borné de Rn , et H10 ./ l’espace de Sobolev
des fonctions qui s’annulent sur , alors la solution u 2 H10 ./ du problème précédent est
également la solution du problème faible :
Z Z
1
8v 2 H0 ./; A.x/ru rv D f .x/v (15.2)
où f 2 L2 ./. Pour que ces deux formulation soient équivalentes, on a supposé que A.x/ est
régulière : i.e. que pour chaque x de , A.x/ est une matrice n n, dont les éléments sont
mesurables et qui est symétrique, définie-positive, bornée, i.e. qu’il existe deux constantes K1
et K2 > 0 telles que :
8v 2 Rn ; K1 v v 6 A.x/v v 6 K2 v v (15.3)
Mentionnons le cas particulier où A.x/ est régulier dans mais pas sur , alors on a deux
conditions d’interface :
Nous savons également que la solution existe et est unique (d’après le théorème de Lax-Milgram,
i.e. d’après le théorème de Riesz-Fréchet en remarquant le produit scalaire qui va bien...)
Intéressons nous maintenant d’un peu plus près à A.x/. Supposons que A.x/ possède une
certaine périodicité ", par exemple, que le matériau constituant le domaine soit constitué de couches
successives de deux matériaux différents, mais avec une périodicité, comme illustré sur la figure 15.1.
Dans un tel cas, nous dirons que les coefficients A dépendent de x=". Cette écriture sous forme
y
0 ε 2ε 3ε x
F IGURE 15.1: Matériau périodique
< d A x
8 du
D f .x/ pour x 2 Œ0; l
dx " dx (15.5)
uD0 pour x D 0 et x D l
:
A est une fonction de R dans R de période 1, et nous la supposerons régulière par morceaux et telle
qu’il existe K1 et K2 > 0 telles que K1 6 A.v/ 6 K2 .
La solution du problème est alors :
Z x Z y
1
u.x/ D A .y="/ f .v/dv C c1 dy C c2 (15.6)
0 0
avec c2 D 0 et
Z l Z y
1
A .y="/ f .v/dvdy
0 0
c1 D Z l
(15.7)
1
a .y="/dy
0
Le nombre d’intervalles pour le calcul approché de c1 est très grand (i.e. très supérieur à l=").
Exprimé autrement, ce résultat devient : le coefficient homogénéisé est 1=h A1 i et non pas hAi
(où hi désigne la moyenne).
où p 2 C.Œ0 I 1/ ; p.x/ > 0; x 2 Œ0 I 1 ; f 2 C.Œ0 I 1/ ; A 2 R, B 2 R et " > 0 le petit paramètre.
En négligeant le terme "2 u00 , nous arrivons à l’équation p.x/ur D f .x/. En d’autres termes,
et comme montré sur la figure 15.2, nous disposons de la courbe correspondant à y D ur D
f .x/=p.x/ (en noir) qui approche la solution du problème donnée par la courbe rouge.
La solution exacte se confond avec ur sur la plus grande partie de l’intervalle Œ0 I 1, mais u
et ur sont fortement différents dans un voisinage des extrémités.
où les uri .x/ sont les termes réguliers, et u0i .x="/ et u1i ..x 1/="/ les termes de couche limite
correspondant aux extrémités x D 0 et x D l.
Pour le calcul on procède donc comme suit :
— on injecte les termes réguliers de l’ansatz dans la formulation du problème ;
— on traite les extrémités pour vérifier les conditions aux limites, ce qui donne les termes de
couche limite ;
— on vérifie que la solution trouvés conviennent bien (i.e. que l’erreur commise est négligeable).
15.2.1 Introduction
D’une manière générale, tous les matériaux, mêmes isotropes, sont hétérogènes en dessous d’une
certaine échelle. Il peut sembler naturel d’utiliser des propriétés homogènes équivalentes correspon-
dant à des propriétés mécaniques effectives. Toutefois, comme nous l’avons déjà vu, ces propriétés
effectives ne s’obtiennent pas par une simple moyenne des propriétés des constituants, mêmes
pondérées par les fractions volumiques.
Les propriétés effectives du milieu homogène équivalent cherché peuvent être obtenues en
résolvant un problème aux limites sur un volume élémentaire d V, à condition que celui-ci soit
suffisamment grand pour être représentatif de la microstructure du matériau hétérogène. Dans le cas
où les constituants présentent une structure périodique, le volume d V peut se réduire à un volume
élémentaire.
Une fois le volume élémentaire déterminé, on le soumet à des sollicitations élémentaires pour
déterminer la réponse résultante. La difficulté réside en fait dans le choix des conditions aux limites
à appliquer au volume élémentaire considéré pour imposer une déformation ou contrainte globale
moyenne donnée. Dans le cas linéaire, on peut prouver l’existence et l’unicité pour les différents
cas de conditions aux limites existants.
Principalement pour les composites stratifiés ou sandwichs, il y a deux niveaux d’homogénéisa-
tion :
— du niveau micromécanique au niveau mésoscopique : Les hétérogénéités de base sont les
fibres et la matrice. On effectue ici une étape d’homogénéisation locale.
— du niveau mésoscopique au niveau macroscopique : Les hétérogénéités de base sont les
différentes couches du stratifié. Ces couches sont considérées comme « homogènes » (étape
précédente). Cette fois, il s’agit d’une homogénéisation dans l’épaisseur du stratifié.
V D Vm C Vf D 1 (15.14)
Elongitudinal D El D Ef Vf C Em Vm (15.15)
C’est la loi des mélanges, qui est bien vérifiée dans la direction des fibres. Il s’agit de la
borne supérieure de Voigt (1887).
— deuxième essai — Il s’effectue dans la direction perpendiculaire aux fibres (compression
transversale)
1 1 Vf Vm
D D C (15.16)
Etransverse Et Ef Em
C’est la loi des mélanges en souplesse. Cette relation n’est pas très bien vérifiée transversale-
ment mais donne une indication sur la borne inférieure, dite de Reuss (1929).
— Module de cisaillement et coefficient de Poisson d’un UD par la loi des mélanges :
1 Vf Vm
lt D f Vf C m Vm et D C (15.17)
Glt Gf Gm
Sur le bord des trous, on peut imposer, soit des conditions de Dirichlet :
uD0 (15.19)
n A.x="/ru D 0 (15.20)
u.2/ D u0 .x; x="/ C "u1 .x; x="/ C "2 u2 .x; x="/ (15.21)
1
r u0 C "0 .r u1 C rx u0 / C ".r u2 C rx u1 / C "2 rx u2 D 0
n A./ " (15.23)
d’où :
8
< n A./r u0
ˆ D0
n A./.r u1 C rx u0 / D0 (15.24)
ˆ
n A./.r u2 C rx u1 / D0
:
On peut prouver que si v0 est solution du problème, alors v0 C cte aussi. Des conditions d’existence
O :
de la solution il vient, tous calculs faits, la matrice des coefficients homogénéisés A
n
@Nj
Z X
O D ŒaO i;j
A avec aO i;j D ai k C aij d (15.25)
QnG0 @k
kD1
où G0 est un trou de la cellule élémentaire Q. Le même type de calcul pourrait être mené avec les
conditions de Dirichlet sur le bord des trous. Cette homogénéisation des milieux poreux est valable
pour l’acoustique comme pour la mécanique.
Résumé — Dans ce chapitre, nous nous intéresserons au cas non stationnaire. Toutefois,
nous n’aborderons pas les éléments finis espace-temps, car il s’agit d’une formulation
gourmande en ressources, d’où sa très faible utilisation (bien que la méthode soit en
elle-même intéressante).
Dans ce chapitre, nous aurons besoin de « dériver numériquement », i.e. de construire
des schémas numériques approchant des dérivées. Le chapitre C en annexe permettra
à certains de se rafraîchir la mémoire en regardant comment on résout les équations
différentielles et aux dérivées partielles... puis numériquement. Nous utiliserons en effet
la méthode de Newmark décrite au paragraphe C.2.5.
R / C Cu.t
Mu.t P / C Ku.t / D F.t / (16.1)
R / C Cq.t
Mq.t P / C Kq.t / D F.t / (16.2)
2t
q t C t D q t C t qP t C .1 2ˇ/qR t C 2ˇ qR t C t
2 (16.3)
qP tC t D q t C t .1
/qR t C
qR tC t
P
que nous avons écrit sous la forme discrétisée, mais qui serait valable pour l’approximation continue
(comme présenté au chapitre au paragraphe C.2.5).
Les différents schémas de Newmark correspondent à des valeurs particulières de ˇ et
.
Kq t C t D R (16.5)
avec :
1 1
KD 2
MC C (16.6)
t 2 t
et :
1 1
R D Ft Kq t C M 2q t qt t C Cq t t (16.7)
2t 2 t
qui ne dépendent bien que de données disponibles. Le calcul se déroule donc comme suit :
— Solution initiale à t D 0 :
— connaissant q0 et qP 0 , on calcule qR 0 en résolvant l’équation « classique » de la dyna-
mique ;
2
— on calcule ensuite q t D q0 t qP 0 C 2t qR 0 .
— On dispose alors de tous les éléments pour, à chaque pas de temps t C t :
— résoudre l’équation modifiée, qui nous fournit q t Ct ;
— puis obtenir qP t et qR t par le schéma de Newmark adopté, ici celui des différences finies
centrées.
On peut faire les remarques suivantes sur ce type de calcul :
— Pour un pas de temps donné, q t ne dépend que des données du temps passé, on a donc une
résolution vectorielle rapide.
— Si les matrices M et C sont diagonales, alors cette méthode est très efficace même pour les
problèmes de grande taille.
— Ce schéma est inconditionnellement stable si t 6 Tmin =, avec Tmin la plus petite période
du système correspondant à l’équation de la dynamique (classique).
— la précision est de l’ordre de 2t .
— L’amortissement numérique est nul.
— Il est possible d’introduire un amortissement numérique pour contrôler les hautes fréquences.
Dans ce cas, il faut considérer le schéma avec ˇ D 0 et
> 1=2. Il n’y a pas automatiquement
stabilité du schéma, celle-ci est à calculer pour chaque schéma.
III ÉLÉMENTS FINIS 16.2 Schéma explicite : différences finies centrées 214
et l’équation de la dynamique discrétisée s’écrit encore sous la forme d’une équation modifiée :
Kq t C t D R (16.9)
qui dépendent également des données au même pas de temps. Le calcul se déroule donc comme
suit :
— Solution initiale à t D 0 :
— connaissant q0 et qP 0 , on calcule qR 0 en résolvant l’équation « classique » de la dyna-
mique ;
— on construit K et, si M, C, K et t sont constants (ce qui est généralement le cas), la
triangulariser.
— À chaque pas de temps t C t :
— calculer R (que l’on appelle le résidu, d’où le choix de la notation) ;
— calculer K et triangulariser si nécessaire ;
— résoudre l’équation modifiée, qui nous fournit q t Ct ;
— puis obtenir qP t et qR t par le schéma de Newmark adopté, ici celui de Newmark implicite.
On peut faire les remarques suivantes sur ce type de calcul :
— Pour un pas de temps donné, q t dépend également des données du même pas de temps, on a
donc une résolution matricielle coûteuse.
— Ce schéma est inconditionnellement stable, et donc comme on peut utiliser de plus grands
pas de temps, on réduit le coût mentionné à la ligne précédente.
— la précision est de l’ordre de 2t , et donc comme on ne peut pas utiliser de trop grands pas de
temps sans réduire la précision...
— L’amortissement numérique est nul.
— Il est possible d’introduire un amortissement numérique pour contrôler les hautes fréquences.
Dans ce cas, il faut considérer le schéma avec
> 1=2 et ˇ D .
C 12 /2 . On obtient encore
un schéma stable.
Endommagement
Explicite
Flambement
Plasticité
Implicite
Elasticité
V
Statique Vibrations Propagation d'ondes Supersonique
1. Lorsque l’on parle de codes éléments finis, on cite les grands noms du domaine... et ils correspondent à des codes
extraordinairement puissants, mais souvent chers.
Or un certain nombre de codes éléments finis professionnels sont disponibles gratuitement (au téléchargement et à
l’utilisation). Les deux plus connus sont C AST 3M (anciennement castem) et C ODE A STER développés et maintenus par
le CEA et par EDF respectivement.
2. Notons que pour ceux qui préfèrent, il existe aussi R HEOLEF, un environnement C++ pour le calcul par éléments
finis, qui reste proche, dans sa philosophie, de ce que propose F REE F EM ++. Nous ne le connaissons pas assez pour
fournir des exemples dans ce document.
III ÉLÉMENTS FINIS 16.5 Exemple : un calcul de propagation avec F REE F EM ++ 216
17 / / p l o t ( Th , w a i t = 1 ) ;
Nous définissons quelques valeurs afin de pouvoir définir le chargement (le chargement défini
par la fonction w1 ne sera pas utilisée dans le calcul effectué). Notons comme il est aisé de définir
des fonctions :
Ensuite, nous définissons l’espace des éléments finis utilisés, et les variables appartenant à ces
espaces. L’analogie avec la « formulation mathématique » saute au yeux... et c’est pourquoi cet
outil nous semble particulièrement intéressant sur le plan pédagogique.
50 / / Espaces elements f i n i s
51 fespace Vh ( Th , P2 ) ;
52 fespace Wh ( Th , P1dc ) ;
53
54 Vh w , wa , wd , wda , wdd , wdda , fi ;
55 Wh sxz , syz ;
56 \ end { Verbatim }
Il s’agit d’un exemple non stationnaire, nous allons donc faire une boucle sur le temps. Dans
F REE F EM ++, nous entrons directement la formulation variationnelle du problème sous forme
« mathématique »... le lien entre pratique des éléments finis et théorie est beaucoup plus clair.
La boucle de temps est basée sur un schéma de Newmark implicite. Le temps « courant » est
le temps t C t et le temps précédent est le temps t pour rester cohérent avec les notations du
chapitre 16.
Il s’agit d’un problème de déplacement hors plan, dont le champ inconnu est traditionnellement
noté w. Ainsi, dans le listing, w représente q t C t , wa q t , wd qP t C t , wda qP t , wdd qR t C t et
wdda qR t . fi ' est la fonction test. La formulation variationnelle, en n’indiçant pas le temps
actuel et en utilisant l’indice a pour le pas de temps précédent (a comme avant), est :
2t
@wP @' @wP @'
ZZ
P '/ D
antiplane.w; P C
w' C
Th 4 Th @x @x @y @y
t t
Z Z
g ' wP a C wR a '
2 Th Th 2
t @wa @' @wa @'
Z
(16.12)
C C
2 Th @x @x @y @y
2
@wP a @' @wP a @'
Z
C t C
4 Th @x @x @y @y
CCondition .wP D 0 sur les lignes 1; 2; 3; 4/