Panorama du cours d’aujourd’hui
Motivation et cadre math´ematique (poly §4.1)
M´ethode de Galerkine (poly §4.2)
´ ements finis en dimension 1 (poly §4.3)
El´
Bref historique de la m´ethode
Introduction aux ´
el´
ements finis
issue des travaux de Courant (∼1930)
d´evelopp´ee par les ing´enieurs en a´eronautique (∼1950)
analys´ee dans les ann´ees ∼1970
utilis´ee dans de multiples applications (cf. cours de M´ecanique)
Alexandre Ern
[email protected]
http://cermics.enpc.fr/cours/CS
Richard Courant (1888–1972)
Calcul Scientifique, 26 mars 2015
A. Ern (ENPC)
26 mars 2015
1 / 24
Probl`eme mod`ele
A. Ern (ENPC)
26 mars 2015
Cadre math´ematique
Membrane ´elastique `a l’´equilibre sous un chargement
f
Formulation variationnelle (cf. cours d’Analyse)
Ω

1
 Chercher u ∈ H0 (Ω) tel que
Z
Z

∇u·∇v =
fv ∀v ∈ H01 (Ω)
0000000000000000000000000
1111111111111111111111111
1111111111111111111111111
0000000000000000000000000
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
∂Ω
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
0000000000000000000000000
1111111111111111111111111
u
Ω
membrane tendue occupant le domaine Ω ⊂ R2
fix´ee sur son pourtour
on applique un chargement vertical f
H01 (Ω) = {v ∈ L2 (Ω); ∇v ∈ L2 (Ω)2 ; v |∂Ω = 0}
´equip´e de la norme
La fonction u est solution du probl`eme aux limites
dans Ω
(bilan des forces)
u=0
sur ∂Ω
(condition limite)
Ω
L’espace fonctionnel pour la solution et les fonctions tests est
u:Ω→R
−∆u = f
(1)
principe des travaux virtuels
v est la fonction test
L’inconnue est le d´eplacement vertical de la membrane `a l’´equilibre
A. Ern (ENPC)
2 / 24
Z
kv kH 1 =
3 / 24
A. Ern (ENPC)
Z
2
|∇v |
v +
Ω
26 mars 2015
2
1/2
= kv k2L2 + k∇v k2L2
1/2
Ω
26 mars 2015
4 / 24
Th´eor`eme de Lax–Milgram
Application `a l’´equilibre d’une membrane
On suppose f ∈ L2 (Ω)
Le probl`eme abstrait
(
Le probl`eme (1) est bien pos´e. En effet,
Chercher u ∈ V tel que
(2)
a(u, v ) = b(v ) ∀v ∈ V
lm2 la forme lin´eaire b est continue
Z
|b(v )| = fv ≤ kf kL2 kv kL2 ≤ kf kL2 kv kH 1
est bien pos´e sous les conditions (suffisantes) suivantes :
lm1 l’espace V ´equip´e de la norme k·kV est un espace de Hilbert
∀v ∈ V ,
β = kf kL2
Ω
lm3 la forme bilin´eaire a est continue
Z
|a(v , w )| = ∇v ·∇w ≤ k∇v kL2 k∇w kL2 ≤ kv kH 1 kw kH 1
lm2 la forme lin´eaire b est continue sur V
∃β,
lm1 H01 (Ω) ´equip´e de la norme k·kH 1 est un espace de Hilbert
|b(v )| ≤ βkv kV
ω=1
Ω
lm3 la forme bilin´eaire a est continue sur V × V
∃ω,
∀v , w ∈ V ,
lm4 grˆ
ace `
a l’in´egalit´e de Poincar´e
|a(v , w )| ≤ ωkv kV kw kV
∃CΩ ,
lm4 la forme bilin´eaire a est coercive sur V
∃α > 0,
∀v ∈ V ,
∀v ∈ H01 (Ω),
kv kL2 ≤ CΩ k∇v kL2
la forme bilin´eaire a est coercive
a(v , v ) ≥ αkv k2V
A. Ern (ENPC)
a(v , v ) = k∇v k2L2 ≥
26 mars 2015
5 / 24
M´ethode de Galerkine
2
1
2 kv kH 1
1+CΩ
α=
1
2
1+CΩ
A. Ern (ENPC)
26 mars 2015
6 / 24
Interpr´etation ´energ´etique
Rappel du probl`eme mod`ele
(
Chercher u ∈ V tel que
On suppose que a est sym´etrique
(2)
a(u, v ) = b(v ) ∀v ∈ V
Rappel partie optimisation : la solution exacte u est l’unique minimiseur
(point critique) sur V de la fonctionnelle d’´energie
Vh de dimension finie et Vh ⊂ V (approximation conforme)
E(v ) :=
on cherche une solution approch´ee uh ∈ Vh
on restreint les fonctions tests `
a vh ∈ Vh
on conserve les formes (bi)lin´eaires a et b (mˆeme physique!)
On obtient le probl`eme discret
(
Chercher uh ∈ Vh tel que
(3)
a(uh , vh ) = b(vh ) ∀vh ∈ Vh
1
a(v , v ) − b(v )
2
Exemple pour la membrane : la solution exacte u minimise sur H01 (Ω)
l’´energie m´ecanique
Z
Z
1
2
E(v ) :=
|∇v | −
fv
2 Ω
Ω
La m´ethode de Galerkine pr´eserve ce principe de moindre ´energie : la solution
approch´ee uh minimise la mˆeme ´energie E, mais uniquement sur Vh
Il s’agit de la m´
ethode de Galerkine
Comme Vh ⊂ V , on a
Boris Grigorievitch Galerkine (1871–1945)
E(uh ) = min E(vh ) ≥ min E(v ) = E(u)
vh ∈Vh
A. Ern (ENPC)
26 mars 2015
7 / 24
A. Ern (ENPC)
v ∈V
26 mars 2015
8 / 24
Le syst`eme lin´eaire (1)
Le syst`eme lin´eaire (2)
On introduit la matrice de rigidit´e A ∈ RN×N et le vecteur chargement
B ∈ RN
Aij = a(ϕj , ϕi )
Bi = b(ϕi )
Rappel du probl`eme discret
(
Chercher uh ∈ Vh tel que
(3)
a(uh , vh ) = b(vh ) ∀vh ∈ Vh
R´esoudre (3) ⇐⇒ R´esoudre AU = B
Un point de vue ´equivalent est de r´esoudre un syst`eme lin´eaire
uh solution de (3) ⇐⇒ a(uh , vh ) = b(vh )
⇐⇒ a(uh , ϕi ) = b(ϕi )
P
N
⇐⇒ a
= Bi
j=1 Uj ϕj , ϕi
Soit (ϕ1 . . . ϕN ) une base de Vh
Au lieu de chercher uh ∈ Vh , on cherche ses composantes dans cette base
uh (x) =
N
X
⇐⇒
U = (Uj )1≤j≤N ∈ RN
Uj ϕj (x)
N
X
a(ϕj , ϕi )Uj = Bi
∀vh ∈ Vh
∀i ∈ {1 . . . N}
∀i ∈ {1 . . . N}
∀i ∈ {1 . . . N}
j=1
⇐⇒
j=1
PN
j=1
Aij Uj = Bi
∀i ∈ {1 . . . N}
⇐⇒ AU = B
A. Ern (ENPC)
26 mars 2015
9 / 24
La matrice A est d´efinie positive
Pour tout X ∈ RN , il vient avec ξ(x) =
hAX , X i =
PN
=
PN
i,j=1
A. Ern (ENPC)
10 / 24
Lemme de C´ea
PN
j=1
Xj ϕj (x),
Estimation de l’erreur ku − uh kV
Aij Xi Xj
Lemme de C´ea :
ku − uh kV ≤
i,j=1 Xi Xj a(ϕj , ϕi )
P
PN
N
=a
j=1 Xj ϕj ,
i=1 Xi ϕi
ω
inf ku − vh kV
α vh ∈Vh
Estimation optimale : constante d´ependant de la physique × valeur minimale
possible pour l’erreur
= a(ξ, ξ) ≥ αkξk2V ≥ 0
par coercivit´
e. D’o`
u
Corollaire imm´ediat : si (par miracle) u ∈ Vh , alors uh = u!
hAX , X i = 0 =⇒ ξ = 0 =⇒ X = 0
Attention ! Le lemme de C´ea ne porte que sur la norme k·kV (ou toute autre
norme ´equivalente)
Corollaire : Le probl`eme discret est bien pos´e
A. Ern (ENPC)
26 mars 2015
26 mars 2015
11 / 24
A. Ern (ENPC)
26 mars 2015
12 / 24
´ ements finis 1D
El´
Preuve du lemme de C´ea
Corde ´elastique `a l’´equilibre sous un chargement
f
Orthogonalit´e de Galerkine : Pour tout yh ∈ Vh , on a
a(u − uh , yh ) = a(u, yh ) − a(uh , yh ) = b(yh ) − b(yh ) = 0
Ω
u
D’o`
u, par coercivit´
e et continuit´
e de a,
corde tendue occupant le domaine Ω = ]0, 1[
fix´ee `
a ses 2 extr´emit´es
on applique un chargement vertical f
αku − uh k2V ≤ a(u − uh , u − uh )
= a(u − uh , u − vh )
∀vh ∈ Vh
≤ ωku − uh kV ku − vh kV
∀vh ∈ Vh
si bien que
ku − uh kV ≤
L’inconnue est le d´eplacement vertical de la corde u : Ω → R
ω
ku − vh kV
α
Formulation variationnelle

1
 Chercher u ∈ H0 (Ω) tel que
Z
Z

u0v 0 =
fv ∀v ∈ H01 (Ω)
et on passe `a l’infimum en vh ∈ Vh
Ω
A. Ern (ENPC)
26 mars 2015
13 / 24
Construction de l’espace Vh
A. Ern (ENPC)
(4)
Ω
26 mars 2015
14 / 24
Construction de l’espace Vh
Il faut imposer la nullit´e au bord
´
Etape
1 : mailler le domaine Ω
un maillage de Ω est d´efini par la donn´ee de N points distincts dans Ω
Il faut imposer la continuit´e sur Ω
Ki
(1)
0 = x0
x1
xi
xi+1
xN
Vh
xN+1 = 1
= {vh ∈ C 0 (Ω); ∀i ∈ {0 . . . N}, vh |Ki ∈ P1 ; vh (0) = vh (1) = 0}
(1)
On a Vh
(N + 1) mailles Ki := [xi , xi+1 ], ∀i ∈ {0 . . . N}
maillage uniforme : h = 1/(N + 1), xi = ih, ∀i ∈ {0 . . . (N + 1)}
⊂ H01 (Ω)
(1)
Pour tout vh ∈ Vh , vh0 est constante par morceaux, obtenue en d´
erivant
maille par maille (formule des sauts, cf. cours d’Analyse)
´
Etape
2 : fixer un comportement (simple) dans chaque maille
Fonctions affines par morceaux
{vh : Ω → R; ∀i ∈ {0 . . . N}, vh |Ki ∈ P1 }
o`
u P1 := {p(x) = ax + b; a, b ∈ R}
A-t-on la conformit´e dans H01 (Ω)?
A. Ern (ENPC)
26 mars 2015
15 / 24
A. Ern (ENPC)
26 mars 2015
16 / 24
Fonctions chapeau
Base de Vh
(1)
(1)
dim(Vh ) = N et {ϕ1 . . . ϕN } est une base de Vh
Fonctions chapeau {ϕ1 . . . ϕN }
La famille est libre : ∀(α1 . . . αN ) ∈ RN ,
N
X
1
ϕ1
ϕi
ϕN
αi ϕi (x) ≡ 0 =⇒
i=1
N
X
αi ϕi (xj ) = 0
=⇒ αj = 0
0
x1
xi
xN
∀j ∈ {1 . . . N}
i=1
∀j ∈ {1 . . . N}
1
(1)
On a ϕi ∈
(1)
Vh
La famille est g´
en´
eratrice : ∀vh ∈ Vh ,
et ϕi (xj ) = δij , ∀i, j ∈ {1 . . . N}

−1

h (x − xi−1 )
ϕi (x) = h−1 (xi+1 − x)


0

−1

h
0
ϕi (x) = −h−1


0
x ∈ Ki−1
x ∈ Ki
sinon
vh (x) =
x ∈ Ki−1
x ∈ Ki
sinon
N
X
vh (xi )ϕi (x)
i=1
ces deux fonctions sont affines par morceaux
elles co¨ıncident en deux points distincts sur chaque maille
A. Ern (ENPC)
26 mars 2015
17 / 24
Assemblage de la matrice de rigidit´e (1)
Aij = a(ϕj , ϕi ) =
26 mars 2015
On calcule le coefficient diagonal
Z
Z
Aii =
ϕ0i (x)ϕ0i (x)dx =
ϕ0i (x)ϕ0j (x)dx
Z
Observation essentielle : Aij 6= 0 seulement si l’intersection des supports de
ϕi et ϕj est de mesure non-nulle
Z
... +
Ki−1
Ki
2
2
1
1
2
... = h ×
+h× −
=
h
h
h
De mˆeme, Ai,i−1 = Ai−1,i = − h1
...
..
.
..
.
•
•
Au final

0
.. 
.


0


•
•

2

−1
1

A=  0
h
 .
 ..
0
Que valent les coefficients diagonaux et extra-diagonaux?
A. Ern (ENPC)
=
(ϕ0i (x))2 dx
Ki−1 ∪Ki
Ω
Ω
La matrice de rigidit´e est tridiagonale

• •
0

• •
•


.
.
A = 0 . . . .

. .
.. •
 ..
0 ... 0
18 / 24
Assemblage de la matrice de rigidit´e (2)
Rappel du terme g´en´erique
Z
A. Ern (ENPC)
−1
..
..
0
2 −1
..
.
.
.
...
−1
0

0
.. 
.


0


2 −1
−1
2
...
..
.
..
.
En 1D, A a les dimensions de l’inverse d’une longueur
26 mars 2015
19 / 24
A. Ern (ENPC)
26 mars 2015
20 / 24
Interpolation
Preuve
Op´erateur d’interpolation (H01 (Ω) ⊂ C 0 (Ω) en 1D)
(1)
(1)
Ih : H01 (Ω) 3 v 7−→ Ih v :=
N
X
Soit v ∈ C ∞ (Ω) ∩ H01 (Ω), soit une maille Ki = [xi , xi+1 ], soit x ∈ Ki
Z
v (xi+1 ) − v (xi )
1 xi+1 0
(1)
v 0 (x) − (Ih v )0 (x) = v 0 (x) −
=
(v (x) − v 0 (t)) dt
h
h xi
Rx
Par in´egalit´e de Cauchy–Schwarz, comme v 0 (x) − v 0 (t) = t v 00 (s)ds,
(1)
v (xi )ϕi ∈ Vh
i=1
(1)
Ih v prend la mˆeme valeur que v en tous les sommets du maillage
v
|v 0 (x) − v 0 (t)| ≤ |x − t|1/2 kv 00 kL2 (Ki ) ≤ h1/2 kv 00 kL2 (Ki )
(1)
Ih v
(1)
Par suite, |v 0 (x) − (Ih v )0 (x)| ≤ h1/2 kv 00 kL2 (Ki ) ; d’o`
u, par CS,
(1)
kv 0 − (Ih v )0 kL2 (Ki ) ≤ h1/2 h1/2 kv 00 kL2 (Ki ) = hkv 00 kL2 (Ki )
x
En sommant sur les mailles, on conclut que, pour tout v ∈ C ∞ (Ω) ∩ H01 (Ω),
Th´eor`eme d’interpolation : ∃CI t.q. ∀v ∈ H 2 (Ω) ∩ H01 (Ω),
(1)
kv 0 − (Ih v )0 kL2 (Ω) ≤ hkv 00 kL2 (Ω)
(1)
kv 0 − (Ih v )0 kL2 ≤ CI hkv 00 kL2
On ´etend l’in´egalit´e `a H 2 (Ω) ∩ H01 (Ω) par densit´e
v s’´ecarte de son interpol´e affine lorsque sa courbure est grande
A. Ern (ENPC)
(CI = 1)
26 mars 2015
21 / 24
Estimation d’erreur
A. Ern (ENPC)
26 mars 2015
22 / 24
En conclusion et en ouverture
Estimation d’erreur : rappel du lemme de C´ea
ku − uh kV ≤
ω
α v inf
h ∈Vh
ku − vh kV
principe : m´
ethode de Galerkine
0
On travaille en norme kv kV := kv kL2 , ´equivalente `a la norme kv kH 1 sur
H01 (Ω) (in´egalit´e de Poincar´e)
matrice de rigidit´e d´
efinie positive
analyse d’erreur : lemme de C´
ea
La solution exacte v´erifie −u 00 = f , elle est donc dans H 2 (Ω), d’o`
u par le
(1)
th´
eor`
eme d’interpolation dans Vh
´el´ement fini de degr´e 1 : convergence `a l’ordre 1 en norme H 1
(1)
Aujourd’hui en TD : exercices sur l’optimisation
ku 0 − (Ih u)0 kL2 ≤ CI hku 00 kL2 = CI hkf kL2
La semaine prochaine : amphi sur les EF en 2D, exercices en 1D
(1)
On d´eduit du lemme de C´ea avec le choix vh = Ih u
0
ku −
uh0 kL2
≤
0
ω
α ku
−
(1)
(Ih u)0 kL2
≤
ω
α CI kf
Dans quinze jours : TP et exercices en 2D
kL2 h
convergence `
a l’ordre 1 en h en norme H 1
(1)
Ih u n’est pas calcul´e, il sert juste `
a prouver l’estimation d’erreur!
A. Ern (ENPC)
26 mars 2015
23 / 24
A. Ern (ENPC)
26 mars 2015
24 / 24