% Bastian Groß
% Universität Trier
%
% 2. Juni 2014
%
% Präkonditioniertes konjugierte Gradientenverfahren

fid=fopen ('Programmieraufgabe_4.txt','w');
fprintf(fid, 'Präkonditioniertes konjugierte Gradientenverfahren\n\n');
fprintf(fid,'  n     mu     i_ges       Fehler               Zeit \n');

for n = [10,100,1000,2000]
e = ones(n,1);
A = spdiags([-e 2*e -e], -1:1, n, n);   %Problematrix 
b=10/n^3*[1:n]';                        %Rechte Seite
x0=zeros(n,1);                          %Startvektor


 
for mu = [0,1,2]
W= 0.5* (diag(diag(A)) + mu*(tril(A)-diag(diag(A))));
Aquer=inv(W)*A*inv(W)';
bquer=inv(W)*b;

if n < 2000
    Af=full (Aquer);
    lambda2=eig(Af);
    figure
    plot (lambda2);
end

%Präkonditionierung
g = Aquer*x0-bquer;
v = -g;
x = x0;
i = 0; %Iterationszähler


tic
while norm(g) > 10^(-12)
    Ahilf = Aquer*v;
    alpha = norm(g)^2/(v'*Ahilf);
    x     = x+ alpha*v;
    gneu  = g+alpha*Ahilf;
    beta  = norm(gneu)^2/norm(g)^2;
    g     = gneu;
    v     = -g+beta*v;
    i     = i+1;
    
end;
t=toc;

fprintf(fid,'%4.0f     %1.0f     %4.0f     %1.12e     %2.6f \n', n, mu, i , norm(A*x-b) , t);

end;

end;
fclose(fid);



