function [v,res] = powerit2(A, nit, tol, v) 
%% [v,err] = powerit(A, nit, tol, v0) 
%% does nit steps of the power method 
%% starging with (input) v
%% 
 for i=1:nit
     v = A*v;
     [t,k] = max(abs(v));
     alpha = v(k);
     if (alpha == 0) 
        error('ERR. division by 0 -- STOP') 
     else 
       v = v / alpha;
       res = norm(A*v - alpha*v) ;
fprintf(1,...
' iter %3d -> res=%12.4e alph=%12.4e\n',...
  i,res,alpha);
     end
     if (res < tol)
        return
     end
 end

