function fx = p23_fx ( n, x )

%*****************************************************************************80
%
%% p23_fx() evaluates the function for problem 23.
%
%    The Chandrasekhar H function is
%
%      F(H)(mu) = H(mu)
%        - 1 / ( 1 - c/2
%        * Integral ( 0 <= nu <= 1 ) [ mu * H(mu) / ( mu + nu ) ] dnu )
%
%    Applying the composite midpoint rule, and setting, for i = 1 to N,
%
%      mu(i) = ( i - 0.5 ) / N
%
%    we have the discrete function
%
%      F(h)(i) = h(i)
%        - 1 / ( 1 - c/(2N)
%        * sum ( 1 <= j <= N ) ( mu(i) * h(j) / ( mu(i) + mu(j) ) )
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    10 October 2020
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Subramanyan Chandrasekhar,
%    Radiative Transfer,
%    Dover, 1960,
%    ISBN13: 978-0486605906,
%    LC: QB461.C46.
%
%  Input:
%
%    integer N, the number of variables and equations.
%
%    real X(N), the evaluation point.
%
%  Output:
%
%    real FX(N), the function value.
%
%  Local:
%
%    real C, the constant, which must be between 0 and 1.
%
  fx = zeros ( n, 1 );

  c = 0.9;

  fx(1:n) = x(1:n);

  for i = 1 : n
    mu(i) = ( 2 * i - 1 ) / ( 2 * n );
  end

  for i = 1 : n

    s = 0.0;
    for j = 1 : n
      s = s + ( mu(i) * x(j) ) / ( mu(i) + mu(j) );
    end

    term = 1.0 - c * s / ( 2 * n );

    fx(i) = fx(i) - 1.0 / term;

  end

  return
end

