function [EVALUE,VECTOR,IFLAG] = invitr ( B , EGUESS , VGUESS )

%*****************************************************************************80
%
%MATLAB version of SUBROUTINE INVITR on page 194 of conte/deboor
% the function calls FACTOR, SUBST
%**************   I N P U T   ************
%  B  the matrix whose eigenvalue/vector is sought
%  EGUESS  a first guess for the eigenvalue
%  VGUESS  a first guess for the eigenvector
%**************  O U T P U T  ************
%  EVALUE  computed approximation to the eigenvalue
%  EVECTOR  computed approximation to the eigenvector
%  IFLAG, an integer
%       = 1  or  -1 (as set in FACTOR) indicates that all is well
%       = 0  indicates that something went wrong; see printed 
%            error message.
%**************  M E T H O D  ************
%  Inverse iteration, as described in the text, is used.

%**************
%  The following  T E R M I N A T I O N  P A R A M E T E R S are set here:
%  a tolerance  E P S L O N  on the difference between successive
%  eigenvalue iterates, and an upper bound  I T E R M X  on the number
%  of iteration steps. 
%
%  Reference:
%
%    Samuel Conte, Carl de Boor,
%    Elementary Numerical Analysis,
%    Third Edition,
%    SIAM, 2017,
%    ISBN: 978-1-611975-19-2.
%
  EPSLON = .000001; ITERMX = 20;
%
% put  B - (EGUESS)*IDENTITY  into  W
%
  [N,~] = size(B);
  W = B;
  diag = (1:N)*(N+1)-N;
  W(diag) = W(diag) - EGUESS;
  [W,IPIVOT,IFLAG] = factor_cd( W );

  if ~IFLAG
    disp(' eigenvalue guess too close; no eigenvalue calculated'),
    return
  end
%
% iteration starts here
%
  EVOLD = 0; records = zeros(ITERMX,N+1);  % for printout of iteratins
  for ITER=1:ITERMX
   VGUESS = VGUESS/norm(VGUESS);% normalize the current vector guess
   VECTOR = subst (W,IPIVOT,VGUESS); % get next vector guess
   EVNEW = VGUESS'*VECTOR;           % compute Rayleigh quotient 
   EVALUE = EGUESS + 1/EVNEW;       % add reciprocal to current eigenvalue guess
   records(ITER,:) = [EVALUE,VECTOR'];
       % stop iteration if current guess is close to previous eigenvalue guess 

   if abs(EVNEW-EVOLD) <= EPSLON*abs(EVNEW)
      disp('   EIGENVALUE        EIGENVECTOR COMPONENTS')
%     format SHORTE
      disp(records(1:ITER,:)), return
   end

   EVOLD = EVNEW; VGUESS = VECTOR;
  end

  IFLAG = 0;
  fprintf('No convergence to within %14.4e  after  %2d  steps', [EPSLON,ITERMX])

  return
end

