chenmh_eigifp.m

Frequently used algorithms for numerical analysis and signal processing

function [x,s,e] = Chenmh_Eigifp(A,B,m,tol)
% m is the dimension of the Krylov subspace.
n = size(A,1); 
x = randn(n,1); 
s = (x'*A*x)/(x'*B*x);
e = norm(A*x-s*B*x)/norm(x);
while (e(end)>tol)
    C = A - s*B;
    V = x/norm(x); H = []; istop = 0; i = 0;
    while (i<m & istop==0)
        [V,H,i,istop] = Chenmh_ArnoldiPro(C,n,i,V,H);
    end
    H_s = H(1:end-1,:); V_s = V(:,1:end-1); 
%     [V_s'*C*V_s H_s]
%     pause
    [U,Sig] = eig(H_s,V_s'*B*V_s);
    [val,ndx] = min(diag(Sig));
    x = V_s*U(:,ndx);
    s = s + val;
    e = [e; norm(A*x-s*B*x)/norm(x)];
    if istop==1 break; end;
end