ANALYSE NUMERIQUE 2

SALEM NAFIRI 1

Ecole Hassania des Travaux Publics 2017

1. Email: [email protected]
Table des matières

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3

1 Dérivation numérique 4
1.1 Dérivation numérique d’ordre 1 . . . . . . . . . . . . . . . . . . . 4
1.2 Evaluation de l’erreur . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3 Dérivation numérique d’ordre 2 . . . . . . . . . . . . . . . . . . . 7
1.4 Dérivation numérique d’ordre supérieur . . . . . . . . . . . . . . 8

2 Intégration numérique 10
2.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.2 Généralités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.3 Formules de quadratures . . . . . . . . . . . . . . . . . . . . . . . 12
2.3.1 Estimation de l’erreur . . . . . . . . . . . . . . . . . . . . 13
2.3.2 Poids d’une formule de quadrature . . . . . . . . . . . . . 14
2.4 Formule du rectangle, du trapèze, de Simpson . . . . . . . . . . . 15
2.4.1 Formule du rectangle . . . . . . . . . . . . . . . . . . . . . 15
2.4.2 Formule du trapèze . . . . . . . . . . . . . . . . . . . . . . 15
2.4.3 Formule de Simpson . . . . . . . . . . . . . . . . . . . . . 16
2.5 Points d’intégration et formules de Gauss-Legendre . . . . . . . . 16

3 Equations différentielles ordinaires (EDOs) 19

3.1 Problème de Cauchy . . . . . . . . . . . . . . . . . . . . . . . . . 19
3.2 Well posedness du problème de Cauchy . . . . . . . . . . . . . . 20
3.3 Approximation du problème de Cauchy . . . . . . . . . . . . . . 22
3.3.1 Schémas d’Euler . . . . . . . . . . . . . . . . . . . . . . . 22
3.3.2 Schéma d’Euler progressif vs Schéma d’Euler rétrograde . 23
3.3.3 Schémas de Runge Kutta . . . . . . . . . . . . . . . . . . 24
3.4 Approximation des EDOs d’ordre supérieur . . . . . . . . . . . . 26

4 Equations aux dérivées partielles 27

4.1 Approximation des équations aux dérivées partielles . . . . . . . 27
4.1.1 Exemples des EDPs classiques . . . . . . . . . . . . . . . 27
4.1.2 Méthode des différences finies . . . . . . . . . . . . . . . . 28
4.1.3 Méthode de Galerkin . . . . . . . . . . . . . . . . . . . . . 29
Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33

2
Introduction

Dans l’élément Analyse numérique 1, nous avons étudiés les sujets suivants :
— Interpolation polynomiale.
— Racine d’équations algébriques.
— Systèmes d’équations linéaires et non linéaires.
Dans l’élement Analyse numérique 2, nous allons suivre notre survol sur
les méthodes numériques élémentaires. En particulier, nous étudions les sujets
suivants :
— Dérivation numérique.
— Intégration numérique.
— Approximations des EDOs.
— Approximations des EDPs.
Nous devrons souligner que ce cours n’est qu’une introduction aux méthodes
numériques élémentaires, qui a pour but de donner à l’élève ingénieur quelques
outils pour étudier certains problèmes de caractère scientifique ou émanents des
sciences de l’ingénieur.
Les notions et concepts du cours d’Analyse numérique 1 sont fondamentaux
pour mieux comprendre le cours d’Analyse numérique 2.
Fait à Casablanca, le 03/02/2017.

3
Chapitre 1

Dérivation numérique

1.1 Dérivation numérique d’ordre 1

La notion de dérivée a vu le jour dans les écrits de Leibniz et de Newton

au 17ième siècle qui le définit comme  le quotient de deux accroissements
évanescents .

Soit f : R −→ R une fonction de classe C 1 et soit x0 ∈ R fixé.

Calculer la dérivée d’une fonction en un point revient à calculer la limite

suivante

f (x0 + h) − f (x0 )
lim = lim tx0 (h).
h→0 h h→0

tx0 (h) est le taux d’accroissement de f au point x0 . Ce calcul de limite revient

graphiquement à rechercher la pente de la tangente à la courbe en ce point.
Ainsi, le nombre dérivé d’une fonction en un point, s’il existe, est égal à la pente
de la tangente à la courbe représentative de la fonction en ce point, voir Figure
1.

4
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

f (x0 +h)−f (x0 )

Figure 1.1: La pente de la tangente en P est lim h .
h→0

Si x0 ∈ R est un réel donné, nous pouvons écrire

f (x0 + h) − f (x0 )
f 0 (x0 ) = lim
h→0 h
f (x0 ) − f (x0 − h)
= lim
h→0 h
f (x0 + h/2) − f (x0 − h/2)
= lim .
h→0 h

De ce qui précède, on déduit que la dérivée d’une fonction en un point x0 , peut

être approchée par différentes taux d’accroissement, on les appelle formules de
différences finies. Dans la suite on va considérer trois types d’approximation de
la dérivée première :
— la formule des différences finies progressive :

f (x0 + h) − f (x0 )
f 0 (x0 ) ≈
h

— la formule des différences finies rétrograde :

f (x0 ) − f (x0 − h)
f 0 (x0 ) ≈
h

— la formule des différences finies centrées :

f (x0 + h2 ) − f (x0 − h2 )
f 0 (x0 ) ≈ .
h

5
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

1.2 Evaluation de l’erreur

Soit h > 0 fixé (petit). Dans ce paragraphe nous allons voir que :
1. La formule des différences finies progressive est d’ordre 1 en h, c.à.d
 

 0 f (x0 + h) − f (x0 ) 
f (x0 ) −  = O(h).
 h 

2. La formule des différences finies rétrograde est d’ordre 1 en h, c.à.d

 

 0 f (x0 ) − f (x0 − h) 
f (x0 ) −  = O(h).
 h 

3. La formule des différences finies centrées est d’ordre 2 en h, c.à.d

 

 0 f (x0 + h2 ) − f (x0 − h2 ) 
f (x0 ) −  = O(h2 ).
 h 

Théorème 1.2.1. Soient f ∈ C 2 et x0 ∈ R, il existe C = C(f, x0 ) > 0, tel que

pour tout 0 < h 6 1 on ait
 

 0 f (x0 + h) − f (x0 ) 
f (x0 ) −  6 Ch.
 h 

Preuve. Soient f : R −→ R une fonction de classe C 2 , x0 ∈ R fixé et soit h > 0.

Le développement de Taylor de f sur l’intervalle [x0 , x0 + h] donne :
h2 00
f (x0 + h) = f (x0 ) + hf 0 (x0 ) + f (ξ), x0 < ξ < x0 + h
2
 

 0 f (x0 + h) − f (x0 )  h  00 
donc f (x0 ) −  = 2 f (ξ)
 h 
 
On ne peut pas choisir C = 12 f 00 (ξ), car ξ dépend de h.
Comme x0 < ξ < x0 + h < x0 + 1, alors
 

 0 f (x0 + h) − f (x0 )  h  
f (x0 ) − 6 max f 00 (ξ)
 
 h  2 x0 6ξ6x0 +1
