From 1077fa4d772832f77a677ce7fb7c2d513b959e3f Mon Sep 17 00:00:00 2001 From: Eric Andersen Date: Thu, 10 May 2001 00:40:28 +0000 Subject: uClibc now has a math library. muahahahaha! -Erik --- libm/float/spencef.c | 135 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 135 insertions(+) create mode 100644 libm/float/spencef.c (limited to 'libm/float/spencef.c') diff --git a/libm/float/spencef.c b/libm/float/spencef.c new file mode 100644 index 000000000..52799babe --- /dev/null +++ b/libm/float/spencef.c @@ -0,0 +1,135 @@ +/* spencef.c + * + * Dilogarithm + * + * + * + * SYNOPSIS: + * + * float x, y, spencef(); + * + * y = spencef( x ); + * + * + * + * DESCRIPTION: + * + * Computes the integral + * + * x + * - + * | | log t + * spence(x) = - | ----- dt + * | | t - 1 + * - + * 1 + * + * for x >= 0. A rational approximation gives the integral in + * the interval (0.5, 1.5). Transformation formulas for 1/x + * and 1-x are employed outside the basic expansion range. + * + * + * + * ACCURACY: + * + * Relative error: + * arithmetic domain # trials peak rms + * IEEE 0,4 30000 4.4e-7 6.3e-8 + * + * + */ + +/* spence.c */ + + +/* +Cephes Math Library Release 2.1: January, 1989 +Copyright 1985, 1987, 1989 by Stephen L. Moshier +Direct inquiries to 30 Frost Street, Cambridge, MA 02140 +*/ + +#include + +static float A[8] = { + 4.65128586073990045278E-5, + 7.31589045238094711071E-3, + 1.33847639578309018650E-1, + 8.79691311754530315341E-1, + 2.71149851196553469920E0, + 4.25697156008121755724E0, + 3.29771340985225106936E0, + 1.00000000000000000126E0, +}; +static float B[8] = { + 6.90990488912553276999E-4, + 2.54043763932544379113E-2, + 2.82974860602568089943E-1, + 1.41172597751831069617E0, + 3.63800533345137075418E0, + 5.03278880143316990390E0, + 3.54771340985225096217E0, + 9.99999999999999998740E-1, +}; + +extern float PIF, MACHEPF; + +/* pi * pi / 6 */ +#define PIFS 1.64493406684822643647 + + +float logf(float), polevlf(float, float *, int); +float spencef(float xx) +{ +float x, w, y, z; +int flag; + +x = xx; +if( x < 0.0 ) + { + mtherr( "spencef", DOMAIN ); + return(0.0); + } + +if( x == 1.0 ) + return( 0.0 ); + +if( x == 0.0 ) + return( PIFS ); + +flag = 0; + +if( x > 2.0 ) + { + x = 1.0/x; + flag |= 2; + } + +if( x > 1.5 ) + { + w = (1.0/x) - 1.0; + flag |= 2; + } + +else if( x < 0.5 ) + { + w = -x; + flag |= 1; + } + +else + w = x - 1.0; + + +y = -w * polevlf( w, A, 7) / polevlf( w, B, 7 ); + +if( flag & 1 ) + y = PIFS - logf(x) * logf(1.0-x) - y; + +if( flag & 2 ) + { + z = logf(x); + y = -0.5 * z * z - y; + } + +return( y ); +} -- cgit v1.2.3