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