= Ch
1
 00 
avec C = max
2 x 6ξ6x
f (ξ).
0 0 +1

Théorème 1.2.2. Soient f ∈ C 2 et x0 ∈ R, il existe C = C(f, x0 ) > 0, tel que

pour tout 0 < h 6 1 on ait
 

 0 f (x0 ) − f (x0 − h) 
f (x0 ) −  6 Ch.
 h 

6
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Preuve. On applique le développement

  de Taylor sur l’intervalle [x0 − h, x0 ].
Dans ce cas C = 12 max f 00 (ξ)
x0 −16ξ6x0

Remarque 1.2.1. • La constante C dans les deux théorèmes ci-dessus peut

dépendre de f et de x0 , mais pas de h.
• Interprétation : Si on choisit f et x0
— l’erreur est divisée par 2 chaque fois que h est divisé par 2,
— l’erreur est divisée par 10 chaque fois que h est divisé par 10.
Théorème 1.2.3. Soient f ∈ C 3 et x0 ∈ R, il existe C = C(f, x0 ) > 0, pour
tout 0 < h 6 1 on a
 

 0 f (x0 + h2 ) − f (x0 − h2 ) 
f (x0 ) −  6 Ch2 .
 h 

Preuve. On applique le développement de Taylor à la fonction f sur l’intervalle


[x0 − h2 , x0 ], puis sur [x0 , x0 + h2 ]. Dans ce cas C = 24
1
1
max 1 f 00 (ξ)
x0 − 2 6ξ6x0 + 2

Remarque 1.2.2. • La constante C dans le théorème précédent peut dépendre

de f et de x0 , mais pas de h.
• Interprétation : Si on choisit f et x0
— l’erreur est divisée par 22 = 4 chaque fois que h est divisé par 2,
— l’erreur est divisée par 102 = 100 chaque fois que h est divisé par 10.
• Le fait que f ∈ C 3 entraine une amélioration dans l’ordre d’estimation
de l’erreur. Il est donc naturel de choisir la formule des différences finies
centrées pour approher la dérivée d’une fonction.

1.3 Dérivation numérique d’ordre 2

Soit f : R −→ R une fonction de classe C 2 et soit x0 ∈ R fixé.
On remarque que f 00 (x0 ) = (f 0 )0 (x0 ), c. à. d, que la dérivée numérique
d’ordre 2, peut être obtenue en utilisant les schémas numériques déjà étudiés
pour la dérivée première. On choisit, à cet effet, la formule des différences finies
centrées pour trouver une approximation de la dérivée d’ordre 2

f 00 (x0 ) = (f 0 )0 (x0 )
f 0 (x0 + h2 ) − f 0 (x0 − h2 )
= lim
h→0 h
or comme
h
h f (x0 + + h2 ) − f (x0 + h
− h2 ) f (x0 + h) − f (x0 )
f 0 (x0 + )≈ 2 2
=
2 h h
h
h f (x0 − + h2 ) − f (x0 − h
− h2 ) f (x0 ) − f (x0 − h)
f 0 (x0 − ) ≈ 2 2
=
2 h h

7
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

donc
f (x0 + h) − 2f (x0 ) + f (x0 − h)
f 00 (x0 ) ≈ .
h2
Le théorème suivant montre que l’erreur commise sur l’approximation de la
dérivée d’ordre 2 est d’ordre 2.
Théorème 1.3.1. ∀f ∈ C 4 , ∀x0 ∈ R, ∃C > 0, ∀0 < h 6 1
 

 00 f (x0 + h) − 2f (x0 ) + f (x0 − h) 
f (x0 ) −  6 Ch2 .
 h2 

Remarque 1.3.1. • La constante C dépend de f et de x0 , mais pas de h.

• Interprétation : Si on choisit f et x0
— l’erreur est divisée par 22 = 4 chaque fois que h est divisé par 2,
— l’erreur est divisée par 102 = 100 chaque fois que h est divisé par 10.
Preuve. La preuve est laisée en exercice.

1.4 Dérivation numérique d’ordre supérieur

Cette partie généralise les résultats déjà énoncés dans les paragraphes précédents.
Pour ce faire, on considère les notations suivantes

Dp f (x0 ) = f (x0 + h) − f (x0 )

Dr f (x0 ) = f (x0 ) − f (x0 − h)
 h  h
Dc f (x0 ) = f x0 + − f x0 − .
2 2
Ainsi, on peut calculer les dérivées numériques d’ordres supérieurs (m > 1), en
utilisant l’observation suivante

Dpm f (x0 ) = Dp (Dpm−1 f )(x0 )

Drm f (x0 ) = Dr (Drm−1 f )(x0 )
Dcm f (x0 ) = Dc (Dcm−1 f )(x0 ).

Le théorème suivant généralise Théorème 2.1 et Théorème 2.2.

Théorème 1.4.1. Soient m ∈ N, f ∈ C m+1 , x0 ∈ R. Alors il existe C > 0 tel
que pour tout 0 < h 6 1
 

 (m) Dpm f (x0 ) 
f (x0 ) −  6 Ch
 hm 
 

 (m) Drm f (x0 ) 
f (x0 ) −  6 Ch.
 hm 

8
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Preuve. Exercice.
Contrairement au théorème précédent, le résultat suivant donne une estima-
tion de l’erreur d’ordre 2 pour toute dérivée d’ordre supérieure.
Théorème 1.4.2. ∀m ∈ N, ∀f ∈ C m+2 , ∀x0 ∈ R, ∃C > 0, ∀0 < h 6 1
 

 (m) Dcm f (x0 ) 
f (x0 ) −  6 Ch2
 hm 

Preuve. Exercice.
Remarque 1.4.1. Les problèmes de diffusion d’espèces, de déformations élastiques,
de propagation d’ondes, d’écoulement de fluides, etc... font intervenir des dérivées
deuxième ou quatrième.

Exercice : Utiliser la formule des différences finies centrées et montrer que

f (x0 + 2h) − 4f (x0 + h) + 6f (x0 ) − 4f (x0 − h) + f (x0 − 2h)

f (4) (x0 ) ≈ .
h4
f (x0 +h)−2f (x0 )+f (x0 −h)
(Indication : utiliser f 00 (x0 ) ≈ h2 .)

9
Chapitre 2

Intégration numérique

2.1 Motivation
Le calcul d’aires, de surfaces, de volumes et longueurs d’arc était une des
préoccupations majeures des scientifiques depuis Euclide et Archimède jusqu’à
Kepler, Galilei, Fermat et Descartes. Puis, avec la découverte du Théorème fon-
damental du calcul différentiel, ces  intégrales  n’étaient qu’une inversion de
formules pour la dérivée. Les Bernoullis et surtout Euler, dans son gigantesque
traité (Inst. Calculi Integralis 1768, 1769, 1770, Opera XI-XIII), se sont donc mis
avec enthousiasme à la recherche de formules analytiques pour les primitives.
Cependant, de nombreuses intégrales résistent encore et toujours à l’envahisseur
R ex R dx
(par exemple dx, )et l’enthousiasme laisse lentement place à un re-
x log x
nouveau de l’intérêt pour les formules de quadrature, qui avait été commencé
par Newton et ses successeurs anglais Cotes et Simpson. De nombreux calculs en
astronomie (perturbations des orbites planétaires) contraignent Gauss (1814) à
intensifier cette théorie. Les premiers ordinateurs, dans les années 40 et 50, ont
surtout été développés pour le calcul d’intégrales et d’équations différentielles.
En effet, ces problèmes sont les plus faciles à programmer.

