7.4 Méthode de Crank-Nicolson
217
La dernière égalité, qui découle de (7.18), fait apparaître, à un facteur
1/h près, l’erreur de la formule du trapèze (4.19). En supposant y ∈ C
3
et en utilisant (4.20), on en déduit que
τ n (h) = −
h
2
12
y
(ξ n ) pour un certain ξ n ∈]t n−1 , t n [.
(7.19)
La méthode de Crank-Nicolson est donc consistante à l’ordre 2, i.e. son
erreur de troncature locale tend vers 0 comme h
2 . En procédant comme
pour la méthode d’Euler explicite, on peut montrer que la méthode de
Crank-Nicolson converge à l’ordre 2 en h.
La méthode de Crank-Nicolson est implémentée dans le Programme
7.3. Les paramètres d’entrée et de sortie sont les mêmes que pour les
méthodes d’Euler.
Programme 7.3. cranknic : méthode de Crank-Nicolson
function [t ,u ]= cranknic ( odefun , tspan , y0 ,Nh , varargin )
% CRANKNIC Résout une équation d i f f ére ntie lle avec la
% méthode de Crank - Nicolson .
% [T ,Y ]= CRANKNIC ( ODEFUN , TSPAN ,Y0 , NH) avec
% TSPAN =[T0 , TF]
% intègre le système d ’ é q u a tions d i f f ér enti ell es
% y ’=f (t ,y ) du temps T0 au temps TF avec la c o n ditio n
% initiale Y0 en u t i l isant la méthode de
% Crank - Nicolson sur une grille de NH i n t e rv alles
% é q u i d istr ibu és. La fonction ODEFUN (T ,Y ) doit
% r e t o urner un vecteur c o r r e spon dant à f (t ,y )
% de même d i m e nsio n que Y .
% Chaque ligne de la solution Y c o r r espo nd
% à un temps du vecteur colonne T.
% [T ,Y ] = CRANKNIC ( ODEFUN , TSPAN , Y0 ,NH , P1 ,P2 ,...)
% passe les p a r amè tres s u p p l éme nta ires P1 , P2 ,.. à
% la fonction ODEFUN de la manière suivante :
% ODEFUN (T ,Y ,P1 , P2 ...).
tt= linspace ( tspan (1) ,tspan (2) , Nh +1);
y = y0 (:); % crée toujours un vecteur colonne
u =y . ’;
global glob_h glob_t glob_y g l o b_ odefu n;
glob_h =( tspan (2) -tspan (1))/ Nh;
glob_y = y;
g l o b_ odef un= odefun ;
if ( exist ( ’ O C T A VE _VE RSIO N’) )
o_ver = O C T A VE _VER SION;
version = str2num ([ o_ver (1) , o_ver (3) , o_ver (5)]);
end
if( ~ exist ( ’ O C T A VE _VE RSIO N’) | version >= 320 )
options = optimset ;
options . Display = ’ off ’;
options . TolFun =1.e -12;
options . M a x F un Evals =10000;
end
for glob_t = tt (2: end)
if ( exist ( ’ O C T A VE _VE RSIO N’) & version < 320 )
w = fsolve ( ’ c r a n kni cfun’, glob_y );