%%%%%%%%%%%%%%%%%%%% Aufgabe 2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function s5_2

    A=[ -9  -14  5  21;
	-15   5 -3  18;
	-15  0 -1 21;
	-9 26 -11 -3;
	3 -27 11 12];
format long; format compact;
[U,S,V]=svd(A);
S=diag(S); e=eps;
if(S(1)<e) rho=0;
else
for i=1:length(S) if (S(i)/S(1) < e) rho = i-1; break; end; end;
end
rank=rho,
n=4;
KA=V(:,rank+1:n),    % b) o.n. basis of kern(A)
     RA=U(:,1:rank), % c) o.n. basis of im(A) 


%%%%%%%%%%%%%%%%%%%%%%%% Aufgabe 3 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5

function s5a3(n,its)

% solving eigenvalue problem with inverse vector iterations
% in : n    meshwidth is h=1/(2n)
%      its  number of iterations

% check input, set defaults

  if (nargin < 2) its=150; end
if (nargin < 1) n=3; end
  if(n<3 | mod(n,1)~=0) error('choose n>2 integer'); end

M=2*n-1; N=n-1;  % interior grid is of size M x N
  NB = N; NA = M*N;  %sizes of B and A
h=1/(2*n);

%generate A
  d=4*ones(NA,1);   % diagonal
du=-ones(NA-1,1); du(NB:NB:NA-n+1)=0; %upper and lower diagonal
  do=-ones((M-1)*NB,1);    % off blocks
A=sparse(diag(d)+diag(du,1)+diag(du,-1)+diag(do,n-1)+diag(do,-(n-1)));


   U=eye(NA,1);   %start for iteration
for i=1:its
  U=A\U;           %inverse iteration
end


lambda=U'*A*U/(U'*U);
test=lambda/(h*pi)^2; % analytical eigenvalue
disp([n,lambda,test]);


% plot the eigenmode
a=0.5; b=1; x=0:h:a; y=0:h:b;
[X,Y] = meshgrid(0:h:a,0:h:b);
Z = reshape(U,N,M)';
Z=[zeros(1,N+2); zeros(M,1) Z zeros(M,1); zeros(1,N+2)];
surf(X,Y,Z);
title(['n=', int2str(n)]);
xlabel('x');ylabel('y');
zlabel('first eigenmode');


 
