Notes sur les méthodes numériques pour la dynamique des fluides - 1° partie, Notes de Économie
Tina_920
Tina_92024 mars 2014

Notes sur les méthodes numériques pour la dynamique des fluides - 1° partie, Notes de Économie

PDF (909 KB)
46 pages
160Numéro de visites
Description
Notes d’économie sur les méthodes numériques pour la dynamique des fluides - 1° partie. Les principaux thèmes abordés sont les suivants: Introduction et classication des EDP, EDP elliptiques linéaires, EDP paraboliques ...
20 points
Points de téléchargement necessaire pour télécharger
ce document
Télécharger le document
Aperçu3 pages / 46

Ceci c'est un aperçu avant impression

3 shown on 46 pages

Télécharger le document

Ceci c'est un aperçu avant impression

3 shown on 46 pages

Télécharger le document

Ceci c'est un aperçu avant impression

3 shown on 46 pages

Télécharger le document

Ceci c'est un aperçu avant impression

3 shown on 46 pages

Télécharger le document

Université de Lorraine Master 2 IMOI Ingénierie Mathématique et Outils Informatiques Option Calcul Scientique

Méthodes numériques pour la dynamique des uides

J.-F. Scheid

Année 2011-2012

2

Table des matières

1 Introduction et classication des EDP 5

1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5

1.2 EDP linéaires du second ordre . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5

1.3 EDP du premier ordre . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

1.4 Problème bien posé . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7

1.5 Classication des EDP du second ordre . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7

2 EDP elliptiques linéaires 11

2.1 Introduction - Propriétés des solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11

2.1.1 Quelques rappels sur l'existence, l'unicité et la régularité des solutions . . . . . . . 11

2.1.2 Principes du maximum . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12

2.2 Diérences nies pour le cas 1D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

2.2.1 Erreur de consistance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

2.2.2 Matrices monotones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

2.2.3 Stabilité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

2.2.4 Convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

2.2.5 Autres conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

2.3 Diérences nies pour le cas 2D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

2.3.1 Un schéma à 5 points pour le Laplacien . . . . . . . . . . . . . . . . . . . . . . . . 23

2.3.2 Un autre schéma à 5 points . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

2.3.3 Un schéma à 9 points pour fonction harmonique . . . . . . . . . . . . . . . . . . . 27

2.3.4 Cas d'un domaine non-rectangulaire . . . . . . . . . . . . . . . . . . . . . . . . . . 28

2.3.5 Opérateur sous forme de divergence . . . . . . . . . . . . . . . . . . . . . . . . . . 33

2.3.6 Autres conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

2.4 Evaluation pratique de l'ordre de convergence d'une méthode . . . . . . . . . . . . . . . . 36

2.5 Méthode de Richardson . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36

2.6 Conditionnement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

3 EDP paraboliques linéaires 39

3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

3.1.1 Existence et unicité des solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

3.1.2 Principes du maximum . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

3.2 Equation de la chaleur en dimension 1 d'espace . . . . . . . . . . . . . . . . . . . . . . . . 40

3.2.1 Schéma d'Euler explicite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41

3.2.2 Schéma d'Euler implicite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44

3.2.3 Schéma de Crank-Nicholson (1947) . . . . . . . . . . . . . . . . . . . . . . . . . . . 48

3.2.4 θ-schéma pour l'équation de la chaleur . . . . . . . . . . . . . . . . . . . . . . . . . 50

3.2.5 Autres schémas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50

3.3 Cas de la dimension 2 d'espace . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

3.3.1 θ-schéma pour l'équation de la chaleur en 2D . . . . . . . . . . . . . . . . . . . . . 51

3.3.2 Directions alternées . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52

1

2 TABLE DES MATIÈRES

4 EDP hyperboliques linéaires 57

4.1 Equation des ondes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

4.1.1 Existence, unicité et propriétés des solutions . . . . . . . . . . . . . . . . . . . . . . 57

4.1.2 Un θ-schéma centré pour l'équation des ondes (1D) . . . . . . . . . . . . . . . . . . 59

4.2 Equation de transport . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

4.2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

4.2.2 Schéma centré . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63

4.2.3 Schéma de Lax . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64

4.2.4 Schéma décentré (upwind) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68

4.2.5 Schéma de Lax-Wendro . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71

4.2.6 Comparaison des diérents schémas . . . . . . . . . . . . . . . . . . . . . . . . . . 75

4.2.7 Quelques schémas pour le 2D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77

5 EDP hyperboliques non-linéaires - Lois de conservation 81

5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81

5.2 Solutions classiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82

5.3 Solutions faibles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83

5.4 Relations de Rankine-Hugoniot . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84

5.5 Solutions d'entropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87

5.6 Problème de Riemann . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

5.6.1 Cas où f est strictement convexe . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

5.6.2 Cas général . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93

5.7 Schémas d'approximations aux Diérences Finies . . . . . . . . . . . . . . . . . . . . . . . 96

5.7.1 Introduction et généralités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

5.7.2 Schéma de Godounov . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97

5.7.3 Schéma d'Engquist-Osher . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

6 Equations de Stokes 101

6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101

6.2 Adimensionalisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102

6.3 Réductions des équations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102

6.4 Discrétisation des équations de Stokes par Diérences Finies . . . . . . . . . . . . . . . . . 103 6.4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103

6.4.2 Schéma MAC pour le problème de Stokes . . . . . . . . . . . . . . . . . . . . . . . 104

6.4.3 Forme matricielle du schéma de MAC . . . . . . . . . . . . . . . . . . . . . . . . . 105

6.5 Résolution du système discrétisé de Stokes . . . . . . . . . . . . . . . . . . . . . . . . . . . 107

6.6 Conditions de Dirichlet non-homogènes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108

6.7 Traitement pratique des conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

6.8 Exemples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111

7 Equations de Navier-Stokes 115

7.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115

7.2 Semi-discrétisation en temps . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115

7.3 Discrétisation totale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116

7.4 Forme matricielle du schéma semi-implicite . . . . . . . . . . . . . . . . . . . . . . . . . . 118

7.5 Résolution du système discrétisé de Navier-Stokes . . . . . . . . . . . . . . . . . . . . . . . 120

7.6 Conditions de Dirichlet non-homogènes - Traitement des conditions limites . . . . . . . . . 120

7.7 Exemples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120

TABLE DES MATIÈRES 3

8 Volumes Finis 123 8.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 8.2 Volumes Finis pour les problèmes elliptiques 1D . . . . . . . . . . . . . . . . . . . . . . . . 123

8.2.1 Maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 8.2.2 Formulation en Volumes Finis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 8.2.3 Système linéaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 8.2.4 Convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 8.2.5 Equation elliptique 1D à coecients discontinus . . . . . . . . . . . . . . . . . . . . 126

8.3 Volumes Finis pour les problèmes elliptiques 2D . . . . . . . . . . . . . . . . . . . . . . . . 127 8.3.1 Maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 8.3.2 Formulation en Volumes Finis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 8.3.3 Exemples de maillages admissibles . . . . . . . . . . . . . . . . . . . . . . . . . . . 130 8.3.4 Système linéaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 8.3.5 Equation elliptique 2D avec coecients discontinus . . . . . . . . . . . . . . . . . . 132

8.4 Volumes Finis pour l'équation de transport . . . . . . . . . . . . . . . . . . . . . . . . . . 133 8.4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133 8.4.2 Maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134 8.4.3 Formulation en Volumes Finis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135 8.4.4 Système linéaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136 8.4.5 Condition de stabilité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138

8.5 Volumes Finis pour les équations de Stokes . . . . . . . . . . . . . . . . . . . . . . . . . . 139 8.5.1 Maillages . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 8.5.2 Formulation en Volumes Finis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139

