EDO 2 Méthodes à un pas
2 Méthodes à un pas
Constituent l’algorithme de base qui permet d’estimer la valeur de la solution à
l’instant ti+1 = ti + h, connaissant seulement ui , celle à ti .
La valeur à estimer peut être approchée par un développement limité de Taylor :
dy h2 d2 y
y(ti + h) = y(ti ) + h (ti ) + 2
(ti ) + · · · (1)
dt 2 dt
Ordre n de la méthode = plus grande puissance de h prise en compte dans
l’approximation.
— Somme des termes négligés = erreur de troncature locale ∝ hn+1
déterministe, augmente si le pas h augmente et si l’ordre de la méthode diminue
— Précision finie des opérations sur les réels ⇒ erreur d’arrondi aléatoire
augmente lorsque les calculs se compliquent, en particulier si le pas h diminue.
Indépendamment du coût (en temps de calcul) des opérations, et des cas où la
fonction est tabulée, ne pas croire que diminuer le pas améliore toujours la
qualité du résultat : un compromis doit être trouvé entre ces deux types d’erreurs.
MNCS 12 2019-2020
EDO 2 Méthodes à un pas 2.1 Méthodes du premier ordre
2.1 Méthodes du premier ordre
2.1.1 Méthode d’Euler progressive (explicite)
Méthode du premier ordre d’intérêt pédagogique, à éviter en pratique
ui+1 = ui + hf (ti , ui ) (2)
Exemple : stabilité
dy y
=− ⇒ solution analytique y = y0 e−t/τ ⇒ yn = y0 (e−h/τ )n
dt τ
h n
ui+1 = ui − ui ⇒ solution numérique un = y0 (1 − h/τ )
τ
Si τ> 0, la solution exacte vérifie y(∞) = 0,
Mais pour l’approximation, un → 0 ⇐⇒ |1 − h/τ | < 1 ⇐⇒ 0 < h < 2τ .
Condition de stabilité : h < 2τ (pas h petit)
Mais, si h > τ , alors (1 − h/τ ) < 0 : alternance de signe de la solution un .
MNCS 13 2019-2020
EDO 2 Méthodes à un pas 2.1 Méthodes du premier ordre
y Méthode d’Euler
Méthode explicite qui ne
nécessite qu’une
seule évaluation
ui+1
de la fonction second
k1 h membre f par pas :
k1 k1 = f (ti , ui )
ui facilement instable
h
ui+1 − ui
= f (ti , ui )
ti ti+1 t h
voir dérivée avant
MNCS 14 2019-2020
EDO 2 Méthodes à un pas 2.1 Méthodes du premier ordre
2.1.2 Méthode d’Euler rétrograde (implicite)
ui+1 = ui + hf (ti+1 , ui+1 ) (3)
Méthode implicite : résolution itérative , plus difficile à mettre en œuvre, sauf si la
forme de f (t, u) permet le calcul analytique de ui+1 à partir de l’équation (3).
Avantage : meilleure stabilité que la méthode progressive explicite.
Exemple : stabilité
dy y
=− ⇒ solution analytique y = y0 e−t/τ ⇒ yn = y0 (e−h/τ )n
dt τ
h ui
ui+1 = ui − ui+1 ⇒ solution numérique ui+1 =
τ 1 + h/τ
y0
un = n
(1 + h/τ )
Si τ > 0, y(∞) = 0, et aussi un → 0 ∀τ > 0, ∀h > 0 solution stable
MNCS 15 2019-2020
EDO 2 Méthodes à un pas 2.1 Méthodes du premier ordre
Mise en œuvre de la méthode d’Euler rétrograde : résolution de l’équation
implicite par itération
ui+1 = ui + hf (ti+1 , ui+1 )
Itérer l’application g pour rechercher son point fixe
v20 = g(v2 ) = ui + hf (t2 , v2 )
Ce point fixe est la solution de l’équation implicite.
— Utilise plusieurs évaluations du second membre, sans calcul de ses dérivées.
— Très peu d’itérations nécessaires
Initialisation par le prédicteur avec Euler progressif
t2 = ti + h
k1 = f (ti , ui )
v2 = ui + hk1
MNCS 16 2019-2020
EDO 2 Méthodes à un pas 2.1 Méthodes du premier ordre
Boucle pour recherche du point fixe de g(v2 ) = v20 = ui + hf (t2 , v2 )
k2 = f (t2 , v2 )
v20 = ui + hk2
δv2 = v20 − v2
arrêt si |δv2 |2 6 α2 |v2 |2
v2 = v20
∂f
La fonction g est contractante si g 0 (v2 ) = h ∂v 61
2
ce qui est vérifié si le pas est assez faible.
∂f
Rappel : la condition de Lipschitz est ∂v
2
6K
Critère d’arrêt : choisir α faible, mais α > ε.
MNCS 17 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
2.2 Méthodes du deuxième ordre
Première idée : augmenter le nombre de termes du développement de Taylor :
rarement utilisé, car nécessite l’évaluation des dérivées partielles de f .
dy d2 y ∂f ∂f dy ∂f ∂f
= f (t, y(t)) ⇒ = + = +f (6)
dt dt2 ∂t ∂y dt ∂t ∂y
Préférer utiliser plusieurs évaluations du second membre f en des points adaptés.
Centrer l’évaluation de la dérivée au point milieu tm = (ti + ti+1 )/2.
h dy 1 h2 d2 y 3
y(ti + h) = y(tm ) + (tm ) + 2
(t m ) + O(h ) (7a)
2 dt 2 4 dt
h dy 1 h2 d2 y 3
y(ti ) = y(tm ) − (tm ) + (t m ) + O(h ) (7b)
2 dt 2 4 dt2
Par différence, (approximation locale parabolique, voir aussi dérivée centrée à 2 termes)
dy
y(ti + h) − y(ti ) = h (tm ) + O(h3 )
dt
MNCS 18 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
2.2.1 Méthode du point milieu
Nécessite l’évaluation du second membre f en 2 points :
en (ti , ui ) et au milieu (ti+1/2 = ti + h/2, ui+1/2 ) d’un pas (hors grille).
h h
ui+1 = ui + hf ti + , ui + f (ti , ui )
2 2
k1 = f (ti , ui ) (8a)
h h
(ui+1/2 calculé via Euler) k2 = f (ti + , u i + k1 ) (8b)
2 2
ui+1 = ui + hk2 (8c)
MNCS 19 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
ui+1 Méthode du
k2 point milieu
k2 h Méthode explicite
k2 qui nécessite deux
évaluations du second
k1 h/2
k1 membre par pas dont
ui
une hors grille.
ti ti + h/2 ti+1 t
MNCS 20 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
2.2.2 Méthode d’Euler modifiée
En appliquant 7a et 7b à la dérivée et en faisant la somme, on peut remplacer la
dérivée au milieu par la moyenne des dérivées aux extrémités de l’intervalle
(voir méthode de quadrature dite des trapèzes) :
dy dy dy
(ti ) + (ti+1 ) = 2 (tm ) + O(h2 )
dt dt dt
D’où une approximation n’utilisant pas la valeur de f au point milieu tm :
h
ui+1 = ui + [f (ti , ui ) + f (ti+1 , ui+1 )]
2
De nouveau, méthode a priori implicite, plus stable, mais plus lourde.
⇒ Contournement du problème en utilisant l’approximation d’Euler explicite (voir 2)
pour évaluer ui+1 intervenant dans f .
h
ui+1 = ui + [f (ti , ui ) + f (ti+1 , ui + hf (ti , ui ))]
2
MNCS 21 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
Bilan : méthode de type prédicteur-correcteur équivalent à
— un demi-pas avec la pente initiale k1
— et un demi-pas avec la pente k2 du point prédit par Euler progressif.
k1 = f (ti , ui ) (9a)
k2 = f (ti+1 , ui + k1 h) (9b)
h
ui+1 = ui + [k1 + k2 ] (9c)
2
Remarques
— deuxième ordre comme point milieu mais sans évaluation hors grille
— la résolution de l’équation implicite peut se faire en itérant la correction jusqu’à
ce qu’elle devienne négligeable.
MNCS 22 2019-2020
EDO 2 Méthodes à un pas 2.2 Méthodes du deuxième ordre
corrigé Méthode d’Euler
ui+1 modifiée
k2 h/2 k2 k2
Méthode explicite qui
prédit
nécessite deux éva-
k1 h/2 luations de la fonction
k1
ui par pas en des points
de la grille.
ti ti + h/2 ti+1 t
MNCS 23 2019-2020
EDO 2 Méthodes à un pas 2.3 Méthodes de Runge Kutta
2.2.3 Méthode de Heun
k1 = f (ti , ui ) (10a)
2 2
k2 = f (ti + h, ui + k1 h) (10b)
3 3
h
ui+1 = ui + [k1 + 3k2 ] (10c)
4
2.3 Méthodes de Runge Kutta
Plus généralement, avec r évaluations de f , on peut atteindre une méthode d’ordre
r si r 6 4. Pour atteindre l’ordre 5, six évaluations sont nécessaires.
=⇒ la méthode de Runge Kutta d’ordre 4 est très utilisée.
Les méthodes de Runge-Kutta sont stables.
MNCS 24 2019-2020
EDO 2 Méthodes à un pas 2.3 Méthodes de Runge Kutta
2.3.1 Méthode de Runge Kutta d’ordre 3
k1 = f (ti , ui ) (11a)
h h
k2 = f (ti + , ui + k1 ) (11b)
2 2
k3 = f (ti + h, ui + (2k2 − k1 )h) (11c)
h
ui+1 = ui + (k1 + 4k2 + k3 ) (11d)
6
MNCS 25 2019-2020
EDO 2 Méthodes à un pas 2.3 Méthodes de Runge Kutta
2.3.2 Méthode de Runge Kutta d’ordre 4
k1 = f (ti , ui ) (12a)
h h
k2 = f (ti + , ui + k1 ) (12b)
2 2
h h
k3 = f (ti + , ui + k2 ) (12c)
2 2
k4 = f (ti + h, ui + k3 h) (12d)
h
ui+1 = ui + (k1 + 2k2 + 2k3 + k4 ) (12e)
6
MNCS 26 2019-2020
EDO 2 Méthodes à un pas 2.4 Erreur absolue en fonction du pas et de l’ordre
2.4 Erreur absolue en fonction du pas et de l’ordre
nombre de pas = L/h =⇒ erreur globale ∼ erreur locale × L/h
TABLE 1 – Erreur de troncature seule
Méthode ordre erreur locale erreur globale
Euler explicite 1 ∝ h2 ∝h
Point milieu – Euler modifiée 2 ∝ h3 ∝ h2
Runge-Kutta 3 3 ∝ h4 ∝ h3
Runge-Kutta 4 4 ∝ h5 ∝ h4
Erreur d’arrondi locale indépendante de h ⇒ erreur d’arrondi globale ∝ 1/h
MNCS 27 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
dy
2.5 Exemple de l’équation logistique = y(1 − y/2)
dt
sol_exacte_logistique
2.0 [ 0.0e+00 , 2.0e+01] h= 2.00e-02 y0=1.00e-01 n= 1001
1.5
logistique
1.0
0.5
0.0
0 5 10 15 20
temps
F IGURE 1 – Solution analytique de l’équation logistique pour t0 = 0, y(t0 ) = 0.1,
a = 1 et k = 2.
MNCS 28 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
fichier [Link] ordre 1 erreur max 0.830E−02 position de l’err max 0.296E+01
0.001
0.000
−0.001
Erreur = Estimé − Exact −0.002
−0.003
−0.004
−0.005
−0.006
−0.007
−0.008
−0.009
0 2 4 6 8 10 12 14 16 18 20
estimé − exact
temps
F IGURE 2 – Erreur dans l’intégration de l’équation logistique avec la méthode d’Euler
pour h = 0, 02. L’allure régulière montre que l’erreur de troncature domine.
Erreur de troncature locale liée à la courbure de la solution
MNCS 29 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
fichier [Link] ordre 4 erreur max 0.954E−06 position de l’err max 0.190E+01
1e−06
8e−07
6e−07
Erreur = Exact −Estimé
4e−07
2e−07
0e+00
−2e−07
−4e−07
−6e−07
−8e−07
−1e−06
0 2 4 6 8 10 12 14 16 18 20
exact − estimé
temps
F IGURE 3 – Erreur dans l’intégration de l’équation logistique avec Runge Kutta
d’ordre 4 pour h = 0, 02. L’allure bruitée est caractéristique de l’erreur d’arrondi
et on retrouve les niveaux de quantification des réels sur 32 bits.
MNCS 30 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
2.5.1 Exemple d’erreur totale maximale en simple précision (32 bits)
Valeur absolue de l’erreur absolue pour y’= y(1−y/2) sur [0,20] y(0)=0.1 simple précision
100
euler sp
milieu sp
rk3 sp
10
−1 rk4 sp pente +1
10−2
erreur maximale
pente +2
−3
10
pente −1
−4
10 pente +3
10−5
pente +4
−6
10
10−7 −4
10 10−3 10−2 10−1 100
pas
MNCS 31 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
2.5.2 Exemple d’erreur totale maximale en double précision (64 bits)
Erreur absolue y’= y(1−y/2) sur [0,20] y(0)=0.1 double précision
100
10−2 pente +1
10−4
pente +2
10−6
erreur maximale
10−8
pente +3
10−10
10−12 pente +4
−14 euler dp
10
milieu dp
rk3 dp
rk4 dp
10−16 −4
10 10−3 10−2 10−1 100
pas
MNCS 32 2019-2020
EDO 2 Méthodes à un pas 2.5 Exemple de l’équation logistique
2.5.3 Comparaison des erreurs maximales simple/double précision
Erreur absolue y’= y(1−y/2) sur [0,20] y(0)=0.1 simple/double précision
100
10−2
10−4
10−6
erreur maximale
10−8
10−10
euler sp
10−12 milieu sp
rk3 sp
rk4 sp
−14 euler dp
10
milieu dp
rk3 dp
rk4 dp
10−16 −4
10 10−3 10−2 10−1 100
pas
MNCS 33 2019-2020
