function  X = subst ( W , IPIVOT, B)

%*****************************************************************************80
%
%MATLAB version of SUBROUTINE SUBST on page 164 of conte/deboor
%**************   I N P U T   ************
%  W, IPIVOT, N  are as on output from  F A C T O R, applied to the
%  matrix  A  of order  N .
%  B  is an N-vector giving the right side of the linear system to be solved
%**************  O U T P U T  ************
%  X  is the N-vector solving A*? = B .
%**************  M E T H O D  ************
% Algorithm 4.4 is used, i.e., the factorization of  A  contained in  W  
% and IPIVOT (as generated in  F A C T O R ) is used to solve  A*X = B
% for  X  by solving two triangular systems.
% 
% This version can handle a matrix  B
%
%  Reference:
%
%    Samuel Conte, Carl de Boor,
%    Elementary Numerical Analysis,
%    Third Edition,
%    SIAM, 2017,
%    ISBN: 978-1-611975-19-2.
%
  [N,~] = size(W);

  if N == 1
    X = B/W; 
    return
  end

  X = B(IPIVOT,:);
  for I=2:N
    X(I,:) = X(I,:) - W(I,1:I-1)*X(1:I-1,:);
  end

  X(N,:) = X(N,:)/W(N,N);
  for I=N-1:-1:1
    X(I,:) = (X(I,:) - W(I,I+1:N)*X(I+1:N,:))/W(I,I);
  end

  return
end