8.6 Volumes Finis pour les équations de Navier-Stokes incompressibles . . . . . . . . . . . . . 142

Références 145

Exercices 147

TP 159

4 TABLE DES MATIÈRES

Chapitre 1

Introduction et classication des EDP

1.1 Introduction

On s'intéresse aux équations aux dérivéees partielles sous la forme générale

F (x, u,Du, · · · , Dαu) = 0, (1.1)

où u est une fonction inconnue des N variables regroupées dans le vecteur x = (x1, · · · , xN ) ; α est un

multi-indice α = (α1, · · · , αN ) et Dα = ∂|α|

∂xα11 · · · ∂x αN N

avec |α| = α1 + α2 + · · ·+ αN .

L'équation (1.1) est du premier ordre si |α| = 1 et du deuxième ordre si |α| = 2.

L'équation est dite :

• linéaire si F est linéaire en u,Du, · · · , Dαu ( u et v solutions ⇒ λu+ βv solution). • semi-linéaire si F est linéaire en Du, · · · , Dαu. • quasi-linéaire si F est linéaire en Dαu. • non-linéaire si F n'est pas linéaire en au moins une dérivée.

1.2 EDP linéaires du second ordre

Donnons pour commencer, quelques exemples fondamentaux d'EDP linéaires du second ordre qu'on étudiera dans les chapitres suivants.

Equation de Laplace.

En dimension 2, on cherche u = u(x, y) vériant

uxx + uyy = f(x, y), pour (x, y) ∈ (0, 1)× (0, 1), (1.2)

avec les conditions limites u(0, y) = φ0(y), u(1, y) = φ1(y), u(x, 0) = ψ0(x), u(x, 1) = ψ1(x),

(1.3)

les fonctions f, φ0, φ1, ψ0, ψ1 étant données.

Plus généralement en dimension n, on cherche u = u(x1, · · · , xn) telle que{ −∆u = f, dans un ouvert Ω ⊂ Rn

u = g, sur ∂Ω

où ∆ désigne l'opérateur de Laplace et g est une fonction donnée sur le bord. Cette équation modélise par exemple (de façon très simpliée...) le déplacement d'une membrane soumise à une force extérieure f avec un déplacement g imposé sur le bord ∂Ω. L'équation avec f ≡ 0 est appelée équation de Poisson et

5

6 CHAPITRE 1. INTRODUCTION ET CLASSIFICATION DES EDP

correspond à la recherche des fonctions harmoniques. L'équation de Poisson peut être obtenue de la façon suivante. En l'absence de force extérieure (f ≡ 0) on peut écrire la loi de conservation du ux sur tout le bord fermé ∂Σ d'un sous-domaine quelconque Σ, c'est-à-dire∫

∂Σ ∇u · n dσ = 0,

où n désigne la normale extérieure au domaine Σ. La formule de Green donne alors∫ Σ

∆u dx = 0.

Ceci étant vrai pour tout sous-domaine Σ (et u susamment régulière), on en déduit ∆u = 0 dans Ω.

Equation de la chaleur.

Pour un domaine Ω ⊂ Rn et T > 0, on cherche une fonction u = u(x1, · · · , xn, t) dépendante du temps t telle que 

∂u

∂t −∆u = f, dans Ω× (0, T )

u = g, sur ∂Ω× (0, T ) u(x, t = 0) = u0(x), x ∈ Ω.

La fonction u représente par exemple la température en un point x du domaine Ω et à l'instant t ∈ [0, T ].

Equation des ondes.

On cherche u = u(x, y) (y représente le temps) vériant

uxx − uyy = f(x, y), pour (x, y) ∈ (0, 1)× (0, T ), (1.4)

avec les conditions limites u(0, y) = φ0(y), u(1, y) = φ1(y), u(x, 0) = ψ0(x), uy(x, 0) = ψ1(x),

(1.5)

ou bien plus généralement, pour un domaine Ω ⊂ Rn on cherche u = u(x1, · · · , xn, t) telle que

∂2u

∂t2 −∆u = f, dans Ω× (0, T )

u = g, sur ∂Ω× (0, T ) u(x, t = 0) = u0(x),

∂u

∂t (x, t = 0) = v0(x), x ∈ Ω.

1.3 EDP du premier ordre

On s'intéressera également aux EDP du premier ordre à partir de quelques exemples fondamentaux.

Equation de transport.

En dimension 1 d'espace, on cherche u = u(x, t) avec x ∈ R et t ∈ R+, vériant

∂u

∂t + c

∂u

∂x = 0. (1.6)

La fonction u représente par exemple une quantité en x et à l'instant t transportée par une vitesse c ∈ R.

Equation de Bürgers.

En dimension 1 d'espace, on cherche u = u(x, t) avec x ∈ R et t ∈ R+, vériant

∂u

∂t +

∂x

( u2

2

) = 0. (1.7)

1.4. PROBLÈME BIEN POSÉ 7

Il s'agit d'un modèle simplié de la dynamique des gaz. Cette équation sert aussi de modèle du bang sonique : le bruit engendré par un avion supersonique, loin de l'avion (près du sol), se concentre dans certaines zones où la pression est gouvernée par l'équation de Bürgers.

Lois de conservation.

Les deux équations précédentes sont des cas particuliers de lois de conservation qui s'écrivent plus générale- ment (toujours en dimension 1 d'espace)

∂u

∂t +

∂x (f(u)) = 0. (1.8)

Ces équations sont quasi-linéaires car linéaire par rapport aux dérivées d'ordre le plus élévé (ut et ux).

1.4 Problème bien posé

On dit qu'un problème est bien posé au sens d'Hadamard si sa solution existe, est unique et dépend continûment des données.

Les conditions aux limites jouent un rôle essentiel quant au caractère bien posé ou non d'un problème. Par exemple, l'équation de Laplace avec une condition de Neumann sur le bord (∂u/∂n = 0 sur ∂Ω) n'est pas bien posé en général (il faut des conditions de compatibilité sur la fonction f). De même l'équation de Laplace (1.2) avec les conditions aux limites (1.5) n'est pas bien posé (pas de continuité par rapport aux données), ni l'équation des ondes (1.4) avec les conditions aux limites (1.3)

1.5 Classication des EDP du second ordre

Comme on le verra plus tard, les solutions d'EDP du second ordre ont des propriétés diérentes selon le type d'équations. Dans un domaine Q ⊂ RN , on considère pour u = u(x), l'équation linéaire aux dérivées partielles du second ordre de la forme

N∑ i,j=1

aij(x) ∂2u

∂xi∂xj +

N∑ i=1

bi(x) ∂u

∂xi + c(x)u = f, (1.9)

avec x = (x1, · · · , xN ) ∈ Q. Les coecients aij sont réels et la matrice A formée des coecients aij est supposée symétrique( 1). La fonction f est donnée.

Pour un point x0 xé dans Q, on désigne par N+ = N+(x0) le nombre de valeurs propres de A(x0) strictement positives, par N− = N−(x0) le nombre de celles strictement négatives et par N0 = N0(x0) le nombre de celles nulles (N+ +N− +N0 = N).

L'équation (1.9) est dite

• elliptique en x0 si N+ = N ou N− = N , • parabolique en x0 si N0 > 0, • hyperbolique en x0 si (N+ = N − 1 et N− = 1) ou bien (N− = N − 1 et N+ = 1),

L'équation (1.9) est dite elliptique, parabolique ou hyperbolique en O ⊂ Q, si elle jouit de la propriété en tout point x ∈ O.

1. L'hypothèse de symétrie de la matrice A est raisonnable au sens où l'on peut écrire ∑ i,j

aij ∂2u

∂xi∂xj =

∑ i,j