2.2 Généralités
Soit f : x ∈ [a, b] −→ f (x) ∈ R une fonction continue donnée. Nous désirons
approcher numériquement la quantité

Z b
I= f (x)dx.
a

Comme montre la figure suivante

10
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Figure 2.1: L’intégrale définie de f sur [a, b].

Pour ce faire, nous commençons par partitionner l’intervalle [a, b] en petits in-
tervalles [xi , xi+1 ], i = 0, 1, 2 · · · , N − 1, c’est à dire nous choisissons des points
xi , i = 0, 1, 2, · · · , N tels que
a = x0 < x1 < x2 < x3 < · · · < xN −1 < xN = b. (2.1)
Soit h = max |xi+1 − xi | le réel positif caractérisant la finesse de la partition.
06i6N −1
Il est clair que, lorsque N augmente, nous pouvons placer les points xi de sorte
à ce que h soit petit. Lorsqu’aucune raison nous incite à choisir des intervalles
de longueurs différentes, nous posons
b−a
h= et xi = a + ih, i = 0, 1, 2 · · · , N.
N
Etant donné la partition (2.1), il est naturel d’écrire :
Z b N
X −1 Z xi+1
f (x)dx = f (x)dx.
a i=0 xi
Rx
Remarque 2.2.1. Dans la suite c’est la quantité xii+1 f (x)dx que nous allons
approcher par ce que l’on appelle formules de quadratures.
On peut toujours introduire une changement de variable et ramener [xi , xi+1 ]
à [−1, 1]. On obtient alors
Z xi+1
xi+1 − xi 1
Z
f (x)dx = gi (t)dt
xi 2 −1
 
t+1
où gi (t) = f xi + (xi+1 − xi ) , t ∈ [−1, 1].
2
Nous sommes maintement en mesure Rde définir la notion de formule de qua-
1
drature pour approcher numériquement −1 gi (t)dt, g étant une fonction conti-
nue donnée sur [−1, 1].

11
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

2.3 Formules de quadratures

Définition 2.3.1. Si g une fonction continue sur [−1, 1], la formule de qua-
drature
M
X
J(g) = ωj g(tj )
j=1

est définie par la donnée de M points −1 6 t1 < t2 < · · · < tM 6 1 appelés

points d’intégration et de M nombre réels ω1 , ω2 , · · · , ωM appelés poids de la
formule de quadrature. Ces M points et ces M poids devront R 1 être cherchés de
façon à ce que J(g) soit une approximation numérique de −1 g(t)dt.
Exemple 2.3.1. Si M = 1, on peut choisir t1 = 0 et donc J(g) = 2g(0).
Si M = 2, on peut choisir t1 = −1 et t2 = 1, donc J(g) = g(−1) + g(1).

Figure 2.2: Cas M=1

Figure 2.3: Cas M=2

Retournons maintenant à notre question de départ : comment approcher numériquement

Rb
a
f (x)dx ?
Les calculs précédents montrent que
b N −1 Z
h X 1
Z
t+1
f (x)dx = f (xi + h )dt
a 2 i=0 −1 2
N −1 M
h XX tj + 1
' ωj f (xi + h )dt = Lh (f ).
2 i=0 j=1 2

12
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Rb
Nous allons donc approcher numériquement a f (x)dx par Lh (f ) qui est dite
formule composite.
En règle générale nous pouvons procéder de la manière suivante pour ap-
Rb
procher la quantité a f (x)dx par la quantité Lh (f ) : on définit une formule
de quadrature par la donnée de M points t1 < t2 < · · · < tM et de M poids
ω1 , ω2 , · · · , ωM ; on partitionne l’intervalle [a, b] en intervalles [xi , xi+1 ], puis on
clacule Lh (f ) par la formule composite.

Remarque 2.3.1. Deux questions standards peuvent être posées :

1. Comment choisir les points d’intégrations t1 < t2 < · · · < tM ainsi que
les poids de la formule de quadrature ω1 , ω2 , · · · , ωM ?
Rb
2. Quelle est l’erreur commise entre a f (x)dx et Lh (f ) ?

2.3.1 Estimation de l’erreur

Avant de montrer comment on construit des formules de quadrature, définissons
une propriété désirable de J(g).

Définition 2.3.2. On dira que la formule de quadrature

M
X
J(g) = ωj g(tj )
j=1

R1
pour calculer numériquement −1
g(t)dt est exacte pour les polynômes de degré
r > 0 si
Z 1
J(p) = p(t)dt
−1

pour tout polynôme p ∈ Pr (L’espace des polynômes de degré inférieur ou égal

à r).

Théorème 2.3.1. Supposons que la formule de quadrature

M
X
J(g) = ωj g(tj )
j=1

R1
pour calculer numériquement −1 g(t)dt soit exacte pour les polynômes de degré
r. Soit f une fonction donnée sur l’intervalle [a, b], soit Lh (f ) la formule com-
posite et soit h un pas uniforme positif. Alors, si f ∈ Cr+1 [a, b], il existe C > 0
indépendante de h telle que
Z
 b 
f (x)dx − Lh (f ) 6 Chr+1 .
 

 a

13
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

2.3.2 Poids d’une formule de quadrature

Dans cette section nous supposons donnés M points d’intégration distincts
de l’intervalle [−1, 1]

−1 6 t1 < t2 < t3 < · · · < tM 6 1

et nous cherchons à déterminer

PM les poids ω1 , ω2 , · · · , ωM de sorte que la formule
de quadrature J(g) = j=1 ωj g(tj ) soit exacte pour les polynômes de degré r
aussi élevé que possible.
Pour réaliser cette objectif, considérons la base de Lagrange ϕ1 , ϕ2 , · · · , ϕM
de PM −1 assosiés aux points t1 , t2 , · · · , tM . Par définition, ϕj est le polynôme
de degré M − 1 défini par
M
Y t − tk
ϕj (t) = ,
tj − tk
k=1
k6=j

pour j = 1, 2, · · · , M . Soit g : t ∈ [−1, 1] −→ g(t) ∈ R une fonction continue

donnée. Son interpolant g̃ de degré M − 1 aux points t1 , · · · , tM est défini par
M
X
g̃(t) = g(tj )ϕj (t).
j=1
R1 R1
Il semble naturelle de remplacer −1
g(t)dt par −1
g̃(t)dt. Puisque
Z 1 M
X Z 1
g̃(t)dt = g(tj ) ϕj (t)dt,
−1 j=1 −1

nous constatons immediatement qu’il suffit de poser

