   function [V, lam] = Lan(A, nev, msteps,orient)
%% function [V, lam] = Lan(A, nev, msteps,orient)
%% ---------computes very approximate eigenvectors
%%                     via lanczos
  [m, n] = size(A);
  v  = randn(m,1);
  t  = norm(v,2);
  v  =  v/t; 
  beta = 0; 
  V = [v]; 
  vold = v;
  orthTol = 1.e-08;
  wn = 0.0 ;
%%-------------------- main loop
 for k=1:msteps
    w = A*v; 
    w = w - beta*vold ; 
    alpha = w'*v; 
    w = w - alpha*v;
%%-------------------- Reorthogonalization: 
if (0) 
    t = V'*w;
%% fprintf(1,' it %4d inPr %e \n',k,norm(t,1));
      w = w - V* t; 
end
    beta = w'*w;
%%-------------------- test for lucky breakd. 
    wn = wn + alpha*alpha;
    if (beta*k < orthTol*wn) 
        break
    end
%%-------------------- params for orth. test
    wn   = wn+2.0*beta;
    beta = sqrt(beta) ;
    vold = v; 
    v = w/ beta;
    V(:,k+1) = v;
%%-------------------- tridiag matrix 
    T(k,k)   = alpha; 
    T(k,k+1) = beta; 
    T(k+1,k) = beta; 
 end
%%--------------------
 [U,lam] = eig(T(1:k,1:k));
 [lam, indx] = sort(diag(lam),orient);
%%-------------------- return only nev eigenpairs
 lam = lam(1:nev);
 V = V(:,1:k)*U(:,indx(1:nev));