a′ij ∂2u

∂xi∂xj +

∑ i,j

a′′ij ∂2u

∂xi∂xj , avec a′ij = (aij +aji)/2 et a

′′ ij = (aij−aji)/2. Or

∑ i,j

a′′ij ∂2u

∂xi∂xj = 0 car si u ∈ C2(Q), on a ∂

2u

∂xi∂xj =

∂2u

∂xj∂xi .

Par conséquent, on obtient ∑ i,j

aij ∂2u

∂xi∂xj =

∑ i,j

a′ij ∂2u

∂xi∂xj et la matrice formée des coecients a′ij est symétrique.

8 CHAPITRE 1. INTRODUCTION ET CLASSIFICATION DES EDP

Dans le cas où N est pair, on en déduit que l'équation (1.9) est elliptique, parabolique ou hyperbolique si on a respectivement, det(A) > 0, det(A) = 0 ou det(A) < 0.

Donnons à présent une interprétation géométrique de la classication précédente. On introduit la fonction F : RN → R dénie par F (x) = (Ax,x) où (·, ·) désigne le produit scalaire dans RN . L'équation xN+1 = F (x1, ..., xN ) est l'équation d'une ellipse, d'une parabole ou d'une hyperbole en dimension N + 1, selon que l'équation (1.9) est respectivement elliptique, parabolique ou hyperbolique. En eet, la matrice A étant diagonalisable, on a la décomposition A = QTDQ où D est la matrice diagonale formée des valeurs propres λi (i = 1, ..., N) de A et Q est la matrice orthogonale dont les colonnes sont formées des vecteurs propres. Ainsi F (x) = (Ax,x) = (DQx, Qx) =

∑N i=1 λiz

2 i avec zi = (Qx)i et en fonction du

signe des valeurs propres, l'équation xN+1 = F (x1, ..., xN ) est celle d'une ellipse, d'une parabole ou d'une hyperbole.

Exemples.

Considérons le cas N = 2. Les courbes de niveaux de F sont des ellipses, des hyperboles ou des droites, selon que l'équation (1.9) est respectivement elliptique, hyperbolique ou parabolique . Si la matrice A est dénie positive (det(A) > 0) alors l'équation (1.9) est elliptique et les courbes de niveaux de F sont

des ellipses. Par exemple avec la matrice A =

( 4 2 2 3

) dont les valeurs propres sont λ1 ' 1.4384472 et

λ2 ' 5.561552, on obtient les ellipses de la gure 1.1. Les axes principaux sont déterminés par les vecteurs propres de A.

−1.224 −0.874 −0.524 −0.174 0.176 0.525 0.875 1.225

−0.941

−0.768

−0.595

−0.422

−0.249

−0.075

0.098

0.271

0.444

0.617

0.790

0.2 0.4 0.6

0.8

Figure 1.1  Courbes de niveaux (Ax,x) = cste avec det(A) > 0 (N+ = N = 2).

La gure 1.2 montre un exemple hyperbolique lorsque la matrice A a toutes ses valeurs propres positives

sauf une strictement négative (det(A) < 0). On a choisi la matrice A =

( 2 3 3 1

) dont les valeurs propres

sont λ1 ' −1.5413813 et λ2 ' 4.5413813.

1.5. CLASSIFICATION DES EDP DU SECOND ORDRE 9

−1.500 −1.071 −0.643 −0.214 0.214 0.643 1.071 1.500

−1.060

−0.848

−0.636

−0.424

−0.212

0.000

0.212

0.424

0.636

0.848

1.060

−0.8

−0.8

−0.6

−0.6

−0.4

−0.4

−0.2

−0.2

0.00.0

0.2

0.2

0.4

0.4

0.6

0.6

0.8

0.8

Figure 1.2  Courbes de niveaux (Ax,x) = cste avec det(A) < 0 (N+ = N − 1 = 1 et N− = 1).

La gure 1.3 montre un exemple parabolique lorsque la matrice A est singulière (det(A) = 0). On a

choisi la matrice A =

( 2 1 1 0.5

) dont les valeurs propres sont λ1 = 0 et λ2 = 2.5.

−1.415 −1.132 −0.849 −0.566 −0.283 0.000 0.283 0.566 0.849 1.132 1.415

−1.0

−0.8

−0.6

−0.4

−0.2

0.0

0.2

0.4

0.6

0.8

1.0

0.1 0.10.2 0.20.3 0.30.4 0.40.5 0.5

0.6

0.6

Figure 1.3  Courbes de niveaux (Ax,x) = cste avec det(A) = 0 (N0 = 1).

Remarque : En général, une équation peut avoir diérents types dans diérentes parties du domaine Q. Par exemple l'équation de Tricomi

yuxx + uyy = f

est elliptique pour y > 0, parabolique pour y = 0 et hyperbolique pour y < 0.

10 CHAPITRE 1. INTRODUCTION ET CLASSIFICATION DES EDP

Chapitre 2

EDP elliptiques linéaires

2.1 Introduction - Propriétés des solutions

Pour un domaine (ouvert connexe) Ω ⊂ Rn borné et de frontière ∂Ω régulière, on considère le problème aux limites suivant pour une fonction u :

Lu := −div (A∇u) + c u = f dans Ω (2.1) u = 0 sur ∂Ω (2.2)

où A est une matrice de taille n× n avec A = A(x) = (aij(x))1≤i,j≤n pour x ∈ Ω. On suppose dans tout ce chapitre (sauf précision contraire) que les coecients aij ∈ C1(Ω). Les fonctions f et c sont données et on suppose également dans tout ce chapitre, que la fonction c ∈ C(Ω) est positive ou nulle i.e. c(x) ≥ 0, ∀x ∈ Ω. L'opérateur L apparaît ici sous forme divergentielle et on doit lire

div (A∇u) = n∑

i,j=1

∂xi

( aij

∂u

∂xj

) .

On dit que l'équation (2.1) est (uniformément) elliptique ou bien que l'opérateur L est elliptique si

∃α > 0 tel que (A(x)ξ, ξ)Rn = n∑

i,j=1

aij(x)ξiξj ≥ α|ξ|2, ∀ξ ∈ Rn, ∀x ∈ Ω. (2.3)

On supposera désormais que la condition d'ellipticité (2.3) est vériée.

Remarque. Si A vérie la condition (2.3) alors l'équation (2.1) est elliptique au sens de la classication faite au chapitre 1 (cf. section 1.5). En eet, soit λ une valeur propre de A et u 6= 0 un vecteur propre associé. On a d'une part (Au, u) ≥ α|u|2 et d'autre part (Au, u) = λ(u, u) = λ|u|2. On en déduit que λ ≥ α > 0 donc toutes les valeurs propres de A sont strictement positives. La condition d'ellipticité (2.3) implique en fait que la matrice A est dénie positive.

2.1.1 Quelques rappels sur l'existence, l'unicité et la régularité des solutions

Pour les quelques résultats d'existence, d'unicité et de principe du maximum rappelés ci-après, on pourra consulter par exemple [1], [3], [5], [4].

Commençons par le cas de la dimension 1 avec Ω = (0, 1). Dans ce cas, le problème s'écrit

− ( au′ )′

+ c u = f dans (0, 1) (2.4)

u(0) = u(1) = 0. (2.5)

• Si f, c ∈ C([0, 1]) avec c ≥ 0 et si a ∈ C1([0, 1]) avec a(x) ≥ α > 0 pour tout x ∈ (0, 1)( 1), alors le problème (2.4)-(2.5) admet une unique solution classique u ∈ C2([0, 1]).

1. Cette condition ne traduit rien d'autre que la condition d'ellipticité (2.3) dans le cas 1D.

11

12 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

