function value = r8_lnrel ( x )

%*****************************************************************************80
%
%% R8_LNREL evaluates log ( 1 + X ) for an R8 argument.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    01 October 2011
%
%  Author:
%
%    Original FORTRAN77 version by Wayne Fullerton.
%    MATLAB version by John Burkardt.
%
%  Reference:
%
%    Wayne Fullerton,
%    Portable Special Function Routines,
%    in Portability of Numerical Software,
%    edited by Wayne Cowell,
%    Lecture Notes in Computer Science, Volume 57,
%    Springer 1977,
%    ISBN: 978-3-540-08446-4,
%    LC: QA297.W65.
%
%  Parameters:
%
%    Input, real X, the argument.
%
%    Output, real VALUE, the value of LOG ( 1 + X ).
%
  persistent alnrcs
  persistent nlnrel
  persistent xmin

  if ( isempty ( nlnrel ) )

    alnrcs = [ ...
      +0.10378693562743769800686267719098E+01, ...
      -0.13364301504908918098766041553133, ...
      +0.19408249135520563357926199374750E-01, ...
      -0.30107551127535777690376537776592E-02, ...
      +0.48694614797154850090456366509137E-03, ...
      -0.81054881893175356066809943008622E-04, ...
      +0.13778847799559524782938251496059E-04, ...
      -0.23802210894358970251369992914935E-05, ...
      +0.41640416213865183476391859901989E-06, ...
      -0.73595828378075994984266837031998E-07, ...
      +0.13117611876241674949152294345011E-07, ...
      -0.23546709317742425136696092330175E-08, ...
      +0.42522773276034997775638052962567E-09, ...
      -0.77190894134840796826108107493300E-10, ...
      +0.14075746481359069909215356472191E-10, ...
      -0.25769072058024680627537078627584E-11, ...
      +0.47342406666294421849154395005938E-12, ...
      -0.87249012674742641745301263292675E-13, ...
      +0.16124614902740551465739833119115E-13, ...
      -0.29875652015665773006710792416815E-14, ...
      +0.55480701209082887983041321697279E-15, ...
      -0.10324619158271569595141333961932E-15, ...
      +0.19250239203049851177878503244868E-16, ...
      -0.35955073465265150011189707844266E-17, ...
      +0.67264542537876857892194574226773E-18, ...
      -0.12602624168735219252082425637546E-18, ...
      +0.23644884408606210044916158955519E-19, ...
      -0.44419377050807936898878389179733E-20, ...
      +0.83546594464034259016241293994666E-21, ...
      -0.15731559416479562574899253521066E-21, ...
      +0.29653128740247422686154369706666E-22, ...
      -0.55949583481815947292156013226666E-23, ...
      +0.10566354268835681048187284138666E-23, ...
      -0.19972483680670204548314999466666E-24, ...
      +0.37782977818839361421049855999999E-25, ...
      -0.71531586889081740345038165333333E-26, ...
      +0.13552488463674213646502024533333E-26, ...
      -0.25694673048487567430079829333333E-27, ...
      +0.48747756066216949076459519999999E-28, ...
      -0.92542112530849715321132373333333E-29, ...
      +0.17578597841760239233269760000000E-29, ...
      -0.33410026677731010351377066666666E-30, ...
      +0.63533936180236187354180266666666E-31 ]';
  
    nlnrel = r8_inits ( alnrcs, 43, 0.1 * r8_mach ( 3 ) );
    xmin = - 1.0 + sqrt ( r8_mach ( 4 ) );

  end

  if ( x <= - 1.0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'R8_LNREL - Fatal error!\n' );
    fprintf ( 1, '  X <= -1.\n' );
    error ( 'R8_LNREL - Fatal error!' )
  elseif ( x < xmin )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'R8_LNREL - Warning!\n' );
    fprintf ( 1, '  Result is less than half precision.\n' );
    fprintf ( 1, '  X is too close to - 1.\n' );
  end

  if ( abs ( x ) <= 0.375 )
    value = x * ( 1.0 - x * r8_csevl ( x / 0.375, alnrcs, nlnrel ) );
  else
    value = log ( 1.0 + x );
  end

  return
end
