Je vous présente ici la méthode du gradient conjugué. Cette méthode numérique permet de résoudre des grands systèmes linéaires dont la matrice est symétrique définie positive. Elle repose sur la recherche de directions successives permettant d’atteindre la solution exacte du système étudié.

Description du problème

Considérons le système suivant:

\[Ax=b,\]

où $A$ est une matrice de taille $n\times n$ symétrique définie positive ($A^\top =A$ et $x^\top Ax>0$, pour tout vecteur $x\in \mathbf{R}^n$ non nul). Soit $x_\star$ la solution exacte de ce système.

Directions conjuguées

Comme la matrice $A$ est symétrique définie positive, on peut définir le produit scalaire suivant sur $\mathbf{R}^n$: $\langle u,v \rangle_A = u^\top A v.$ Deux éléments $u,v\in \mathbf{R}^n$ sont dit $A$-conjugués si:

\[u^\top A v = 0\]

La méthode du gradient conjugué consiste à  construire une suite $(p_k)_{k\in\mathbf{N}^*}$ de $n$ vecteurs $A$-conjugués. Dès lors, la suite $p_1,p_2,\ldots,p_n$ forme une base de $\mathbf{R}^n$. La solution exacte $x_\star$ peut donc se décomposer comme suit:

\[x_\star = \alpha_1 p_1 + \cdots + \alpha_n p_n\]

où $\alpha_k = \dfrac{p_k^\top b}{p_k^\top A p_k},k=1,\ldots,n.$

Construction des directions conjuguées

La solution exacte $x_\star$ peut-àªtre également vu comme l’unique minimisant de la fonctionnelle

\[J(x) = \frac12 x^\top A x - b^\top x, \quad x\in\mathbf{R}^n\]

On a donc clairement

\(\nabla J(x) = A x - b, \quad x\in\mathbf{R}^n\) d’où \(\nabla J(x_\star)=0_{\mathbf{R}^n}\)

On définit le résidu du système d’équation comme suit

\[r_k = b - Ax_k=-\nabla J(x_k)\]

$r_k$ représente donc la direction du gradient de la fonctionnelle $J$ en $x_k$ (à  un signe près).

La nouvelle direction de descente $p_{k+1}$ suit donc celle du résidu modulo, sa $A$-conjugaison avec $p_k$, on a alors:

\[p_{k+1} = r_k - \frac{p_k^\top A r_k}{p_k^\top A p_k} p_k\]

C’est le choix du coefficient $\dfrac{p_k^\top A r_k}{p_k^\top A p_k}$ qui assure la $A$-conjugaison des directions $p_k$. Pour vous en assurez calculer $(Ap_{k+1},p_k)$, cette quantité est nulle !

Algorithme du gradient conjugué

Calcul du résidu $r_0=b-Ax_0$ pour un vecteur $x_0$ quelconque. On fixe $p_0=r_0$.

Pour $k=1,2,\ldots$

\[q_k=A p_k\]

If $k=1$

\[p_1=r_{0}\]

Else

\[\beta_{k-1}= \frac{r_{k-1}^\top r_{k-1}}{r_{k-2}^\top r_{k-2}}\] \[p_k=r_{k-1} + \beta_{k-1}p_{k-1}\]

EndIf

\[\alpha_k=\frac{r_{k-1}^\top r_{k-1}}{p_k^\top q_k}\] \[x_k=x_{k-1}+\alpha_k p_k\] \[r_k=r_{k-1}-\alpha_k q_k\]

FinPour