Considérons à présent le cas de la dimension quelconque avec Ω ⊂ Rn un domaine borné et régulier. On rappelle tout d'abord un résultat d'existence de solution faible( 1).

• Si f ∈ L2(Ω), c ∈ C(Ω) avec c ≥ 0 et si les coecients aij ∈ C(Ω) vérient (2.3), alors le problème (2.1)-(2.2) admet une unique solution faible u ∈ H10 (Ω) (Lax-Milgram) et ‖u‖H1 ≤ C‖f‖L2 où C est une constante indépendante de u et f .

Si de plus, aij ∈ C1(Ω) alors u ∈ H2(Ω)( 2) et u vérie les équations (2.1)-(2.2) au sens presque partout.

Donnons enn un résultat de régularité dans les espaces de Hölder. Pour 0 < α < 1, on dénit les espaces

C0,α(Ω) = {u ∈ C0(Ω), sup x 6=y

|u(x)− u(y)| |x− y|α

<∞}

et pour m ∈ N, Cm,α(Ω) = {u ∈ Cm(Ω), Dβu ∈ C0,α(Ω) avec |β| = m}.

• Si aij ∈ Ck+1,α(Ω) vérient (2.3) et f , c ∈ Ck,α(Ω) avec c ≥ 0, alors il existe une unique solution u ∈ Ck+2,α(Ω).

Ces résultats indiquent un eet régularisant d'un opérateur elliptique, par rapport aux données. Par exemple pour l'opérateur Laplacien (A = Id), si on prend un second membre f aléatoire (entre -100 et 100) sur le pavé unité Ω = (0, 1)× (0, 1) on obtient une solution plus régulière, comme l'illustre la gure 2.1. Si à présent on prend cette solution comme second membre de l'équation de Laplace, on obtient une solution encore plus régulière (cf. Fig. 2.1).

Remarques sur la régularité des solutions.

L'étude de la convergence d'approximation de solution doit tenir compte de la régularité de la solution exacte. Si la solution est faible (L2, H1, · · · ), il faut regarder la convergence dans L2, H1, · · · ; si la solution est régulière (C2, · · · ), alors on peut regarder la convergence dans L2, H1, H2, L∞, C1, C2, · · · .

Par ailleurs, même pour des problèmes elliptiques simples, la solution n'est pas toujours régulière ... jusqu'au bord. Par exemple, considérons le problème

−∆u = 1 dans Ω = (0, 1)× (0, 1) (2.6) u = 0 sur ∂Ω (2.7)

Ce problème n'a pas de solution dans C2(Ω). En eet, si une telle solution existait, on aurait u(x, 0) = 0

pour tout x ∈ [0, 1] et donc ∂ 2u

∂x2 (0, 0) = 0. De même, on aurait u(0, y) = 0 pour tout y ∈ [0, 1] et donc

∂2u

∂y2 (0, 0) = 0. Ainsi, on obtiendrait ∆u(0, 0) = 0, ce qui contredirait (2.6) en faisant tendre (x, y) vers

(0, 0).

2.1.2 Principes du maximum

On donne à présent des résultats de principe du maximum pour des solutions classiques.

On suppose que aij ∈ C1(Ω) vérient (2.3) et c ∈ C(Ω) avec c ≥ 0. Soit u ∈ C2(Ω) ∩ C0(Ω) vériant Lu = f dans Ω.

1. u est solution faible de (2.1),(2.2) si u ∈ H10 (Ω) et vérie ∫

Ω A∇u · ∇v +

∫ Ω c uv =

∫ Ω fv, ∀v ∈ H10 (Ω).

2. On a supposé Ω régulier ; en général, on n'a pas la régularité u ∈ H2(Ω) si Ω présente des coins rentrants par exemple.

2.1. INTRODUCTION - PROPRIÉTÉS DES SOLUTIONS 13

0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

0.6

0.8

1

−150

−100

−50

0

50

100

150

Terme source initial f .

0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

0.6

0.8

1

−0.5

0

0.5

Solution de -∆u1 = f .

0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

0.6

0.8

1

−4

−3

−2

−1

0

1

2

3

4

x 10 −3

Solution de -∆u2 = u1.

Figure 2.1  Eet régularisant du Laplacien.

14 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

i) • On prend c ≡ 0. Si f ≤ 0 (resp. f ≥ 0) alors u atteint son maximum (resp. minimum) sur ∂Ω.

ii) • Si c ≡ 0 et f ≡ 0 alors inf ∂Ω u ≤ u ≤ sup

∂Ω u. Ceci montre qu'avec c ≡ 0 et f ≡ 0, si on choisit en

particulier u|∂Ω = 0 alors on obtient u ≡ 0.

iii) • Si f ≥ 0 et u|∂Ω ≥ 0 alors u ≥ 0 dans Ω.

iv) • (Principe de Hopf) Soit f ≤ 0 et soit x0 ∈ ∂Ω tel que u(x0) > u(x), ∀x ∈ Ω. On suppose que u est dérivable en x0. Si c ≡ 0 ou bien si u(x0) = 0, alors

∂u

∂n (x0) > 0,

où n désigne la normale extérieure à ∂Ω.

Démonstration (directe) de iii) dans le cas c > 0.

Soit x0 ∈ Ω tel que u(x0) = min Ω u, i.e. le point de Ω où le minimum de u est atteint.

- Si x0 ∈ ∂Ω alors u(x) ≥ u(x0) ≥ 0 pour tout x ∈ Ω et la démonstration est terminée.

- Si x0 ∈ Ω alors ∇u(x0) = 0 et ∂2u

∂x2i (x0) ≥ 0 pour i = 1, · · · , n. On va alors montrer que

∑ i,j

aij(x0) ∂2u

∂xixj (x0) ≥ 0. (2.8)

On peut tout d'abord supposé, sans perte de généralité, que la matrice A est symétrique (cf. note de bas de page de la Section 1.5, Chap. 1). La matrice A est donc diagonalisable par

D = CACT , (2.9)

avec D = diag(λ1, · · · , λn) où les λi sont les valeurs propres (positives) de A et C = C(x0) est la matrice orthogonale (C−1 = CT ) formée des vecteurs propres associés. On considère alors la fonction ũ(y) = u(x) avec y = Cx. La fonction ũ est dénie sur l'ouvert O = {y ∈ Rn | y = Cx, x ∈ Ω} et atteint son minimum au point intérieur y0 = Cx0 ∈ O( 1). Par conséquent, on a, pour j = 1, · · · , n,

∂2ũ

∂y2j (y0) ≥ 0, avec y0 = Cx0. (2.10)

On note (cij) et (dij) les coecients de C et D. En utilisant le fait que ∂yi ∂xj

= cij et la décomposition

(2.9), on montre que

Lu = ∑ i,j

aij ∂2u

∂xixj = ∑ i,j

dij ∂2ũ

∂yiyj = ∑ i

λi ∂2ũ

∂y2i . (2.11)

