4.3 M´ ethodes it´ eratives stationnaires et instationnaires
147
on se donne x
(0)
∈ R
n et on pose r
(0) = b − Ax
(0) , z
(0) = P
−1 r
(0) et p
(0) =
z
(0) , la k-i` eme it´ eration s’´ ecrit alors
α k =
z
(k) T r
(k)
p (k) T Ap (k)
,
x
(k+1) = x
(k) + α k p
(k) ,
r
(k+1) = r
(k)
− α k Ap
(k) ,
Pz
(k+1) = r
(k+1) ,
β k =
z
(k+1) T r
(k+1)
z (k) T r (k)
,
p
(k+1) = z
(k+1) + β k p
(k) .
Le coˆ ut du calcul est plus ´ elev´ e que pour GC puisqu’on doit r´ esoudre ` a chaque
it´ eration le syst` eme lin´ eaire Pz
(k+1) = r
(k+1) . Pour ce syst` eme, on peut utiliser les pr´ econditionneurs sym´ etriques vus `
a la Section 4.3.2. L’estimation de
l’erreur est la mˆ eme que pour la m´ ethode non pr´ econditionn´ ee, ` a condition de
remplacer A par P
−1 A.
Nous donnons dans le Programme 20 une impl´ ementation MATLAB du
gradient conjugu´ e pr´ econditionn´ e. Pour une description des param` etres en
entr´ ee et en sortie, voir le Programme 19.
Programme 20 - conjgrad : M´ ethode du gradient conjugu´ e pr´ econditionn´ e
function [x,relres,iter,flag]=conjgrad(A,b,x,P,nmax,tol)
%CONJGRAD M´ ethode du gradient conjugu´ e
% [X,RELRES,ITER,FLAG]=CONJGRAD(A,B,X0,P,NMAX,TOL) tente
% de r´ esoudre le syst` eme A*X=B avec la m´ ethode du gradient conjugu´ e. TOL est
% la tol´ erance de la m´ ethode. NMAX est le nombre maximum d’it´ erations.
% X0 est la donn´ ee initiale. P est un pr´ econditionneur. RELRES est le residu
% relatif. Si FLAG vaut 1, alors RELRES > TOL. ITER est l’it´ eration
% `
a laquelle la solution X a ´ et´ e calcul´ ee.
flag=0; iter=0; bnrm2=norm(b);
if bnrm2==0, bnrm2=1; end
r=b-A*x; relres=norm(r)/bnrm2;
if relres for iter = 1:nmax
z=P\r; rho=r’*z;
if iter>1
beta=rho/rho1;
p=z+beta*p;
else
p=z;