%Méthode du gradient conjugué
function [X,n]=gc(A,B,X0,nmax,tol)

X=X0; R=B-A*X; P=R;
r2=R'*R; %le test porte sur la norme au carré (qui intervient ailleurs)
for n=0:nmax
  if (r2<=tol^2), break, end
  AP=A*P; %on va calculer un produit de matrices en moins
  alpha=r2/(AP'*P);
  X=X+alpha*P;
  R=R-alpha*AP; %on a moins d'opérations qu'avec "R=B-A*X"
  s2=R'*R; %on calcule la norme au carré du résidu suivant
  P=R+(s2/r2)*P;
  r2=s2;
end