On conclut à (2.8) par (2.11), (2.10) et le fait que toutes les valeurs propres λi de A soient positives (condition d'ellipticité (2.3)).

On termine alors la démonstration en écrivant que −div (A∇u) (x0) + c(x0)u(x0) = f(x0) ; puisque

−div (A∇u) (x0) = − ∑

i,j aij(x0) ∂2u

∂xixj (x0) ≤ 0 et f(x0) ≥ 0, on en déduit que c(x0)u(x0) ≥ f(x0) ≥ 0

et donc que u(x0) ≥ 0 (on a supposé c > 0), ce qui montre que u ≥ 0 dans Ω. 

1. O = f−1(Ω) est un ouvert comme image réciproque de l'ouvert Ω par l'application continue f : y 7→ C−1y.

2.2. DIFFÉRENCES FINIES POUR LE CAS 1D 15

2.2 Diérences nies pour le cas 1D

On choisit pour simplier Ω = (0, 1). Le problème consiste à trouver une fonction u = u(x) qui vérie

(P )

{ Lu := − (au′)′ + c u = f dans (0, 1)

u(0) = u(1) = 0.

On se donne les fonctions f, c ∈ C([0, 1]) avec c ≥ 0 et a ∈ C1([0, 1]) avec a(x) ≥ α > 0 pour tout x ∈ [0, 1], de sorte que (P ) admet une unique solution u ∈ C2([0, 1]).

Principe des Diérences Finies (DF).

On discrétise le segment [0, 1] et on approche les dérivées aux points de discrétisation, par des opéra- teurs aux diérences. Plus précisément, on se donne un entier N à partir duquel on dénit le pas de discrétisation h = 1/(N + 1) ; on introduit les N + 2 points xi = i h pour i = 0, · · · , N + 1 qui forment alors une subdivision régulière de l'intervalle [0, 1]. On dénit Ωh le maillage intérieur à (0, 1) par

Ωh = {xi, 1 ≤ i ≤ N} (2.12)

et Ωh le maillage complet de [0, 1] par

Ωh = {xi, 0 ≤ i ≤ N + 1}. (2.13)

On désigne par Vh (resp. V h) l'ensemble des fonctions (discrètes) dénies sur le maillage Ωh (resp. Ωh), qu'on identie à RN (resp. RN+2) c'est-à-dire qu'un élément vh de Vh (resp. V h) est identié à un vecteur vh de RN (resp. RN+2) dont les composantes sont (vh(x1), · · · , vh(xN )) (resp. (vh(x0), · · · , vh(xN+1))). Les composantes de vh seront désormais notées vi = vh(xi).

On dénit l'opérateur aux diérences centrées Dh : V h/2 7→ Vh/2 par

Dhv(x) = v(x+ h/2)− v(x− h/2)

h , x ∈ Ωh (2.14)

pour une fonction v ∈ V h/2. Comme on le verra plus tard, Dhv est une approximation de la dérivée de v aux points du maillage. Par application itérée de l'opérateur aux diérences Dh, on obtient

Dh (aDhuh) (x) = 1

h2 [a(x+ h/2)uh(x+ h)− (a(x+ h/2) + a(x− h/2))uh(x) + a(x− h/2)uh(x− h)] ,

(2.15) qui est déni pour uh ∈ V h et x ∈ Ωh. On remarquera que seules interviennent les valeurs de la fonction uh aux points xi.

On considère alors le problème approché qui consiste à trouver une fonction uh ∈ V h telle que

(Ph)

{ Lhuh := −Dh (aDhuh) + c uh = f dans Ωh

uh(0) = uh(1) = 0.

Les valeurs des fonctions c et f aux pointx xi seront notées ci et fi.

Cas particulier a ≡ 1. Le problème approché (Ph) s'écrit dans ce cas, de la façon suivante. Trouver les ui tels que −(ui+1 − 2ui + ui−1)h2 + ciui = fi pour i = 1, · · · , Nu0 = uN+1 = 0. (2.16)

On peut écrire le système précédent sous forme matricielle. Les inconnues sont regroupées dans le vecteur uh = (u1, · · · , uN )> ∈ RN qui vérie

Ahuh = b, (2.17)

16 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

où la matrice Ah de taille N ×N est donnée par

Ah = 1

h2

 2 −1 0 −1 2 −1

. . . . . .

. . .

−1 2 −1 0 −1 2

+ Ch avec Ch = diag(c1, · · · , cN ) (2.18) et le second membre b = (f1, · · · , fN )> ∈ RN .

Cas général a 6≡ 1. Le problème approché (Ph) peut s'écrire de la façon suivante. Trouver les ui tels que{ − 1 h2 [ ai+1/2ui+1 −

( ai+1/2 + ai−1/2

) ui + ai−1/2ui−1

] + ciui = fi pour i = 1, · · · , N

u0 = uN+1 = 0, (2.19)

avec ai+1/2 = a(xi+h/2) et ai−1/2 = a(xi−h/2). La matrice Ah du système linéaire (2.17) correspondant est donnée dans le cas général par

Ah = 1

h2

 α1 β1 0 β1 α2 β2

. . . . . .

. . .

βN−2 αN−1 βN−1 0 βN−1 αN

+ Ch (2.20) avec αi = ai+1/2 + ai−1/2, βi = −ai+1/2 et Ch = diag(c1, · · · , cN ). Cette matrice est symétrique et dénie positive (voir Proposition 2.7), de sorte que le système linéaire (2.17) (et donc le problème approché) admet une unique solution uh.

2.2.1 Erreur de consistance

L'erreur de consistance correspond à l'erreur commise lorsqu'on remplace l'opérateur diérentiel L : v 7→ − (a v′)′+cv par l'opérateur aux diérences Lh : v 7→ −Dh (aDhv)+cv. L'erreur de consistance de (Ph) par rapport à (P ) est dénie par

Rhv(x) = Lv(x)− Lhv(x), (2.21)

pour v ∈ C2([0, 1]) et x ∈ Ωh (i.e. pour x = xi, i = 1, · · · , N). Le problème (Ph) est dit consistant par rapport à (P ) si pour chaque v, on a Rhv = Lv − Lhv → 0, quand h → 0. La convergence précédente doit être précisée. On s'intéresse ici à la norme discrète du maximum, dénie pour des fonctions v ∈ Vh par

‖v‖h,∞ = max 1≤i≤N

|v(xi)|. (2.22)

On remarquera que dans la norme discrète, on ne prend pas en compte les valeurs de la fonction aux points extrèmes x0 = 0 et xN+1 = 1. Mais, comme on le verra, ceci n'est pas génant pour évaluer la diérence entre les solutions exacte et approchée car toutes deux coïncident et sont nulles en 0 et 1 (Dirichlet homogène). Enn, on rappelle qu'on identie une fonction v ∈ Vh avec un vecteur v ∈ RN et on a ‖v‖h,∞ = ‖v‖∞ := max1≤i≤N |vi| où les vi sont les composantes de v (avec vi = v(xi)). Proposition 2.1 Le problème (Ph) est consistant par rapport à (P ) et on a

i) Pour a ≡ 1, si v ∈ C4([0, 1]) alors on a

‖Rhv‖h,∞ = ‖Lv − Lhv‖h,∞ ≤ 1

12 ‖v(4)‖∞ h

2. (2.23)

ii) Si a ∈ C3([0, 1]) et v ∈ C4([0, 1]) alors il existe une constante C = C(‖a‖C3 , ‖v‖C4) > 0 indépen- dante de h telle que

‖Rhv‖h,∞ = ‖Lv − Lhv‖h,∞ ≤ C h 2. (2.24)

2.2. DIFFÉRENCES FINIES POUR LE CAS 1D 17

Démonstration. La preuve est donnée dans le cas a ≡ 1 (traiter le cas général a ∈ C3 en exercice). Remarquons tout d'abord que

Dh (Dhv) = v(x+ h)− 2v(x) + v(x− h)

h2 .

Les développements de Taylor-MacLaurin de v ∈ C4([0, 1]) donnent

v(x+ h) = v(x) + hv′(x) + h2

2 v′′(x) +

h3

3! v(3)(x) +

h4

4! v(4)(x+ θ1h) avec 0 < θ1 < 1.

v(x− h) = v(x)− hv′(x) + h 2

2 v′′(x)− h

3

3! v(3)(x) +

h4

4! v(4)(x+ θ2h) avec 0 < θ2 < 1.

Des deux équations précédentes, on obtient

v(x+ h)− 2v(x) + v(x− h) = h2v′′(x) + h 4

4!