Z 1
ωj = ϕj (t)dt
−1
PM R1
pour que J(g) = j=1 ωj g(tj ) soit une approximation de −1
g(t)dt.
Nous obtenons ainsi le théorème suivant
Théorème 2.3.2. Soit t1 , t2 , · · · , tM M points distincts de l’intervalle [−1, 1]
et soit ϕ1 , ϕ2 , · · · , ϕM la base de Lagrange de PM −1 associée à ces M points.
Alors la formule de quadrature
M
X
J(g) = ωj g(tj )
j=1

est exacte pour les polynôme de degré M − 1 si et seulement si

Z 1
ωj = ϕj (t)dt, j = 1, 2, · · · , M.
−1

14
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

2.4 Formule du rectangle, du trapèze, de Simp-

son
2.4.1 Formule du rectangle
La formule du rectangle correspond au cas M = 1, c’est à dire J(g) = 2g(0),
ce qui donne
N −1
X h
Lh (f ) = h f (xi + ).
i=0
2

Selon Théorème 2.3.2, cette formule de quadrature est exacte pour des po-
lynômes de degré zéro, mais en fait elle est meilleure ; elle est exacte pour des
polynôme p ∈ P1 . Ainsi, l’estimation du Théorème 2.3.1 devient
Z
 b 
f (x)dx − Lh (f ) 6 Ch2 .
 

 a

L’interprétation géométrique de la formule de rectangle est la suivante : on

somme les aires des rectangles dont la base est le segment [xi + xi+1 ] et dont la
hauteur est f (ξi ), où ξi est le point milieu de [xi , xi+1 ].

Figure 2.4: Formule du rectangle sur [-1,1].

2.4.2 Formule du trapèze

La formule du trapèze correspond au cas M = 2, c’est à dire J(g) = g(−1) +
g(1), ce qui donne
N −1  
h X
Lh (f ) = f (xi ) + f (xi+1 ) .
2 i=0

Selon Théorème 2.3.2, Rcette formule de quadrature est exacte pour des po-
1
lynômes de degré 1, i.e., −1 p(t)dt = J(p), ∀p ∈ P1 . L’estimation du Théorème
2.3.1 devient alors
Z
 b 
f (x)dx − Lh (f ) 6 Ch2 .
 

 a

15
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Figure 2.5: Formule du trapèze sur [-1,1]

2.4.3 Formule de Simpson

La formule de Simpson correspond au cas M = 3, c’est à dire
J(g) = ω1 g(−1) + ω2 g(0) + ω3 g(1).
R1
D’après Théorème 2.3.2, on a ωj = −1 ϕj (t)dt, j = 1, 2, 3. Un calcul simple
donne
1 4 1
ω1 = , ω2 = , ω3 = .
3 3 3
On déduit alors que
N −1  
h X
Lh (f ) = f (xi ) + 4f (xi + h/2) + f (xi+1 ) .
6 i=0

Cette formule de quadrature est exacte pour des polynômes de degré 2 et 3.

L’estimation du Théorème 2.3.1 devient alors
Z
 b 
f (x)dx − Lh (f ) 6 Ch4 .
 

 a

Remarque 2.4.1. La formule de Simpson donne une erreur d’ordre h4 . C’est

une formule souvent utilisée dans la pratique car Lh (f ) converge très rapidement
Rb
vers a f (x)dx lorsque h → 0.

2.5 Points d’intégration et formules de Gauss-

Legendre
L’idée des formules de Gauss-Legendre est de placer au mieux les P points
M
d’intégrations t1 , t2 , · · · , tM de sorte que la formule de quadrature J(p) = j=1 ωj p(tj )
R1
soit égale à −1 p(t)dt pour des polynôme de degré r aussi grand que possible.
En effet, rappelons que si f est une fonction donnée et si Lh (f ) est l’approxi-
Rb
mation de a f (x)dx, alors en vertu du Théorème 2.3.1, plus r est grand et plus
Rb
l’erreur entre a f (x)dx et Lh (f ) tend rapidement vers zéro avec h.
Commençons par la définition suivante

16
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Définition 2.5.1. Le polynôme de Legendre de degré M est défini par

1 dM 2
LM (t) = (t − 1)M .
2M M ! dtM
Ainsi donc nous avons, si t ∈ R :

3t2 − 1
L0 (t) = 1, L1 (t) = t, L2 (t) = , ···
2
Les polynômes de Legendre L0 , L1 , L2 , · · · , vérifient de nombreuses propriétés.
Entre autres, nous démontrons certaines propriétés qui nous serons utiles dans
la suite.
Théorème 2.5.1. Les polynômes de Legendre L0 , L1 , L2 , · · · vérifient les pro-
priétés suivantes
i) L0 , L1 , L2 , · · · , LM forment une base de PM .
R1
ii) Si i 6= j alors −1 Li (t)Lj (t)dt = 0 (Propriété d’orthogonalité)
iii) LM a exactement M zéros réels distincts tous compris dans l’intervalle
ouvert ]-1,1[. Ces zéros sont appelés les points de Gauss.
Définition 2.5.2. Nous dirons que la formule de quadrature
M
X
J(g) = ωj g(tj )
j=1

est la formule de Gauss-Legendre à M points si

i) Les points d’intégration t1 , t2 , · · · , tM sont les M zéros du polynôme de Le-
gendre LM , c’est à dire les M points de Gauss.
R1
ii) Les poids sont déinis par la relation ωj = −1 ϕj (t)dt, j = 1, 2, · · · , M , où
ϕ1 , ϕ2 , · · · , ϕM est la base de Lagrange de PM −1 associée aux M points
de Gauss.
Nous sommes maintenant en mesure d’énoncer le résultat suivant
Théorème 2.5.2. La formule de Gauss-Legendre à M points (M entier > 1)
est exacte pour les polynômes de degré r = 2M − 1.
Tenant compte des Théorèmes 2.3.1 et 2.5.2, nous avons
Z
 b 
f (x)dx − Lh (f ) 6 Ch2M .
 

 a

où ici f ∈ C 2M et C une constante idépendante de h.

Exemple 2.5.1. Formule de Gauss-Legendre à un seul point
On a L1 (t) = t et donc le seul zéro de L1 est donné par t1 = 0. Nous
retrouvons dans ce cas-là, la formule du rectangle qui est d’ordre h2 .

17
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Exemple 2.5.2. Formule de Gauss-Legendre à deux points

Nous avons L2 (t) = 12 (3t2 − 1) et donc les deux zéros de L2 sont donnés par

1 1
t1 = − √ et t2 = √ .
3 3
La base de Lagrange ϕ1 , ϕ2 de P2 associée aux points t1 , t2 est définie par
√ √
1 − 3t 1 + 3t
ϕ1 (t) = et ϕ2 (t) =
2 2
et ainsi

ω1 = ω2 = 1.

La formule de Gauss-Legendre à deux points s’écrit donc

√
√

J(g) = g − 1/ 3 + g 1/ 3 .

Ce qui entraine
N −1
( √ ! √ !)
X xi+1 − xi 3−1 3+1
Lh (f ) = f xi + √ (xi+1 − xi ) + f xi + √ (xi+1 − xi ) .
i=0
2 2 3 2 3

Si f est de classe C 4 les théorèmes 2.3.1 et 2.5.2, nous assurent l’existence

d’une constante C indépendante de h telle que
Z
 b 
