v =  randn(n,1); v = v/norm(v);
 K = zeros(n,k);
 K(:,1) = v; 
 resid = zeros(k,1);
 for j = 1:k,
    w = A*v;
    v = w/norm(w);
    K(:,j) = v;
 end
 format long e
 [V,R] = qr(K);
 H = V'*(A*V);
 F = A*V - V*H;
 format short e
 Hval = H

 for j = 1:k,
    resid(j) = norm(F(:,j));
 end
 plot(resid,'*')
 hold
 plot(resid)
 hold