( v(4)(x+ θ1h) + v

(4)(x+ θ2h) )

d'où l'on déduit (2.23). 

2.2.2 Matrices monotones

Avant de s'intéresser à la notion de stabilité, on rappelle quelques résultats utiles par la suite, sur les matrices monotones et les M-matrices.

Dénition 2.1 Une matrice A ∈ RN×N est dite monotone si A est inversible et si A−1 ≥ 0 i.e. tous les coecients de A−1 sont positifs ou nuls.

La terminologie est en fait justiée par le résultat suivant.

Proposition 2.2 Une matrice A est monotone si et seulement si (Ax ≥ 0⇒ x ≥ 0).

Démonstration. - Supposons tout d'abord que A soit une matrice monotone et considérons x ∈ RN tel que Ax ≥ 0.

Montrons que x ≥ 0. La matrice A étant inversible, on écrit x = A−1Ax. Comme A−1 ≥ 0 et que Ax ≥ 0, on en déduit que x ≥ 0.

- Supposons maintenant que la propriété (Ax ≥ 0 ⇒ x ≥ 0) soit vraie. Montrons d'abord que A est inversible. Soit x ∈ RN tel que Ax = 0. On déduit de l'hypothèse que x ≥ 0. De plus, on a A(−x) = 0 d'où l'on déduit également que −x ≥ 0. On a ainsi établi que x = 0. Il reste à montrer que A−1 ≥ 0. Soit alors y ∈ RN tel que y ≥ 0 et on pose x = A−1y et par conséquent y = Ax ≥ 0. L'hypothèse implique que x ≥ 0 et donc que A−1y ≥ 0. Par des choix convenables de y (y = (1, 0, · · · , 0)>,y = (0, 1, 0, · · · , 0)>, · · · ,y = (0, · · · , 0, 1)>), on en déduit que A−1 ≥ 0. 

Dénition 2.2 Une matrice A ∈ RN×N est une M-matrice si elle est monotone et si aij ≤ 0 pour i 6= j.

Voici quelques critères pratiques de détermination de M-matrices.

Proposition 2.3 Soit une matrice A ∈ RN×N telle que i) aij ≤ 0 pour tous i 6= j. ii)

∑ 1≤j≤N

aij > 0 pour i = 1, · · · , N .

Alors A est une M-matrice.

Démonstration. On a clairement aii > 0, pour tout i. On introduit alors la matrice diagonale D = diag (aii) qui est

inversible avec D−1 = diag (1/aii) et par conséquent D −1 ≥ 0. On décompose alors A sous la forme

A = D(I−M) où I ∈ RN×N est la matrice identité etM = (mij) avec mij = − aij aii

pour i 6= j et mii = 0. Par l'hypothèse i), on a M ≥ 0. On va montrer que I −M est inversible puis monotone. Pour cela, on va utiliser un résultat bien connu sur le rayon spectral de matrice (cf. [2]).

18 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

Lemme 2.1 Soit M ∈ RN×N une matrice dont on note λi(M), i = 1, · · · , N les valeurs propres et ρ(M) = max

i |λi(M)| le rayon spectral( 1).

i) Pour toute norme matricielle ‖ · ‖ subordonnée( 2),on a ρ(M) ≤ ‖M‖.

ii) Si ρ(M) < 1 alors I −M est inversible et (I −M)−1 = ∞∑ k=0

Mk.

En utilisant le point i) de ce lemme, on a ρ(M) ≤ ‖M‖∞ = sup 1≤i≤N

N∑ j=1

|mij |. Or, d'après les hypothèses i)

et ii) de la Proposition, il vient

N∑ j=1

|mij | = ∑ j 6=i

∣∣∣∣−aijaii ∣∣∣∣ = − 1aii∑

j 6=i aij < 1,

ce qui implique d'après le point ii) du lemme, que la matrice I−M est inversible et (I −M)−1 = ∞∑ k=0

Mk.

Comme M ≥ 0, on a aussi Mk ≥ 0 et donc (I −M)−1 ≥ 0. Ainsi, A = D(I −M) est inversible et A−1 = (I −M)−1D−1 ≥ 0. Les coecients aij étant tous négatifs ou nuls par hypothèse, la matrice A est bien une M-matrice. 

On peut relacher la contrainte de positivité stricte de la somme des coecients et obtenir le résultat suivant.

Proposition 2.4 Soit une matrice A ∈ RN×N telle que i) aij ≤ 0 pour tous i 6= j. ii)

∑ 1≤j≤N

aij ≥ 0 pour i = 1, · · · , N .

iii) A est inversible.

Alors A est une M-matrice.

Démonstration. - Montrons d'abord que A + εI est monotone, quelque soit ε > 0. On pose A + εI =

( aεij

) i,j

avec

aεij = aij + εδij . Pour ε > 0, on a a ε ij ≤ 0 pour i 6= j et

∑ j

aεij = ∑ j

aij + ε ≥ ε > 0. Par conséquent,

d'après la proposition 2.3, la matrice A+ εI est monotone. - Montrons à présent que A est monotone par passage à la limite ε → 0. La matrice A est inversible

donc on peut écrire A+ εI = A ( I + εA−1

) . De plus, on a ρ(εA) = ερ(A) < 1 pour ε susamment petit

donc I + εA−1 est inversible et ( I + εA−1

)−1 =

∞∑ k=0

( −εA−1

)k . Par ailleurs, A+ εI est inversible avec

(A+ εI)−1 = ( I + εA−1

)−1 A−1 =

( ∞∑ k=0

(−ε)k(A−1)k ) A−1 =

( I +

∞∑ k=1

(−ε)k(A−1)k ) A−1.

Ainsi, (A+ εI)−1 −A−1 = ∞∑ k=1

(−ε)k(A−1)k A−1 et donc

∥∥∥(A+ εI)−1 −A−1∥∥∥ ≤ ( ∞∑ k=1

εk‖A−1‖k ) ‖A−1‖ = ε‖A

−1‖ 1− ε‖A−1‖

‖A−1‖ → 0, quand ε→ 0,

1. Si M est symétrique, on a ‖M‖2 := sup x 6=0

‖Mx‖2 ‖x‖2

= ρ(M).

2. ‖M‖ = sup x6=0

‖Mx‖ ‖x‖

2.2. DIFFÉRENCES FINIES POUR LE CAS 1D 19

pour n'importe quelle norme matricielle ‖·‖ subordonnée. Par conséquent chaque coecient de (A+ εI)−1 qui est positif ou nul, tend vers le coecient correspondant de A−1 qui est donc également positif ou nul. Ainsi A−1 ≥ 0 et A est une M-matrice. 

2.2.3 Stabilité

Soient u la solution de (P ) (avec Lu = f dans Ω) et uh la solution de (Ph) (avec Lhuh = f dans Ωh). On a Lu−Lhuh = 0 dans Ωh et donc 0 = Lu−Lhu+Lh(u−uh) = Rhu+Lh(u−uh) où Rhu est l'erreur de consistance (encore appelée résidu). On obtient ainsi

Rhu = Lh(uh − u) dans Ωh. (2.25)

On rappelle (cf. Proposition 2.1) que ‖Rhu‖h,∞ ≤ Ch2 si la solution u est susamment régulière. La question que l'on se pose à présent est de savoir si la diérence u − uh entre les solutions exacte et approchée est petite, lorsque le résidu Rhu est petit (avec h petit) ? C'est la notion de stabilité.

Plus précisement, on dit que le problème approché (Ph) est stable si pour toute fonction f , la solution correspondante uh de (Ph) avec f , est bornée par la donnée f . Les normes utilisées sont liées à la régularité des solutions c'est-à-dire aux espaces dans lesquels on cherche la solution et dans lesquels les fonctions sont données. La notion de stabilité dépend donc des normes utilisées et comme on le verra plus tard, un problème peut être stable pour une norme et ne pas l'être pour une autre norme.

