From 7ce331c01ce6eb7b3f5c715a38a24359da9c6ee2 Mon Sep 17 00:00:00 2001 From: Eric Andersen Date: Thu, 22 Nov 2001 14:04:29 +0000 Subject: Totally rework the math library, this time based on the MacOs X math library (which is itself based on the math lib from FreeBSD). -Erik --- libm/double/ellpk.c | 234 ---------------------------------------------------- 1 file changed, 234 deletions(-) delete mode 100644 libm/double/ellpk.c (limited to 'libm/double/ellpk.c') diff --git a/libm/double/ellpk.c b/libm/double/ellpk.c deleted file mode 100644 index 8b36690e2..000000000 --- a/libm/double/ellpk.c +++ /dev/null @@ -1,234 +0,0 @@ -/* ellpk.c - * - * Complete elliptic integral of the first kind - * - * - * - * SYNOPSIS: - * - * double m1, y, ellpk(); - * - * y = ellpk( m1 ); - * - * - * - * DESCRIPTION: - * - * Approximates the integral - * - * - * - * pi/2 - * - - * | | - * | dt - * K(m) = | ------------------ - * | 2 - * | | sqrt( 1 - m sin t ) - * - - * 0 - * - * where m = 1 - m1, using the approximation - * - * P(x) - log x Q(x). - * - * The argument m1 is used rather than m so that the logarithmic - * singularity at m = 1 will be shifted to the origin; this - * preserves maximum accuracy. - * - * K(0) = pi/2. - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC 0,1 16000 3.5e-17 1.1e-17 - * IEEE 0,1 30000 2.5e-16 6.8e-17 - * - * ERROR MESSAGES: - * - * message condition value returned - * ellpk domain x<0, x>1 0.0 - * - */ - -/* ellpk.c */ - - -/* -Cephes Math Library, Release 2.8: June, 2000 -Copyright 1984, 1987, 2000 by Stephen L. Moshier -*/ - -#include - -#ifdef DEC -static unsigned short P[] = -{ -0035020,0127576,0040430,0051544, -0036025,0070136,0042703,0153716, -0036402,0122614,0062555,0077777, -0036441,0102130,0072334,0025172, -0036341,0043320,0117242,0172076, -0036312,0146456,0077242,0154141, -0036420,0003467,0013727,0035407, -0036564,0137263,0110651,0020237, -0036775,0001330,0144056,0020305, -0037305,0144137,0157521,0141734, -0040261,0071027,0173721,0147572 -}; -static unsigned short Q[] = -{ -0034366,0130371,0103453,0077633, -0035557,0122745,0173515,0113016, -0036302,0124470,0167304,0074473, -0036575,0132403,0117226,0117576, -0036703,0156271,0047124,0147733, -0036766,0137465,0002053,0157312, -0037031,0014423,0154274,0176515, -0037107,0177747,0143216,0016145, -0037217,0177777,0172621,0074000, -0037377,0177777,0177776,0156435, -0040000,0000000,0000000,0000000 -}; -static unsigned short ac1[] = {0040261,0071027,0173721,0147572}; -#define C1 (*(double *)ac1) -#endif - -#ifdef IBMPC -static unsigned short P[] = -{ -0x0a6d,0xc823,0x15ef,0x3f22, -0x7afa,0xc8b8,0xae0b,0x3f62, -0xb000,0x8cad,0x54b1,0x3f80, -0x854f,0x0e9b,0x308b,0x3f84, -0x5e88,0x13d4,0x28da,0x3f7c, -0x5b0c,0xcfd4,0x59a5,0x3f79, -0xe761,0xe2fa,0x00e6,0x3f82, -0x2414,0x7235,0x97d6,0x3f8e, -0xc419,0x1905,0xa05b,0x3f9f, -0x387c,0xfbea,0xb90b,0x3fb8, -0x39ef,0xfefa,0x2e42,0x3ff6 -}; -static unsigned short Q[] = -{ -0x6ff3,0x30e5,0xd61f,0x3efe, -0xb2c2,0xbee9,0xf4bc,0x3f4d, -0x8f27,0x1dd8,0x5527,0x3f78, -0xd3f0,0x73d2,0xb6a0,0x3f8f, -0x99fb,0x29ca,0x7b97,0x3f98, -0x7bd9,0xa085,0xd7e6,0x3f9e, -0x9faa,0x7b17,0x2322,0x3fa3, -0xc38d,0xf8d1,0xfffc,0x3fa8, -0x2f00,0xfeb2,0xffff,0x3fb1, -0xdba4,0xffff,0xffff,0x3fbf, -0x0000,0x0000,0x0000,0x3fe0 -}; -static unsigned short ac1[] = {0x39ef,0xfefa,0x2e42,0x3ff6}; -#define C1 (*(double *)ac1) -#endif - -#ifdef MIEEE -static unsigned short P[] = -{ -0x3f22,0x15ef,0xc823,0x0a6d, -0x3f62,0xae0b,0xc8b8,0x7afa, -0x3f80,0x54b1,0x8cad,0xb000, -0x3f84,0x308b,0x0e9b,0x854f, -0x3f7c,0x28da,0x13d4,0x5e88, -0x3f79,0x59a5,0xcfd4,0x5b0c, -0x3f82,0x00e6,0xe2fa,0xe761, -0x3f8e,0x97d6,0x7235,0x2414, -0x3f9f,0xa05b,0x1905,0xc419, -0x3fb8,0xb90b,0xfbea,0x387c, -0x3ff6,0x2e42,0xfefa,0x39ef -}; -static unsigned short Q[] = -{ -0x3efe,0xd61f,0x30e5,0x6ff3, -0x3f4d,0xf4bc,0xbee9,0xb2c2, -0x3f78,0x5527,0x1dd8,0x8f27, -0x3f8f,0xb6a0,0x73d2,0xd3f0, -0x3f98,0x7b97,0x29ca,0x99fb, -0x3f9e,0xd7e6,0xa085,0x7bd9, -0x3fa3,0x2322,0x7b17,0x9faa, -0x3fa8,0xfffc,0xf8d1,0xc38d, -0x3fb1,0xffff,0xfeb2,0x2f00, -0x3fbf,0xffff,0xffff,0xdba4, -0x3fe0,0x0000,0x0000,0x0000 -}; -static unsigned short ac1[] = { -0x3ff6,0x2e42,0xfefa,0x39ef -}; -#define C1 (*(double *)ac1) -#endif - -#ifdef UNK -static double P[] = -{ - 1.37982864606273237150E-4, - 2.28025724005875567385E-3, - 7.97404013220415179367E-3, - 9.85821379021226008714E-3, - 6.87489687449949877925E-3, - 6.18901033637687613229E-3, - 8.79078273952743772254E-3, - 1.49380448916805252718E-2, - 3.08851465246711995998E-2, - 9.65735902811690126535E-2, - 1.38629436111989062502E0 -}; - -static double Q[] = -{ - 2.94078955048598507511E-5, - 9.14184723865917226571E-4, - 5.94058303753167793257E-3, - 1.54850516649762399335E-2, - 2.39089602715924892727E-2, - 3.01204715227604046988E-2, - 3.73774314173823228969E-2, - 4.88280347570998239232E-2, - 7.03124996963957469739E-2, - 1.24999999999870820058E-1, - 4.99999999999999999821E-1 -}; -static double C1 = 1.3862943611198906188E0; /* log(4) */ -#endif - -#ifdef ANSIPROT -extern double polevl ( double, void *, int ); -extern double p1evl ( double, void *, int ); -extern double log ( double ); -#else -double polevl(), p1evl(), log(); -#endif -extern double MACHEP, MAXNUM; - -double ellpk(x) -double x; -{ - -if( (x < 0.0) || (x > 1.0) ) - { - mtherr( "ellpk", DOMAIN ); - return( 0.0 ); - } - -if( x > MACHEP ) - { - return( polevl(x,P,10) - log(x) * polevl(x,Q,10) ); - } -else - { - if( x == 0.0 ) - { - mtherr( "ellpk", SING ); - return( MAXNUM ); - } - else - { - return( C1 - 0.5 * log(x) ); - } - } -} -- cgit v1.2.3