function [ wval, dwdx ] = r8poly_lagrange_factor ( npol, xpol, xval )

%*****************************************************************************80
%
%% r8poly_lagrange_factor() evaluates the polynomial Lagrange factor at a point.
%
%  Discussion:
%
%    The power sum form of a polynomial is:
%
%      p(x) = c(0) + c(1) * x + ... + c(n-1) * x^(n-1) + c(n) * x^(n)
%
%    Because MATLAB doesn't allow 0 indexing, we use an adjusted formula:
%
%      p(x) = c(1) + c(2) * x + ... + c(n) * x^(n-1) + c(n+1) * x^(n)
%
%    W(X) = Product ( 1 <= I <= NPOL ) ( X - XPOL(I) )
%
%    Suppose F(X) is at least N times continuously differentiable in the
%    interval [A,B].  Pick NPOL distinct points XPOL(I) in [A,B] and compute
%    the interpolating polynomial P(X) of order NPOL ( and degree NPOL-1)
%    which passes through all the points ( XPOL(I), F(XPOL(I)) ).
%    Then in the interval [A,B], the maximum error
%
%      abs ( F(X) - P(X) )
%
%    is bounded by:
%
%      C * FNMAX * W(X)
%
%    where
%
%      C is a constant,
%      FNMAX is the maximum value of the NPOL-th derivative of F in [A,B],
%      W(X) is the Lagrange factor.
%
%    Thus, the value of W(X) is useful as part of an estimated bound
%    for the interpolation error.
%
%    Note that the Chebyshev abscissas have the property that they minimize
%    the value of W(X) over the interval [A,B].  Hence, if the abscissas may
%    be chosen arbitrarily, the Chebyshev abscissas have this advantage over
%    other choices.
%
%    For a set of points XPOL(I), 1 <= I <= NPOL, the IPOL-th Lagrange basis
%    polynomial L(IPOL)(X), has the property:
%
%      L(IPOL)( XPOL(J) ) = delta ( IPOL, J )
%
%    and may be expressed as:
%
%      L(IPOL)(X) = W(X) / ( ( X - XPOL(IPOL) ) * W'(XPOL(IPOL)) )
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    26 April 2005
%
%  Author:
%
%    John Burkardt
%
%  Input:
%
%    integer NPOL, the number of abscissas.
%    NPOL must be at least 1.
%
%    real XPOL(NPOL), the abscissas, which should 
%    be distinct.
%
%    real XVAL, the evaluation point.
%
%  Output:
%
%    real WVAL, the value of the Lagrange factor at XVAL.
%
%    real DWDX, the derivative of W with respect to XVAL.
%
  wval = prod ( xval - xpol(1:npol) );

  dwdx = 0.0;

  for i = 1 : npol

    term = 1.0;

    for j = 1 : npol
      if ( i ~= j )
        term = term * ( xval - xpol(j) );
      end
    end

    dwdx = dwdx + term;

  end

  return
end
