9.7 Chapitre 7
339
Prouvons à présent que τn+1 est d’ordre deux en h. Supposons que y ∈
C
3 ([t0, T [). Pour alléger les notations, on pose yn = y(tn) pour n ≥ 0. On a
τn+1 =
yn+1 − yn
h
−
1
2
[f (tn, yn) + f (tn+1, yn + hf (tn, yn))]
=
yn+1 − yn
h
−
1
2
y
(tn) −
1
2
f (tn+1, yn + hy
(tn)).
D’après l’expression de l’erreur de la formule du trapèze (4.20), il existe ξn ∈
]tn, tn+1[ tel que
yn+1 − yn =
t n+1
tn
y
(t)dt =
h
2
y
(tn) + y
(tn+1)
−
h
3
12
y
(ξn),
donc
τn+1 =
1
2
y
(tn+1) − f (tn+1, yn + hy
(tn)) −
h
2
6
y
(ξn)
=
1
2
f (tn+1, yn+1) − f (tn+1, yn + hy
(tn)) −
h
2
6
y
(ξn)
.
Ensuite, on utilise que la fonction f est lipschitzienne par rapport à sa
deuxième variable (voir Proposition 7.1), donc
|τn+1| ≤
L
2
|yn+1 − yn − hy
(tn)| +
h
2
12
|y
(ξn)|.
Enfin, on applique la formule de Taylor
yn+1 = yn + hy
(tn) +
h
2
2
y
(ηn),
η n ∈]tn, tn+1[,
et on obtient
|τn+1| ≤
L
4
h
2 |y
(ηn)| +
h
2
12
|y
(ξn)| ≤ Ch
2 .
La méthode de Heun est implémentée dans le Programme 9.2. En utilisant
ce programme, on peut vérifier l’ordre de convergence comme dans la Solution
7.1. Avec les instructions suivantes, on voit que la méthode de Heun est d’ordre
deux par rapport à h :
p=log(abs(e(1:end-1)./e(2:end)))/log(2); p(1:2:end)
ans =
1.7642
1.9398
1.9851
1.9963
1.9991
Programme 9.2. rk2 : méthode de Heun (ou RK2)
function [tt , u ]= rk2( odefun , tspan , y0 ,Nh , varargin )
tt= linspace ( tspan (1) ,tspan (2) , Nh +1);
h =( tspan (2) - tspan (1))/ Nh ; hh= h *0.5;
u = y0;
for t = tt (1: end -1)
y = u ( end ,:);
k1= feval ( odefun ,t ,y , varargin {:});
t1 = t + h; y = y + h*k1;
k2= feval ( odefun , t1 ,y , varargin {:});
u = [ u; u ( end ,:) + hh *( k1+ k2 )];
end