f (x)dx − Lh (f ) 6 Ch4 .
 

 a

La formule de Gauss-Legendre à deux points converge donc au même ordre que

la formule de Simpson.

18
Chapitre 3

Equations différentielles
ordinaires (EDOs)

Soit f : R × R+ −→ R
.
(x, t) 7−→ f (x, t)
On suppose que f est une fonction assez régulière, par exemple continue.
Notre objectif dans cette section est de résoudre le problème suivant.

3.1 Problème de Cauchy

u : R+ −→ R
Etant donné u0 (valeur initiale), trouver de classe
t 7−→ u(t)
C 1 telle que
u0 (t) = f (u(t), t),

t>0
(PC)
u(0) = u0 .
Le problème (PC) s’appelle problème de Cauchy.
Exemple 3.1.1. On prend f (x, t) = 3x − 3t et u0 = α (un nombre qcq). Le
problème de Cauchy devient donc
 0
u (t) = 3u(t) − 3t, t > 0
u(0) = α.
La formule de la variation de la constante
Z t
u(t) = αe3t + e3(t−s) .(−3s)ds
0
Z t
= αe − 3e3t
3t
se−3s) ds
0
1 1
= (α − )e3t + t + .
3 3

19
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

√
Exemple 3.1.2. Pour f (x, t) = 3 x, u0 = 0. Le problème de Cauchy devient
 0 p
u (t) = 3 u(t), t > 0
u(0) = 0.

Exercice : On peut p vérifier que les fonctions u définies par

u(t) = 0 et u(t) = ± 8t3 /27, t > 0 sont toutes trois des solutions du problème
de Cauchy. Cet exemple nous montre que le problème de Cauchy n’a pas forcément
une solution unique.
Exemple 3.1.3. Pour f (x, t) = x3 , u0 = 1. Le problème de Cauchy devient
 0
u (t) = u3 (t), t > 0
u(0) = 1.

Exercice : On peut vérifier que la solution u est donnée par

