function [ x, r ] = inpow ( A, xold )
%
% function [ x, r ] = inpow ( A, xold )
%
%  Discussion:
%
%    INPOW carries out one step of the inverse power method.
%    For efficiency, the matrix A should only be PLU factored ONCE.
%
%    We opt here for convenience, but that means we recompute the
%    same PLU factors every time we call.
%
%  Modified:
%
%    09 April 2000
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real A(N,N), the matrix.
%
%    Input, real XOLD(N), an estimate for the eigenvector.
%
%    Output, real X(N), the value of inv(A)*xold, normalized to have unit L2 norm.
%
%    Output, real R, the Rayleigh estimate of the eigenvalue.
%

%
%  Factor the matrix.
%
[ P, L, U ] = ge_pp ( A );
%
%  Multiply        x = inv ( A ) * xold  
%  <=>  
%  Solve:      A * x =             xold
%
x = plu_solve ( P, L, U, xold );

r = ( x' * A * x ) / ( x' * x );

x = x / norm ( x );
%
%  X and -X are both eigenvectors.  To make it easier to detect
%  convergence, pick the one whose maximum is positive.
%
x_max = max ( x );
if ( x_max < 0 )
  x = - x;
end

