% Script: krylovvar_example.m
% try different variants of Krylov subspace methods

% close all
figure
k = 0;
nn = [20 30 40 50];
tol = 10^-9;
for n = nn
  % Poisson poisson
   A = n^2*gallery('poisson',n);
   b = -ones(n^2,1);
  % Wathen matrix
  % A = n*gallery('wathen',n,n);
  % nnn = length(A);
  % b = -ones(nnn,1);
  % start solving!
  [x1,f1,rr1,it1,rv1] = gmres(A,b,30,tol);
  [x2,f2,rr2,it2,rv2] = cgs(A,b,tol,100);
  [x3,f3,rr3,it3,rv3] = bicgstab(A,b,tol,100);
  [x4,f4,rr4,it4,rv4] = qmr(A,b,tol,100);
  f1, f2, f3, f4
  k = k+1;
  subplot(2,2,k)
  semilogy(1:length(rv1),rv1/norm(b),'-k')
  hold on
  semilogy(1:length(rv2),rv2/norm(b),'--k')
  semilogy(1:length(rv3),rv3/norm(b),'-.k')
  semilogy(1:length(rv4),rv4/norm(b),'^k')
  axis('tight')
  hold off
  title(sprintf('n = %2.0f',n));
  if ((k == 3) | (k==4)) 
       xlabel('iteration');
  end
  if ((k == 3) | (k==1)) 
     ylabel('||b-Ax||/||b||');
  end
  set(gca,'YTick',[10^-9 10^-6 10^-3 10^0 10^3])
end