On va montrer que le problème (Ph) est stable pour la norme discrète du maximum. Commençons d'abord par un principe du maximum discret.

Proposition 2.5 (Principe du maximum discret) Si f ≥ 0, alors la solution uh du problème (Ph) vérie uh ≥ 0 dans Ωh.

Démonstration. On donne une démonstration directe dans le cas où a ≡ 1 et c > 0. On retrouvera ce résultat plus tard, en montrant que la matrice du système linéaire est monotone. La solution uh vérie

−Dh(Dhuh)(x) + c(x)uh(x) = − uh(x− h)− 2uh(x) + uh(x+ h)

h2 + c(x)uh(x) = f(x),

pour x ∈ Ωh, c'est-à-dire aux points intérieurs x = xi, i = 1, ..., N de l'intervalle (0, 1). Soit xi0 le point de Ωh où le minimum de uh ∈ V h est atteint i.e. uh(xi0) ≤ uh(xi) pour tout i = 0, · · · , N + 1. Si i0 ∈ {1, · · · , N}, alors −Dh (Dhuh) (xi0) ≤ 0 et donc c(xi0)uh(xi0) ≥ f(xi0) ≥ 0, d'où uh(xi0) ≥ 0. Si i0 = 0 ou i0 = N + 1, alors uh(xi0) = 0. Dans tous les cas, on a 0 ≤ uh(xi), ∀i ∈ {0, · · · , N + 1}, d'où la conclusion. 

Comme nous allons le voir, ce résultat est en fait une conséquence directe d'une propriété plus générale de monotonie de la matrice Ah du problème approché. Les résultats sur les matrices monotones vont surtout nous servir à établir la stabilité du problème (Ph). On va examiner tout d'abord le cas homogène où a ≡ 1 puis le cas hétérogène où a est quelconque (a ∈ C3([0, 1])).

Stabilité dans le cas homogène a ≡ 1.

On suppose pour simplier que c ≡ 0. La stabilité s'établit à partir des propriétés de la matrice Ah ∈ RN×N du système linéaire (2.17). On rappelle (avec c ≡ 0) que

Ah = 1

h2

 2 −1 0 −1 2 −1

. . . . . .

. . .

−1 2 −1 0 −1 2

 . (2.26)

20 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

Proposition 2.6 La matrice Ah (avec c ≡ 0) vérie les propriétés suivantes. i) Ah est une M-matrice, symétrique dénie positive.

ii) (Stabilité) ‖A−1h ‖∞ = maxi ∑ j

|bij | ≤ 1

8 où les bij sont les coecients de la matrice A

−1 h .

Remarques :

• D'après la propriété i) précédente, la matrice Ah est monotone, ce qui permet de retrouver le principe du maximum discret établi à la Proposition 2.5.

• La proprièté ii) correspond bien à la notion de stabilité introduite au début du paragraphe. En eet, la solution uh du problème approché est identiée au vecteur uh qui vérie Ahuh = b. Par conséquent, on a ‖uh‖h,∞ = ‖uh‖∞ = ‖A

−1 h b‖∞ ≤ ‖A

−1 h ‖∞‖b‖∞ ≤

1 8‖f‖∞.

Démonstration de la Proposition. i) Il sut de vérier que (Ahx, x) = 1

h2

( x21 + x

2 N +

N∑ i=2

(xi − xi−1)2 )

pour tout x = (x1, · · · , xN )> et d'utiliser la Proposition 2.4 sur les M-matrices. ii) On a ‖A−1h ‖∞ = maxi

∑ j

|bij | = max i

∑ j

bij car A −1 h ≥ 0 et par conséquent ‖A

−1 h ‖∞ = ‖A

−1 h e‖∞

avec e = (1, 1, · · · , 1)> ∈ RN . On considère alors le problème suivant

Lu := −u′′ = 1 dans (0, 1) u(0) = u(1) = 0,

qui admet pour solution u(x) = 1

2 x(1 − x). On considère alors la solution uh du problème approché

associé au problème précédent, c'est-à-dire uh ∈ V h telle que Lhuh := −Dh(Dhuh) = 1 dans Ωh et uh(0) = uh(1) = 0. On a Rhu = Lhu − Lu = 0 car d'après la Proposition 2.1, il existe C > 0 telle que ‖Rhu‖h,∞ ≤ C‖u(4)‖∞ h2 et la fonction u est telle que u(4) ≡ 0 car u est un polynôme de degré 2. Ainsi, le résidu étant nul, on obtient d'après (2.25), Lhu = Lhuh = 1 dans Ωh et par conséquent uh(xi) = u(xi) pour tout i = 0, · · · , N+1, par unicité de la solution du problème approché. Finalement, on remarque que matriciellement, la solution approchée vérie le système linéaire Ahuh = e et par conséquent, on obtient

‖A−1h e‖∞ = ‖uh‖∞ = ‖uh‖h,∞ = ‖u‖h,∞ ≤ sup x∈[0,1]

|u(x)| = 1 8 . 

Stabilité dans le cas hétérogène.

On suppose seulement que a ∈ C3([0, 1]) mais toujours (pour simplier) que c ≡ 0. La matrice Ah est donnée cette fois par (2.20) (avec Ch = 0). En suivant la même démarche que dans le cas homogène précédent, on peut montrer les propriétés suivantes sur cette matrice.

Proposition 2.7

i) Ah est une M-matrice, symétrique dénie positive.

ii) (Stabilité) Si h est susamment petit, alors il existe une constante C > 0 indépendante de h telle que ‖A−1h ‖∞ ≤ C.

2.2.4 Convergence

La convergence est établie dans le cas homogène a ≡ 1 et avec c ≡ 0 (pour simplier).

Théorème 2.1 Soient f ∈ C2([0, 1]), a ≡ 1 et c ≡ 0. On note u la solution de (P ) et uh la solution de (Ph). Alors il existe une constante C > 0 indépendante de h telle que

‖u− uh‖h,∞ ≤ C‖f (2)‖∞ h

2.

2.2. DIFFÉRENCES FINIES POUR LE CAS 1D 21

Remarque. La norme discrète ‖ · ‖h,∞ ne porte que sur les points intérieurs xi, i = 1, · · · , N . Mais les fonctions u et uh coïncident et sont nulles en x0 = 0 et xN+1 = 1.

Démonstration. On a Lu = f dans Ω, Lhuh = f dans Ωh et Rhu = Lu − Lhu = Lh(uh − u) (cf. (2.25)). Matriciellement, on peut écrire Ahuh = b et Ah(u − uh) = Rhu avec u = (u(x1), · · · , u(xN ))>, Rhu = (Rhu(x1), · · · , Rhu(xN ))>. On en déduit que u − uh = A−1h Rhu. Puisque f ∈ C

2([0, 1]), alors u ∈ C4([0, 1]) et on a

‖u− uh‖h,∞ = ‖u− uh‖∞ ≤ ‖A −1 h ‖∞‖Rhu‖h,∞ ≤ C‖u

(4)‖∞h 2 = C‖f (2)‖∞h

2,

d'après les Propositions 2.1 et 2.6. 

2.2.5 Autres conditions limites

On peut choisir d'autres types de conditions limites pour le problème (P ).

Conditions de Dirichlet non-homogènes.

Plus généralement, des conditions de Dirichlet non-homogènes

u(0) = g0, u(1) = g1,

peuvent être imposées avec g0, g1 données. Il faut alors modier le second membre b du système linéaire Ahuh = b. Par exemple dans le cas a ≡ 1, la matrice Ah de taille N ×N est toujours donnée par (2.18) mais le second membre devient

