10.7 M´ ethodes pr´ edicteur-correcteur
405
Programme 81 - predcor : Sch´ ema pr´ edicteur-correcteur
function [t,u]=predcor(a,b,at,bt,h,f,t0,u0,tf,pece,m)
%PREDCOR Sch´ ema pr´ edicteur-correcteur.
% [T,U]=PREDCOR(A,B,AT,BT,H,FUN,T0,U0,TF,PECE,M) r´ esout le probl` eme
% de Cauchy Y’=FUN(T,Y) pour T dans ]T0,TF[ en utilisant un pr´ edicteur
% avec des coefficients AT et BT, et un correcteur avec A et B. H est
% le pas de temps. Si PECE=1, on utilise la m´ ethode P(EC)ˆmE, sinon
% la m´ ethode P(EC)ˆm
p = max(length(a),length(b)-1);
pt = max(length(at),length(bt));
q = max(p,pt); if length(u0) t = [t0:h:t0+(q-1)*h]; u = u0; y = u0; fe = eval(f);
k = q;
for t = t0+q*h:h:tf
ut=sum(at.*u(k:-1:k-pt+1))+h*sum(bt.*fe(k:-1:k-pt+1));
y=ut; foy=eval(f);
uv=sum(a.*u(k:-1:k-p+1))+h*sum(b(2:p+1).*fe(k:-1:k-p+1));
k = k+1;
for j=1:m
fy=foy; up=uv+h*b(1)*fy; y=up; foy=eval(f);
end
if pece==1
fe=[fe,foy];
else
fe=[fe,fy];
end
u=[u,up];
end
t=[t0:h:tf];
return
Exemple 10.9 V´ erifions la performance du sch´ ema P (EC)
m E sur le probl` eme
de Cauchy y
(t) = e
−y(t) pour t ∈ [0, 1] avec y(0) = 1. La solution exacte est
y(t) = log(1 + t). Dans les tests num´ eriques, le correcteur est le sch´ ema d’AdamsMoulton du 3` eme ordre (AM3), et le pr´ edicteur est soit le sch´ ema d’Euler explicite
(AB1) soit le sch´ ema d’Adams-Bashforth du second ordre (AB2). La Figure 10.8
montre que la paire AB2-AM3 a une convergence du troisi` eme ordre (m = 1),
et AB1-AM3 est du premier ordre (m = 1). En prenant m = 2, on r´ ecup` ere la
convergence d’ordre 3 du correcteur pour la paire AB1-AM3.
•
Comme pour la stabilit´ e absolue, le polynˆ ome caract´ eristique de P (EC)
m
s’´ ecrit
Π P (EC) m (r) = b −1 r
p ( ρ(r) − hλ σ(r)) +
H
m (1 − H)
1 − H m (˜ ρ(r) σ(r) −
ρ(r)˜ σ(r))