1
u(t) = √ , t ∈ [0, 1/2[.
1 − 2t
Cet exemple montre que le problème de Cauchy n’a pas toujours une solution
globale sur R, puique ici la solution u explose lorsque t → 1/2.
Remarque 3.1.1. Les trois exemples ci-dessus montrent que l’étude mathématique
de l’existence et de l’unicité des solutions du problème de Cauchy peut être une
affaire délicate.

3.2 Well posedness du problème de Cauchy

Dans ce paragraphe nous donnons les hypothèses pour lesquelles le problème
de Cauchy (PC) admet une solution unique. Nous nous contentons de donner
(sans démonstration) un résultat d’existence et d’unicité global au problème de
Cauchy. Avant d’énoncer le théorème qui répond à la question de well posedness
du problème de Cauchy (PC), on rappelle la définition donnée par Hadamard
Définition 3.2.1. Le problème de Cauchy (PC) est dit bien posé, au sens de
Hadamard si :
(i) la solution existe,
(ii) la solution est unique,
(iii) la solution dépend de façon continue des données.
Théorème 3.2.1. Soit f : R×R+ −→ R continue, vérifiant l’hypothèse suivante
(H) il existe ` : t ∈ R+ −→ `(t) ∈ R continue, telle que pour tout x, y ∈ R et
pour tout t ∈ R+ on ait

f (x, t) − f (y, t) (x − y) 6 `(t)|x − y|2 .

Alors le problème de Cauchy (PC) admet une solution globale unique.

20
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Exemple 3.2.1. Soit f : R × R+ −→ R de classe C 1 telle que, il existe K ∈ R :

∂f
(x, t) 6 K, x ∈ R, t > 0.
∂x
Regardons si l’hypothèse (H) est vérifiée. Le théorème des acroissements finis
entraine qu’il existe ξ ∈]x, y[
∂f
f (x, t) − f (y, t) = (ξ, t)(x − y).
∂x
Ainsi nous obtenons
∂f
(ξ, t)(x − y)2 6 K|x − y|2 .


f (x, t) − f (y, t) (x − y) =
∂x
Si on pose `(t) = K, t ∈ R+ , l’hypothèse (H) est bien satisfaite. Le théorème
3.1 permet de conclure que (PC) admet une solution globale unique.
Considérons le problème de Cauchy suivant
(
t2
u0 (t) = −u3 (t) + e− 2 , t > 0
u(0) = 1.
Question : Ce problème, admet il une solution globale unique ?
t2
En effet, si on pose f (x, t) = −x3 + e− 2 , alors
∂f
(x, t) = −3x2 6 0.
∂x
On déduit d’après ce qui précède que (PC) admet une solution globale unique.
Corollaire 3.2.1 (Cauchy-Lipschitz). Si les deux hypothèse suivantes sont vérifiées
i) f : R × R+ −→ R continue.
ii) ∃L ∈ R, ∀x, y ∈ R, ∀t ∈ R+ , on ait
 
f (x, t) − f (y, t) 6 L|x − y|.
Alors (PC) admet une solution globale unique.
t2
Exemple 3.2.2. f (x, t) = |x| + sin(x) + e− 2 , u0 = 1.
   
f (x, t) − f (y, t) = |x| − |y| + sin(x) − sin(y)
  
6 |x| − |y| + | sin(x) − sin(y)
Z 
   x 
6 x − y  +  cosθdθ

 y 
6 2|x − y|.
Le corollaire précédent entraine qu’il existe une solution globale unique pour le
problème de Cauchy (PC).
Remarque 3.2.1. Pour les deux exemples précédents, Exemple 3.1.4 et Exemple
3.1.5, il n’est pas possible de donner une expression explicite de la solution u(t) ;
il est donc nécessaire d’utiliser une méthode pour obtenir des approximations des
valeurs u(t) pour différents t ∈ R+ .

21
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

3.3 Approximation du problème de Cauchy

3.3.1 Schémas d’Euler
Pour établir un schéma d’approximation du problème de Cauchy
 0
u (t) = f (u(t), t), t > 0
(PC)
u(0) = u0

nous supposons que le problème de Cauchy (PC) est bien posé, et nous com-
mençons par partitionner l’axe (temporel) Ot, c. à. d. , nous choisissons des
points t0 , t1 , · · · , tels que

0 = t0 < t1 < t2 < · · · < tn < tn+1 < · · ·

En posant, hn = tn+1 − tn , nous pouvons approcher u0 (tn ) ou u0 (tn+1 ) par

u(tn+1 ) − u(tn )
.
hn
Nous notons dans la suite un ' u(tn ) (un est une approximation de u(tn )).
Ces deux approches, nous suggèrent les schémas suivants :
• Schéma d’Euler Progressif :
 n+1
 u − un
= f (un , tn ), n = 0, 1, 2 · · ·
 0 hn
u = u0 .

• Schéma d’Euler rétrograde :

 n+1
 u − un
= f (un+1 , tn+1 ), n = 0, 1, 2 · · ·
 0 hn
u = u0 .

Ces deux schémas nous permettent de calculer un+1 à partir de un , et il est

donc possible de déterminer successivement u1 , u2 , u3 , · · · , en partant de u0 .
— Le schéma d’Euler progressif est un schéma explicite, car il permet d’ex-
pliciter un+1 en fonction de un :

un+1 = un + hn f (un , tn ).

— Le schéma d’Euler rétrograde est un schéma implicite, car il ne permet

pas d’expliciter directement un+1 en fonction de un lorsque la fonction f
n’est pas triviale

un+1 − hn f (un+1 , tn+1 ) = un .

Dans ce cas pour calculer un+1 , on pose

g(x) = x − hn f (x, tn+1 ) − un

22
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

et nous cherchons un zéro de g(x) en prenant par exemple une méthode

de Newton. On pose x0 = un et xm+1 = xm − gg(x m)
0 (x ) , m = 0, 1, 2, · · · On
m
a
∂f
g 0 (x) = 1 − hn (x, tn+1 ).
∂x
Nous obtenons donc dans ce cas le schéma suivant

x = un ,

 0

xm − hn f (xm , tn+1 ) − un
xm+1 = xm − , m = 0, 1, 2 · · ·
∂f
1 − hn (xm , tn+1 )


∂x
De plus, on peut montrer que xm → un+1 (g(un+1 ) = 0)

3.3.2 Schéma d’Euler progressif vs Schéma d’Euler rétrograde

On prend f (x, t) = −βx, β > 0 donné, alors le problème de Cauchy (PC)
devient
 0
u (t) = −βu(t), t > 0
u(0) = u0 .

— La solution associée à ce problème de Cauchy est donnée par u(t) =

e−βt u0 .
— u(·) décroit exponentiellement vers zéro.
Le schéma d’Euler explicite/progressif associé à ce problème de Cauchy
 n+1
 u − un
= −βun , n = 0, 1, 2, · · ·
 0 hn
u = u0 ,

ce qui est équivalent à

un+1 = un − βhn un , (hn = h, tn = nh)

n
= (1 − βh)u
= (1 − βh)n u0 , n = 0, 1, 2, · · ·

— Bien que u(t) −→ 0, un 6−→ 0 pour u0 6= 0 et |1 − βh| > 1. Il faut

t→+∞ n→+∞
donc chercher une condition pour éviter ce phénomène.
• Si u0 = 0 =⇒ un = 0.

2
• Si u0 6= 0 et |1 − βh| < 1 =⇒ h < β.
2
La condition h < β s’appelle condition de stabilité.

23
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Si on considère maintenant le schéma d’Euler rétrograde

 n+1
 u − un
= −βun+1 , n = 0, 1, 2, · · ·
 0 hn
u = u0 ,

alors

un+1 (1 + βhn ) = un (hn = h, tn = nh)

 1 n
un = u0 , n = 0, 1, 2, · · ·
1 + βh

— Dans ce cas un −→ 0 pour tout h > 0. On dit que le schéma d’Euler

n→+∞
implicite/rétrograde est inconditionnellement stable.

Remarque 3.3.1. Les schémas d’Euler sont d’ordre 1 en h c. à. d., ∃C(h

) > 0n
telle que |u(t) − un | < Ch.

3.3.3 Schémas de Runge Kutta

Dans cette section nous introduisons des nouveaux schémas dits shémas de
Runge-Kutta après deux scientifiques allemands : un mathématicien et physicien
Carl Runge (1856-1927) et un mathématicien Martin Kutta (1867-1944).
Pour bien expliquer la différence entre les schémas d’Euler et les schémas de
Runge-Kutta, on considère notre problème de Cauchy (PC), que l’on suppose
bien posé. On a

u0 (t) = f (u(t), t)
Z tn+1 Z tn+1
0
u (t)dt = f (u(t), t)dt
tn tn
Z tn+1
u(tn+1 ) − u(tn ) = f (u(t), t)dt.
tn

R tn+1
La méthode de Runge Kutta repose sur l’approximation numérique de tn
f (u(t), t)dt.
• Si on utilise la méthode/la formule des trapèzes

1
hn f (un , tn ) + f (un+1 , tn+1 ) ,


u(tn+1 ) − u(tn ) = n = 0, 1, 2, · · ·
2
— Ce schéma est implicite.
— On peut aussi l’obtenir en faisant la moyenne entre le schéma d’Euler
progressif et le schéma d’Euler rétrograde.
— On peut montrer que ce schéma est d’ordre 2 en h, c. à. d.,

|u(t) − un | = O(h2 ).

24
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

— Pour éviter le calcul implicite de un+1 , on peut utiliser une prédiction

d’Euler progressive et remplacer un+1 par ũn+1
ũn+1 = un + hn f (un , tn ).
Nous construisons ainsi un nouveau schéma dit schéma de Heun.
p1 = f (un , tn )
p2 = f (un + hn p1 , tn+1 )
hn
un+1 = un + (p1 + p2 ).
2
— La méthode de Heun fait partie des méthodes de Runge Kutta d’ordre 2,
explicites.
Rt
• Si cette fois ci on approche tnn+1 f (u(t), t)dt par la méthode/la formule
des rectangles
1
un+1 − un = hn f (un+ 2 , tn+ 12 )
— Là encore, on utilise une prédiction d’Euler progressive
1 hn
ũn+ 2 = un + f (un , tn ).
2
— Nous obtenons la méthode d’Euler modifiée
p1 = f (un , tn )
hn hn
p2 = f (un + p1 , tn + )
2 2
un+1 = un + hn p2 .

• Pour obtenir le schéma de Runge Kutta R t classique, dit aussi schéma de

Runge Kutta d’ordre 4, il faut approcher tnn+1 f (u(t), t)dt par la méthode/la
formule de Simpson, puis utiliser la prédiction d’Euler.
On trouve
p1 = f (un , tn )
hn hn
p2 = f (un + p1 , tn + )
2 2
hn hn
p3 = f (un + p2 , tn + )
2 2
p4 = f (un + hn p3 , tn+1 )
hn
un+1 = un + (p1 + 2p2 + 2p3 + p4 ).
6
Remarque 3.3.2. La méthode de Runge Kutta classique est clairement une
méthode explicite. Deplus, c’est une méthode d’ordre 4 en h, c. à. d.,
|u(t) − un | = O(h4 ).

25
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

3.4 Approximation des EDOs d’ordre supérieur

Soit f : R2 × R+ −→ R
une fonction continue, et soient u0 , v0
(x, y, t) 7−→ f (x, t)
donnés. On cherche u de classe C 2 telle que
 00
u (t) = f (u(t), u0 (t), t), t > 0
(PD2)
u(0) = u0 , u0 (0) = v0 .

Il s’agit d’un problème différentiel d’ordre 2. Pour résoudre numériquement le

problème (PD2), on pose v(t) = u0 (t) et donc (PD2) devient
 0
 u (t) = v(t),
v 0 (t) = f (u(t), v(t), t), t > 0 (PD2)
u(0) = u0 , v(0) = v0 .


Remarque 3.4.1. — Pour résoudre ce système, on peut appliquer les méthodes

déjà obtenues dans les sections précédentes pour les systèmes d’ordre 1.
— Il existe des méthodes spécifiquement adaptées aux systèmes différentielles
d’ordre 2. Il s’agit ici de la méthode de Newmark.
On suppose que f (x, y, t) = g(x, t), c. à. d. ,
 00
u (t) = g(u(t), t), t > 0
u(0) = u0 , u0 (0) = v0 .

Le schéma numérique associé à ce problème

 n+1
u − 2un + un−1

 2
= g(un , tn ), t>0
0
h
 u1 = u0 ,

u = u0 + hv0 + 12 h2 g(u0 , 0).

Ce schéma est d’ordre 2 en h, mais comme elle est explicite, elle est stable sous
une condition de stabilité.

26
Chapitre 4

Equations aux dérivées

partielles

4.1 Approximation des équations aux dérivées

partielles
Les équations aux dérivées partielles (EDPs) apparaissent dans tous les do-
maines des mathématiques appliquées, pour décrire l’évolution en temps de cer-
tains quantités qui dépendent du temps. Puisque ces EDPs peuvent avoir des
formes complexes, il est souvent difficile, et parfois impossible de les résoudre
analytiquement. Pour cela, nous avons besoin de recourir aux schémas numériques,
pour obtenir une approximation de la solution exacte.
Dans ce cours nous nous intéressons aux méthodes d’approximations par
rapport à l’espace, c. à. d. , aux procédures numériques qui approchent les
opérateurs différentiels qui figurent dans ces EDPs.

Remarque 4.1.1. Il maintenant bien clair que avant de passer à la résolution

numérique d’un système d’équations différentielles, il faut s’assurer que ce système
est bien posé. Pour les EDOs il existe le théorème de Cauchy-Lipschitz pour
l’existence et l’unicité du problème de Cauchy. Pour les EDPs la situation est
différente, cependant, il existe plusieurs méthodes et résultats qui traitent le
problème de well-posedness des EDPs. Je réfère le lecteur intéressé par ces ques-
tions à [4, 7], voir aussi [3].

4.1.1 Exemples des EDPs classiques

Exemple 4.1.1 (Equation de la chaleur).

 ∂t w(t, x) = ∂xx w(t, x), x ∈ (0, π)
w(t, 0) = w(t, π) = 0
w(0, x) = w0 (x).


27
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Exemple 4.1.2 (Equation de transport).


 ∂t w(t, x) = ∂x w(t, x), x ∈ (0, π)
w(t, 0) = 0
w(0, x) = w0 (x).


Exemple 4.1.3 (Equation des ondes).


 ∂tt w(t, x) = ∂xx w(t, x), x ∈ (0, π)
w(t, 0) = ∂x w(t, 0) = 0
w(0, x) = w0 (x).


L’idée principale pour une discrétisation simple des EDPs est la suivante
1. On discrétise l’opérateur différentiel spacial par rapport à la variable spa-
ciale x.
2. Après obtention du système d’EDOs, on fait une nouvelle discrétisation
par rapport au temps ou bien on utilise les méthodes de discrétisation
temporelle.

Remarque 4.1.2. Dans 1., on parle de semi-discrétisation par rapport à l’es-

pace x ou de la méthode des lignes.
Dans la suite nous allons présenter deux façons d’approcher un opérateur
différentiel spatial et nous allons introduire deux classes des méthodes de discrétisation
spatiale, à savoir, la méthode des différences finies et la méthode de Ga-
lerkin.

4.1.2 Méthode des différences finies

Soit Ω = (a, b) un domaine spatial, et soit (xj )j une subdivision de Ω. Les
points de la grille xj = a + j∆x, j = 0, 1, 2, · · · , N , avec ∆x = b−a
N .
On rappelle que

wj+1 (t) − wj (t)

∂x w(t, xj ) ≈
∆x
ou
wj (t) − wj−1 (t)
∂x w(t, xj ) ≈
∆x
wj+1 (t) − 2wj (t) + wj−1
∂xx w(t, xj ) ≈
(∆x)2

pour j = 1, 2, · · · , N − 1, wj (t) ≈ w(t, xj ).

• Discrétisation de l’équation de la chaleur :

d wj+1 (t) − 2wj (t) + wj−1

wj (t) = , j = 1, 2, · · · , N − 1.
dt (∆x)2

28
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Les cas de j = 0 et j = N sont déjà donnés par les conditions au bord, de

Dirichlet, w(t, 0) = w(t, π) = 0 par w0 (t) = wN (t) = 0, ∀t > 0. On obtient
ainsi un système d’équations différentielles ordinaires qui peut se ramener sous
la forme suivante
d
(
W (t) = M W (t), t > 0
dt
W (0) = W0

t
avec W (t) = w1 (t), w2 (t), · · · , wN −1 (t) ∈ RN −1 et
−2 1
 
0
 1 −2 . . .
 

1  
.. .. ..

M=  ∈ RN −1 × RN −1 .

(∆x)2 
 . . . 
 .. 
 . −2 1 
0 1 −2
M est une matrice creuse, tridiagonale, symétrique, définie négative. On peut
l’écrire aussi
1
M= tridiag(1, −2, 1).
(∆x)2
• Discrétisation de l’équation de transport :

d wj+1 (t) − wj (t)

wj (t) = , j = 1, 2, · · · , N
dt (∆x)2
avec w0 (t) = 0 ∀t > 0.

t 1
En introduisant W (t) = w1 (t), w2 (t), · · · , wN (t) ∈ RN et M = ∆x tridiag(−1, 1, 0) ∈
RN × RN , on obtient le système d’EDOs suivant
d
(
W (t) = M W (t), t > 0
dt
W (0) = W0

4.1.3 Méthode de Galerkin

Contrairement à la méthode des différences finies qui approche la solution
exacte en certains points de la grille, la méthode de Galerkin utilise une com-
binaison linéaire de quelques fonctions de base. La procédure de la méthode de
Galerkin est la suivante
— Soit (H, < ·, · >) un espace de Hilbert muni d’un produit scalaire.
— A : D(A) ⊂ H −→ H, f ∈ H. On considère le problème suivant

∃?u ∈ D(A)
(P)
f = Au in H.
Par exemple Au = u00 sur H = L2 (0, π).

29
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

— Soit Hm ⊂ D(A) ⊂ H un sous espace vectoriel, tel que dim(Hm ) = mn

et soient {ϕ1 , ϕ2 , · · · , ϕm } ⊂ Hm , une base de Hm .
Pm
Pour tout um ∈ Hm : um = k=1 ck ϕk , ck ∈ C, k = 1, 2, · · · , m.
L’idée principal de la méthode de Galerkin est que la solution exacte u ∈ H est
approchée par une suite (um ) ⊂ Hm quand m → ∞.
Revenons au problème (P) : ∃?u ∈ D(A) telle que f = Au in H.
On prend le produit scalaire des deux membres par ϕj , j = 1, 2, · · · , m.

< f, ϕj >=< Au, ϕj >, j = 1, 2, · · · , m.

On définit la méthode de Galerkin, en remplaçant u par um dans l’expression

précédente et on projète um dans la base de Hm .

< f, ϕj > =< Aum , ϕj >

m
X
= ck < Aϕk , ϕj >, j = 1, 2, · · · , m.
k=1
   
< f, ϕ1 > c1
 ·   · 
On pose Φ =  , C =   et
 ·   · 
< f, ϕm > cm

· · ·
 
< Aϕ1 , ϕ1 > < Aϕm , ϕ1 >
 .. 
 < Aϕ1 , ϕ2 >
 · . · 

Am =
 .. .. .. .

 · . . . · 
 .. 
 · . · 
< Aϕ1 , ϕm > · < Aϕm , ϕm >

Ce qui donne Am C = Φ.
L’algorithme de la méthode Galerkin est le suivant
1. Choisir un sous espace vectoriel Hm ⊂ D(A) ⊂ H, avec dimHm = m et
une base {ϕ1 , ϕ2 , · · · , ϕm } ⊂ Hm .
2. Résoudre Am C = Φ dont la solution est C = (c1 , c2 , · · · , cm )t .
Pm
3. L’approximation um de u est construite comme um = k=1 ck ϕk .
Sous certaines hypothèses on peut montrer que um −→ u dans H.
m→+∞
Question : Comment choisir les fonctions de base {ϕj }, j = 1, 2, · · · , m ? En
pratique, il existe deux possibilités : la base {ϕj } est choisie comme étant
(i) les vecteurs propres de l’opérateur A.
(ii) les fonctions splines/chapeaux, fonctions continues par morceaux.
Autrement dit, on peut utiliser soit
→ la méthode spectrale,

30
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

→ la méthode des éléments finis.

Exemple 4.1.4 (Equation de Poisson). On considère l’équation de Poisson en
dimension 1, sur L2 (0, π).

 d2
w(x) = f (x), x ∈ (0, π)
dx2
 w(0) = w(π) = 0.

d2
On a A = , D(A) = H 2 (0, π) ∩ H01 (0, π).
dx2
Pour résoudre cette équation, en utilise la méthode de Galerkin. On commence
par choisir Hm ⊂ D(A), une base {ϕ1 , ϕ2 , · · · , ϕm } ⊂ Hm telle que dimHm =
m.
Xm
∀wm ∈ Hm , wm (x) = ck ϕk (x)
k=1

ce qui implique
m
X
ck < Aϕk , ϕj >L2 =< f, ϕj >L2
k=1
m
X Z π Z π


ck Aϕk (x)ϕj (x)dx = f (x)ϕj (x)dx
k=1 0 0
m
X Z π Z π
(−ck ) ϕ0k (x)ϕ0j (x)dx = f (x)ϕj (x)dx (On a utilisé les conditions au bord)
k=1 0 0

Méthode spectrale
On choisit Hm = vect{sin(jx), j = 1, 2, · · · , m}. Donc
m
X Z π m
X Z π
(−ck ) ϕ0k (x)ϕ0j (x)dx = (−ck jk) cos(kx) cos(jx)dx
k=1 0 k=1 0
m
X π 1

= (−ck jk)δjk cos(a) cos(b) = [cos(a − b) + cos(a + b)]
2 2
k=1
π
= −cj j 2
Z π 2
= f (x)ϕj (x)dx.
0

Par suite
Z π
2
cj = − 2 f (x)ϕj (x)dx
j π 0
m Z π
2X 1
wm = − sin(kx) f (s)sin(ks)ds
π k2 0
k=1

31
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

= −jk π2 δjk .


Dans ce cas Am jk

Méthode des éléments finis

Une approche est de choisir les fonctions de base {ϕj }, j = 1, 2, · · · , m
comme étant les fonctions splines/chapeaux suivantes :

 0 si x < (j − 1)∆x = xj−1 ,
 x − (j − 1) si xj−1 < x < xj ,



ϕj (x) = ∆x x
 (j + 1) − si xj < x < xj+1 ,
∆x



0 si x > xj+1 .


Ici ∆x = π
m, m ∈ N∗ .


0 si x < (j − 1)∆x = xj−1 ,
d 1/∆x

si xj−1 < x < xj ,
ϕj (x) =
dx 
 −1/∆x si xj < x < xj+1 ,
0 si x > xj+1 .


D’après les calculs précédents

m Z π
X  −1 2 −1 
(−ck ) ϕ0k (x)ϕ0j (x)dx = ∆x − cj−1 − cj − c j+1
0 (δx)2 (δx)2 (δx)2
k=1
Z π
1
= (cj−1 − 2cj + cj+1 ) = f (x)ϕj (x)dx, j = 1, · · · , m.
∆x 0
1
Dans ce cas Am = ∆x tridiag(1, −2, 1).
Exercices :
1. Utiliser la formule de Simpson et montrer que le schéma de Runge Kutta
est d’ordre 4.
2. Soit le système suivant
 0
u (t) = −k 2 u(t), t > 0
(S)
u(0) = 1.

a) Calculer la solution u.
b) Utiliser le schéma d’Euler explicite et le schéma d’Euler implicite pour
avoir une approximation de la solution de (S).
c) Quel schéma préserve les propriétés de la solution exacte de (S).
3. On considère l’équation de la chaleur sur Ω = (0, π) × (0, π), t > 0.