b = ( f1 +

g0 h2 , f2, · · · , fN−1, fN +

g1 h2

)> . (2.27)

Conditions de Neumann.

On peut imposer aussi une condition de Neumann en choisissant par exemple

a(1)u′(1) = g1 (2.28)

avec g1 donnée, au lieu de u(1) = 0. En revanche, on conserve (pour simplier) la condition u(0) = 0 . Il y a maintenant N + 1 inconnues uh = (u1, · · · , uN , uN+1) ∈ RN+1. Pour les points intérieurs x1, · · · , xN , on utilise la même discrétisation que précédemment. Dans le cas où c ≡ 0 (pour simplier...), ceci fournit les N équations (cf. (2.19)) :

− 1 h2 [ ai+1/2ui+1 −

( ai+1/2 + ai−1/2

) ui + ai−1/2ui−1

] = fi pour i = 1, · · · , N. (2.29)

Mais, ayant une inconnue en plus (uN+1), il nous faut une équation supplémentaire. Cette dernière est obtenue par discrétisation de la condition limite (2.28). Il faut porter une attention particulière à cette discrétisation si l'on ne veut pas perdre l'ordre de convergence en h2. Une façon de faire est d'utiliser des diérences nies centrées qui font intervenir des points virtuels. Par développement de Taylor, on a

v(xN+1) = v(xN+1 − h/2) + v(xN+1 + h/2)

2 − h

2

16

( v′′(θ1) + v

′′(θ2) ) ,

pour θ1 ∈ (xN+1 − h/2, xN+1) et θ2 ∈ (xN+1, xN+1 + h/2). En prenant v = au′, on obtient

a(xN+1)u ′(xN+1) = g1 =

a(xN+1 − h/2)u′(xN+1 − h/2) + a(xN+1 + h/2)u′(xN+1 + h/2) 2

+O(h2).

Or u′(xN+1 ± h/2) = Dhu(xN+1 ± h/2) + O(h2). L'approximation considérée est alors la suivante (qui conserve l'ordre h2).

2g1 = a(xN+1 − h/2)Dhuh(xN+1 − h/2) + a(xN+1 + h/2)Dhuh(xN+1 + h/2)

22 CHAPITRE 2. EDP ELLIPTIQUES LINÉAIRES

soit

2g1 = aN+1/2 uN+1 − uN

h + aN+3/2

uN+2 − uN+1 h

. (2.30)

L'équation aux diérences nies Lhuh = f , s'écrit au noeud xN+1 :

− 1 h2 ( aN+3/2uN+2 − (aN+3/2 + aN+1/2)uN+1 + aN+1/2uN

) = fN+1.

On substitue alors la relation (2.30) donnant aN+3/2(uN+2 − uN+1) dans l'équation précédente pour obtenir

− 1 h2 ( 2g1h− aN+1/2(uN+1 − uN )− aN+1/2uN+1 + aN+1/2uN

) = fN+1,

ce qui donne

− 1 h2 ( aN+1/2uN − aN+1/2uN+1

) = fN+1

2 + g1 h . (2.31)

On peut à présent, donner le système linéaire Ahuh = b correspondant aux équations (2.29) et (2.31). La matrice Ah est de taille (N + 1)× (N + 1) et vaut

Ah = 1

h2

 α1 β1 0 β1 α2 β2

. . . . . .

. . .

βN−1 αN βN 0 βN −βN

 , (2.32)

avec αi = ai+1/2 + ai−1/2 et βi = −ai+1/2. Le second membre est donné par

b =

( f1, · · · , fN ,

fN+1 2

+ g1 h

)> ∈ RN+1. (2.33)

Remarque. On a introduit un point fantôme xN+3/2 = xN+1 + h/2 et le ux au ′ en ce point est calculé

par extrapolation linéaire du ux aux points xN+1/2 et xN+1. Cette technique parfois appelée mirror imaging permet de conserver un ordre global en h2, ce qui n'est pas le cas si on prend par exemple (ce qu'on aurait tendance à faire...) des diérences nies non-centrées pour discrétiser les ux au bord. Ainsi l'approximation

a(xN+1)u ′(xN+1) = g1 ' a(xN+1)D−h u(xN+1)

avec l'opérateur décentré D−h v(x) = v(x)− v(x− h)

h , ne fournit qu'un ordre de convergence en h car

v′(x) = D−h v(x) +O(h).

2.3 Diérences nies pour le cas 2D

Pour un ouvert borné Ω de R2, on considère le problème suivant : trouver une fonction u = u(x, y) qui vérie

(P )

{ Lu := −div (A(x)∇u) + c u = f dans Ω

u|∂Ω = 0 sur ∂Ω.

On suppose que la fonction c = c(x, y) est positive ou nulle et que l'équation est elliptique c'est-à-dire qu'il existe α > 0 tel que (A(x)ξ, ξ)Rn ≥ α|ξ|2, ∀ξ ∈ Rn, ∀x ∈ Ω. On suppose enn que (P ) admet une unique solution u ∈ C2(Ω) ∩ C(Ω).

2.3. DIFFÉRENCES FINIES POUR LE CAS 2D 23

2.3.1 Un schéma à 5 points pour le Laplacien

Commençons par le cas du Laplacien, c'est-à-dire que l'on prend

L = −∆ + c Id (2.34)

et plaçons nous pour commencer dans un domaine rectangulaire en choisissant

Ω = (0, a)× (0, b). (2.35)

Soient deux entiers N et M qui permettent de dénir les paramètres de discrétisation hx = a/(N + 1) et hy = b/(M + 1). On notera h = max(hx, hy). Enn, on introduit les points Pij = (xi, yj) de Ω avec

xi = i hx, i = 0, · · · , N + 1 (2.36) yj = j hy, j = 0, · · · ,M + 1. (2.37)

Ces points appelés aussi noeuds dénissent un maillage du domaine Ω. Plus précisement, on dénit Ωh le maillage intérieur à Ω par

Ωh = {Pij = (xi, yj), i = 1, · · · , N, j = 1, · · · ,M} (2.38)

et Ωh le maillage complet de Ω par

Ωh = {Pij = (xi, yj), i = 0, · · · , N + 1, j = 0, · · · ,M + 1}. (2.39)

On désignera par Γh les points du maillage complet Ωh qui sont sur le bord ∂Ω, c'est-à-dire

Γh = {Pij = (xi, yj), avec i = 0 ou N + 1, ou j = 0 ou M + 1}. (2.40)

Comme dans le cas monodimensionnel, on désignera par Vh (resp. V h) l'ensemble des fonctions (discrètes) dénies sur le maillage Ωh (resp. Ωh). On identiera un élément vh de Vh (resp. V h) à un vecteur vh de RN×M (resp. R(N+2)×(M+2)).

Le problème approché s'écrit alors : trouver uh ∈ V h telle que

(Ph)

{ −∆huh + c uh = f dans Ωh

uh = 0 sur Γh.

L'opérateur discret ∆h est donné par un schéma à 5 points. Pour un point intérieur P du maillage Ωh et une fonction v ∈ V h, on a

∆hv(P ) = 1

h2x [v(W )− 2v(P ) + v(E)] + 1

h2y [v(S)− 2v(P ) + v(N)] , (2.41)

où les points E,W,N, S sont les voisins de P , c'est-à-dire les 4 points du maillage Ωh les plus proches de P comme indiqué sur la gure 2.2.

s

s

s

ss

N

E

S

W

hy

hx hx

hy

P

Figure 2.2  Schéma à 5 points pour le Laplacien.

commentaires (0)

Aucun commentaire n'a été pas fait

Écrire ton premier commentaire

Ceci c'est un aperçu avant impression

3 shown on 46 pages

Télécharger le document