5.11 Méthode du gradient conjugué
169
solution exacte du système. Soit P la matrice tridiagonale dont les coefficients
diagonaux sont égaux à 2 et les coefficients sur- et sous-diagonaux sont égaux
à −1. Les matrices A et P sont toutes les deux symétriques définies positives.
On utilise le Programme 5.2 pour tester la méthode de Richardson dynamique
préconditionnée par P. On fixe tol=1.e-05, nmax=5000, x0=zeros(100,1). La
méthode converge en 43 itérations. Le même Programme 5.2 avec P=’G’ montre
que pour la méthode de Gauss-Seidel, 1658 itérations sont nécessaires pour
satisfaire le même critère d’arrêt.
5.11 Méthode du gradient conjugué
Dans une méthode itérative du type (5.58), la nouvelle itérée x
(k+1) est
obtenue en ajoutant à l’ancienne x
(k) un vecteur appelé direction de
descente, qui est soit le résidu r
(k) soit le résidu préconditionné z
(k) .
On peut se demander s’il ne serait pas possible de construire d’autres
directions de descente, p
(k) , qui permettraient de converger plus vite.
Quand la matrice A∈ R
n×n est symétrique définie positive, la méthode du gradient conjugué (en abrégé CG) utilise une suite de directions
de descente constituée par des vecteurs A-orthogonaux (on dit aussi Aconjugués), c’est-à-dire vérifiant ∀k ≥ 1,
(Ap
(j) )
T
p
(k+1) = 0,
j = 0, 1, . . ., k.
(5.61)
Pour tout vecteur x
(0) , après avoir posé r
(0) = b − Ax
(0) et p
(0) = r
(0) ,
la méthode du gradient conjugué s’écrit
pour k = 0, 1, . . .
α k =
p
(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) ,
β k =
(Ap
(k) )
T
r
(k+1)
(Ap (k) ) T p (k) ,
p
(k+1) = r
(k+1)
− β k p
(k)
(5.62)
Le paramètre α k permet de minimiser l’erreur e
(k+1)
A le long de la
direction de descente p
(k) , et β k est choisi pour que la nouvelle direction
p
(k+1) soit A-conjuguée avec p
(k) , c’est-à-dire (Ap
(k) )
T
p
(k+1) = 0. En
fait, on peut montrer par récurrence que si la dernière relation est satisfaite alors toutes les relations d’orthogonalité (5.61) pour j = 0, . . ., k −1
sont également satisfaites. On pourra trouver les détails de la construction de cette méthode dans [QSS07, Chapitre 4] ou [Saa03] par exemple.
169
solution exacte du système. Soit P la matrice tridiagonale dont les coefficients
diagonaux sont égaux à 2 et les coefficients sur- et sous-diagonaux sont égaux
à −1. Les matrices A et P sont toutes les deux symétriques définies positives.
On utilise le Programme 5.2 pour tester la méthode de Richardson dynamique
préconditionnée par P. On fixe tol=1.e-05, nmax=5000, x0=zeros(100,1). La
méthode converge en 43 itérations. Le même Programme 5.2 avec P=’G’ montre
que pour la méthode de Gauss-Seidel, 1658 itérations sont nécessaires pour
satisfaire le même critère d’arrêt.
5.11 Méthode du gradient conjugué
Dans une méthode itérative du type (5.58), la nouvelle itérée x
(k+1) est
obtenue en ajoutant à l’ancienne x
(k) un vecteur appelé direction de
descente, qui est soit le résidu r
(k) soit le résidu préconditionné z
(k) .
On peut se demander s’il ne serait pas possible de construire d’autres
directions de descente, p
(k) , qui permettraient de converger plus vite.
Quand la matrice A∈ R
n×n est symétrique définie positive, la méthode du gradient conjugué (en abrégé CG) utilise une suite de directions
de descente constituée par des vecteurs A-orthogonaux (on dit aussi Aconjugués), c’est-à-dire vérifiant ∀k ≥ 1,
(Ap
(j) )
T
p
(k+1) = 0,
j = 0, 1, . . ., k.
(5.61)
Pour tout vecteur x
(0) , après avoir posé r
(0) = b − Ax
(0) et p
(0) = r
(0) ,
la méthode du gradient conjugué s’écrit
pour k = 0, 1, . . .
α k =
p
(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) ,
β k =
(Ap
(k) )
T
r
(k+1)
(Ap (k) ) T p (k) ,
p
(k+1) = r
(k+1)
− β k p
(k)
(5.62)
Le paramètre α k permet de minimiser l’erreur e
(k+1)
A le long de la
direction de descente p
(k) , et β k est choisi pour que la nouvelle direction
p
(k+1) soit A-conjuguée avec p
(k) , c’est-à-dire (Ap
(k) )
T
p
(k+1) = 0. En
fait, on peut montrer par récurrence que si la dernière relation est satisfaite alors toutes les relations d’orthogonalité (5.61) pour j = 0, . . ., k −1
sont également satisfaites. On pourra trouver les détails de la construction de cette méthode dans [QSS07, Chapitre 4] ou [Saa03] par exemple.