∂t w(t, x, y) = ∂xx w(t, x, y) + ∂yy w(t, x, y)

avec les conditions au bord : w(t, x, y) = 0 sur ∂Ω.

On pose L = ∂xx + ∂yy .
Utiliser la méthode des différences finies pour approcher l’opérateur L.

32
ANALYSE NUMERIQUE 2 SALEM NAFIRI 2017

Conclusion
Introduire les méthodes numériques n’est pas une tâche facile en général.
Puisque leurs utilisations nécessite la connaissance préalable d’outils de base
tels que l’analyse de l’erreur, résolution numérique des équations linéaires et
non linéaires, intégration numérique, etc...Pour plus de détails sur ces outils de
base, voir [1, 5, 6]. Une grande partie de ces notes se sont inspirées des livres [2]
et [6]. Quelques résultats que nous n’avons pas démontrés peuvent être trouvés
également dans [6].

33
Bibliographie

[1] G. Allaire. Analyse numérique et optimisation. Edité par les Editions de

l’Ecole Polytechnique, 2012.
[2] A. Batkai, B. Farkas, P. Csomos, A. Ostermann, Operator Semigroups for
Numerical Analysis, 15th Internet Seminar 2011/12.
[3] H. Brezis, Analyse fonctionnelle : Théorie et applications. Paris : Masson,
1983. (Mathematiques Appliquées pour la Maitrise.)
[4] K. Engel and R. Nagel, One-parameter semigroups for linear evolution
equations, Encyclopedia of Mathematics and its Applications. Springer-
Verlag, New York (2000).
[5] A. Fortin, Analyse numérique pour ingénieurs, Editions de l’école polytech-
nique de Montréal, 1995.
[6] J. Rappaz, M. Picasso, Introduction à l’analyse numérique, PPUR
(Presses Polytechniques et Universitaires Romandes), Lausanne, 1998.
http ://dmawww.epfl.ch/rappaz.mosaic
[7] A. Pazy, Semigroups of linear operators and applications to partial diffe-
rential equations 44 of Applied Math. Sciences. Springer-Verlag, New York
(1983).

34