Implements inverse iteration to find an eigenvalue of A and its corresponding eigenvector Input arguments: A, square matrix s, estimate of eigenvalue tol, tolerance Output arguments: v, eigenvector of A c, corresponding eigenvalue k, number of iterations
0001 function [ v,c,k ] = inverse_iteration( A,s,tol ) 0002 % Implements inverse iteration to find an eigenvalue of A 0003 % and its corresponding eigenvector 0004 % Input arguments: 0005 % A, square matrix 0006 % s, estimate of eigenvalue 0007 % tol, tolerance 0008 % Output arguments: 0009 % v, eigenvector of A 0010 % c, corresponding eigenvalue 0011 % k, number of iterations 0012 0013 [n,m]=size(A); % finding the size of A 0014 if n~= m; 0015 error('input is not a square matrix'); 0016 elseif tol<=0; 0017 error('tolerance should be positive'); 0018 end 0019 0020 k=0; % initializing k 0021 0022 while k>=0; 0023 [L,U,P]=LU_pivot(A-s*eye(n)); % calculate the LU factorization 0024 % with pivoting of (A-sI) 0025 d=diag(U); 0026 if min(abs(d))<=1e-12; % if any diagonal elements of U are close to zero 0027 c=s; % accept s as the eigenvalue 0028 v=U; % set output v to U 0029 disp(['eigenvalue found; v lies in the null space of '... 0030 'given matrix']); 0031 return; % return function 0032 else 0033 if k==0; 0034 [~,i]=min(abs(d)); % i is the index of the smallest entry 0035 % along the diagonal of U 0036 ei=zeros(n,1); ei(i)=1; % initialize ith standard basis vector 0037 x0=L*ei; % initialize x0 0038 x1=Backward(U,ei); % calculate x1 0039 x1=P*x1; x0=P*x0; % apply pivoting 0040 else 0041 y=Forward(L,x0); % calculate x1 from (A-sI)x1=x0 0042 x1=Backward(U,y); 0043 x1=P*x1; % apply pivoting 0044 end 0045 0046 mu=x1'*(A-s*eye(n))*x1/(x1'*x1); % calculate mu which minimizes 0047 % |x0-mu*x1| 0048 z=x0-mu*x1; % if |z|<tol|x1| we terminate procedure 0049 if sqrt(z'*z)<tol*sqrt(x1'*x1); 0050 c=mu+s; % accept mu+s as eigenvalue 0051 v=x1/sqrt(x1'*x1); % normalize x1 and accept as eigenvector 0052 return; % return function 0053 else 0054 x0=x1/sqrt(x1'*x1); % normalize x 0055 s=s+mu; % adjust s 0056 k=k+1; % increment k 0057 end 0058 end 0059 end 0060 end 0061 0062 0063