%=== Nettoyage
clc;
clear all;
close all;

%=== Parametres
tolerance=1E-10;    %pour le test de convergence sur la norme du residu
niter_max=100000;

% %=== Creation d'un nouveau systeme
% n=100;              %dimension du systeme
% A=rand(n,n);        %matrice carree aleatoire a valeurs dans [-1;1]
% A=A*A';             %matrice symétrique définie positive
% A=A+2*eye(n,n);     %amélioration du conditionnement
% b=rand(n,1);        %vecteur des inconnues
% nomA=strcat('A',num2str(n));
% nomb=strcat('b',num2str(n));
% save(strcat(nomA,'.mat'),'A');
% save(strcat(nomb,'.mat'),'b');
% pfile=fopen(nomA,'wb');%ecriture de A
% fwrite(pfile,n);
% fwrite(pfile,A,'double');
% fclose(pfile);
% pfile=fopen(nomb,'wb');%ecriture de b
% fwrite(pfile,n);
% fwrite(pfile,b,'double');
% fclose(pfile);


%=== Chargement d'un systeme existant
nomA='A100';
nomb='b100';
% nomA='A200';
% nomb='b200';
pfile=fopen(nomA,'rb');%lecture A
n=fread(pfile,1,'uint8');
dataA=zeros(n*n,1);
dataA=fread(pfile,n*n,'double');
fclose(pfile);
A=zeros(n,n);
for i=1:n
    A(:,i)=dataA((i-1)*n+1:1:i*n,1);
end
pfile=fopen(nomb,'rb');%lecture b
n=fread(pfile,1,'uint8');
b=zeros(n,1);
b=fread(pfile,n*n,'double');
fclose(pfile);


%=== Resolutions
x0=zeros(n,1);  %estimation initiale de la solution
[xSD,residusSD]=        SD(A,b,x0,niter_max,tolerance);
%[xCG1,residusCG1]=      CGk(A,b,x0,tolerance,1);
%[xCG2,residusCG2]=      CGk(A,b,x0,tolerance,2);
%[xCGall,residusCGall]=  CGall(A,b,x0,tolerance);

%=== Affichage de l evolution des normes de residus en fct des iterations
figure(1)
hold on
plot([1:1:length(residusSD)]',residusSD,'-b');
%plot([0:1:n]',residusCG1,'-r');
%plot([0:1:n]',residusCG2,'-c');
%plot([0:1:n]',residusCGall,'-g');
grid;
set(gca,'xscale','log');
set(gca,'yscale','log');
set(gca,'yMinorGrid','on');
legend('SD');
xlabel('nb d iterations');
ylabel('norme de b-A*x');
hold off