156
4 M´ ethodes it´ eratives pour la r´ esolution des syst` emes lin´ eaires
Programme 23 - gmres : M´ ethode GMRES pour la r´ esolution des syst` emes
lin´ eaires
function [x,iter]=gmres(A,b,x0,m,tol)
%GMRES M´ ethode GMRES.
% [X,ITER]=GMRES(A,B,X0,M,TOL) tente de r´ esoudre le syst` eme
% A*X=B avec la m´ ethode GMRES. TOL est la tol´ erance de la m´ ethode.
% M est la taille maximale de l’espace de Krylov. X0 est la donn´ ee
% initiale. ITER est l’it´ eration ` a laquelle la solution X a ´ et´ e calcul´ ee.
r0=b-A*x0; nr0=norm(r0,2);
if nr0 ˜= 0
v1=r0/nr0; V=[v1]; H=[]; iter=0; residual=1;
while iter <= m-1 & residual > tol,
[iter,V,H] = GSarnoldi(A,m,iter,V,H);
[nr,nc]=size(H); y=(H’*H) \ (H’*nr0*[1;zeros(nr-1,1)]);
x=x0+V(:,1:nc)*y; residual = norm(b-A*x,2)/nr0;
end
else
x=x0;
end
Pour am´ eliorer l’efficacit´ e de l’impl´ ementation de GMRES, il est n´ ecessaire
de d´ efinir un crit` ere d’arrˆ et qui ne requiert pas le calcul explicite du r´ esidu
` a chaque pas. Ceci est possible si on r´ esout de fa¸ con appropri´ ee le syst` eme
associ´ e ` a la matrice de Hessenberg
H k .
En pratique,
H k est transform´ e en une matrice triangulaire sup´ erieure
R k ∈ R
(k+1)×k avec r k+1,k = 0 telle que Q
T
k R k =
H k , o` u Q k est le r´ esultat du
produit de k rotations de Givens (voir Section 5.6.3). On peut alors montrer
que, Q k ´ etant orthogonale, minimiser
(0)
2 e 1 −
H k z
(k)
2 est ´ equivalent `
a
minimiser f k − R k z
(k)
2 , avec f k = Q k r
(0)
2 e 1 . On peut aussi montrer que
la valeur absolue de la k + 1-i` eme composante de f k est ´ egale ` a la norme
euclidienne du r´ esidu `
a l’it´ eration k.
Tout comme la m´ ethode d’Arnoldi, GMRES est coˆ uteuse en calcul et en
m´ emoire ` a moins que la convergence ne survienne qu’apr` es peu d’it´ erations.
Pour cette raison, on dispose `
a nouveau de deux variantes de l’algorithme :
la premi` ere, GMRES(m), bas´ ee sur un red´ emarrage apr` es m it´ erations, la
seconde, Quasi-GMRES ou QGMRES, sur l’arrˆ et du proc´ ed´ e d’orthogonalisation d’Arnoldi. Dans les deux cas, on perd la propri´ et´ e de GMRES d’obtenir
la solution exacte en un nombre fini d’it´ erations.
Remarque 4.4 (m´ ethodes de projection) Les it´ erations de Krylov peuvent ˆ etre vues comme des m´ ethodes de projection. En notant Y k et L k deux
sous-espaces quelconques de R
n de dimension m, on appelle m´ ethode de projection un proc´ ed´ e qui construit une solution approch´ ee x
(k) ` a l’´ etape k, en
4 M´ ethodes it´ eratives pour la r´ esolution des syst` emes lin´ eaires
Programme 23 - gmres : M´ ethode GMRES pour la r´ esolution des syst` emes
lin´ eaires
function [x,iter]=gmres(A,b,x0,m,tol)
%GMRES M´ ethode GMRES.
% [X,ITER]=GMRES(A,B,X0,M,TOL) tente de r´ esoudre le syst` eme
% A*X=B avec la m´ ethode GMRES. TOL est la tol´ erance de la m´ ethode.
% M est la taille maximale de l’espace de Krylov. X0 est la donn´ ee
% initiale. ITER est l’it´ eration ` a laquelle la solution X a ´ et´ e calcul´ ee.
r0=b-A*x0; nr0=norm(r0,2);
if nr0 ˜= 0
v1=r0/nr0; V=[v1]; H=[]; iter=0; residual=1;
while iter <= m-1 & residual > tol,
[iter,V,H] = GSarnoldi(A,m,iter,V,H);
[nr,nc]=size(H); y=(H’*H) \ (H’*nr0*[1;zeros(nr-1,1)]);
x=x0+V(:,1:nc)*y; residual = norm(b-A*x,2)/nr0;
end
else
x=x0;
end
Pour am´ eliorer l’efficacit´ e de l’impl´ ementation de GMRES, il est n´ ecessaire
de d´ efinir un crit` ere d’arrˆ et qui ne requiert pas le calcul explicite du r´ esidu
` a chaque pas. Ceci est possible si on r´ esout de fa¸ con appropri´ ee le syst` eme
associ´ e ` a la matrice de Hessenberg
H k .
En pratique,
H k est transform´ e en une matrice triangulaire sup´ erieure
R k ∈ R
(k+1)×k avec r k+1,k = 0 telle que Q
T
k R k =
H k , o` u Q k est le r´ esultat du
produit de k rotations de Givens (voir Section 5.6.3). On peut alors montrer
que, Q k ´ etant orthogonale, minimiser
(0)
2 e 1 −
H k z
(k)
2 est ´ equivalent `
a
minimiser f k − R k z
(k)
2 , avec f k = Q k r
(0)
2 e 1 . On peut aussi montrer que
la valeur absolue de la k + 1-i` eme composante de f k est ´ egale ` a la norme
euclidienne du r´ esidu `
a l’it´ eration k.
Tout comme la m´ ethode d’Arnoldi, GMRES est coˆ uteuse en calcul et en
m´ emoire ` a moins que la convergence ne survienne qu’apr` es peu d’it´ erations.
Pour cette raison, on dispose `
a nouveau de deux variantes de l’algorithme :
la premi` ere, GMRES(m), bas´ ee sur un red´ emarrage apr` es m it´ erations, la
seconde, Quasi-GMRES ou QGMRES, sur l’arrˆ et du proc´ ed´ e d’orthogonalisation d’Arnoldi. Dans les deux cas, on perd la propri´ et´ e de GMRES d’obtenir
la solution exacte en un nombre fini d’it´ erations.
Remarque 4.4 (m´ ethodes de projection) Les it´ erations de Krylov peuvent ˆ etre vues comme des m´ ethodes de projection. En notant Y k et L k deux
sous-espaces quelconques de R
n de dimension m, on appelle m´ ethode de projection un proc´ ed´ e qui construit une solution approch´ ee x
(k) ` a l’´ etape k, en