disp('Ci-contre une illustration des trois méthodes.')
disp('On dessine les approximations successives de la solution, avec :')
nmax=11; tol=1.e-6;
A=[2,1; 1,2], B=[1;2], X0=zeros(2,1);
//Méthode du gradient à pas constant
alpha=0.5; X=X0; R=B-A*X; p1=[];
for n=0:nmax
p1=[p1,X];
if (norm(R)<=tol), break, end
X=X+alpha*R;
R=B-A*X;
end
p1=[X0*ones(1,nmax+1-size(p1,2)),p1];
//Méthode du gradient à pas variable
X=X0; R=B-A*X; p2=[];
for n=0:nmax
p2=[p2,X];
if (norm(R)<=tol), break, end
AR=A*R;
alpha=(R'*R)/(AR'*R);
X=X+alpha*R;
R=R-alpha*AR;
end
p2=[X0*ones(1,nmax+1-size(p2,2)),p2];
//Méthode du gradient conjugué
X=X0; R=B-A*X; P=R; p3=[];
r2=R'*R;
for n=0:nmax
p3=[p3,X];
if (r2<=tol^2), break, end
AP=A*P;
alpha=r2/(AP'*P);
X=X+alpha*P;
R=R-alpha*AP;
s2=R'*R;
P=R+(s2/r2)*P;
r2=s2;
end
p3=[X0*ones(1,nmax+1-size(p3,2)),p3];
xset("window",2), xbasc()
plot2d([p1(1,:)',p2(1,:)',p3(1,:)'],[p1(2,:)',p2(2,:)',p3(2,:)'],[2,3,5],'131',..
leg="gradient à pas constant@gradient à pas variable@gradient conjugué",..
rect=[0 0 0.5 1],nax=[1 15 1 10])
xtitle('Résolution de AX=B avec A=[2,1; 1,2], B=[1;2] et X0=[0;0]','x','y')