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

1 files changed, 0 insertions, 250 deletions

diff --git a/libm/double/log10.c b/libm/double/log10.c
deleted file mode 100644
index 7dc72e253..000000000
--- a/libm/double/log10.c
+++ /dev/null
@@ -1,250 +0,0 @@
-/*							log10.c
- *
- *	Common logarithm
- *
- *
- *
- * SYNOPSIS:
- *
- * double x, y, log10();
- *
- * y = log10( x );
- *
- *
- *
- * DESCRIPTION:
- *
- * Returns logarithm to the base 10 of x.
- *
- * The argument is separated into its exponent and fractional
- * parts.  The logarithm of the fraction is approximated by
- *
- *     log(1+x) = x - 0.5 x**2 + x**3 P(x)/Q(x).
- *
- *
- *
- * ACCURACY:
- *
- *                      Relative error:
- * arithmetic   domain     # trials      peak         rms
- *    IEEE      0.5, 2.0     30000      1.5e-16     5.0e-17
- *    IEEE      0, MAXNUM    30000      1.4e-16     4.8e-17
- *    DEC       1, MAXNUM    50000      2.5e-17     6.0e-18
- *
- * In the tests over the interval [1, MAXNUM], the logarithms
- * of the random arguments were uniformly distributed over
- * [0, MAXLOG].
- *
- * ERROR MESSAGES:
- *
- * log10 singularity:  x = 0; returns -INFINITY
- * log10 domain:       x < 0; returns NAN
- */
-
-/*
-Cephes Math Library Release 2.8:  June, 2000
-Copyright 1984, 1995, 2000 by Stephen L. Moshier
-*/
-
-#include <math.h>
-static char fname[] = {"log10"};
-
-/* Coefficients for log(1+x) = x - x**2/2 + x**3 P(x)/Q(x)
- * 1/sqrt(2) <= x < sqrt(2)
- */
-#ifdef UNK
-static double P[] = {
-  4.58482948458143443514E-5,
-  4.98531067254050724270E-1,
-  6.56312093769992875930E0,
-  2.97877425097986925891E1,
-  6.06127134467767258030E1,
-  5.67349287391754285487E1,
-  1.98892446572874072159E1
-};
-static double Q[] = {
-/* 1.00000000000000000000E0, */
-  1.50314182634250003249E1,
-  8.27410449222435217021E1,
-  2.20664384982121929218E2,
-  3.07254189979530058263E2,
-  2.14955586696422947765E2,
-  5.96677339718622216300E1
-};
-#endif
-
-#ifdef DEC
-static unsigned short P[] = {
-0034500,0046473,0051374,0135174,
-0037777,0037566,0145712,0150321,
-0040722,0002426,0031543,0123107,
-0041356,0046513,0170752,0004346,
-0041562,0071553,0023536,0163343,
-0041542,0170221,0024316,0114216,
-0041237,0016454,0046611,0104602
-};
-static unsigned short Q[] = {
-/*0040200,0000000,0000000,0000000,*/
-0041160,0100260,0067736,0102424,
-0041645,0075552,0036563,0147072,
-0042134,0125025,0021132,0025320,
-0042231,0120211,0046030,0103271,
-0042126,0172241,0052151,0120426,
-0041556,0125702,0072116,0047103
-};
-#endif
-
-#ifdef IBMPC
-static unsigned short P[] = {
-0x974f,0x6a5f,0x09a7,0x3f08,
-0x5a1a,0xd979,0xe7ee,0x3fdf,
-0x74c9,0xc66c,0x40a2,0x401a,
-0x411d,0x7e3d,0xc9a9,0x403d,
-0xdcdc,0x64eb,0x4e6d,0x404e,
-0xd312,0x2519,0x5e12,0x404c,
-0x3130,0x89b1,0xe3a5,0x4033
-};
-static unsigned short Q[] = {
-/*0x0000,0x0000,0x0000,0x3ff0,*/
-0xd0a2,0x0dfb,0x1016,0x402e,
-0x79c7,0x47ae,0xaf6d,0x4054,
-0x455a,0xa44b,0x9542,0x406b,
-0x10d7,0x2983,0x3411,0x4073,
-0x3423,0x2a8d,0xde94,0x406a,
-0xc9c8,0x4e89,0xd578,0x404d
-};
-#endif
-
-#ifdef MIEEE
-static unsigned short P[] = {
-0x3f08,0x09a7,0x6a5f,0x974f,
-0x3fdf,0xe7ee,0xd979,0x5a1a,
-0x401a,0x40a2,0xc66c,0x74c9,
-0x403d,0xc9a9,0x7e3d,0x411d,
-0x404e,0x4e6d,0x64eb,0xdcdc,
-0x404c,0x5e12,0x2519,0xd312,
-0x4033,0xe3a5,0x89b1,0x3130
-};
-static unsigned short Q[] = {
-0x402e,0x1016,0x0dfb,0xd0a2,
-0x4054,0xaf6d,0x47ae,0x79c7,
-0x406b,0x9542,0xa44b,0x455a,
-0x4073,0x3411,0x2983,0x10d7,
-0x406a,0xde94,0x2a8d,0x3423,
-0x404d,0xd578,0x4e89,0xc9c8
-};
-#endif
-
-#define SQRTH 0.70710678118654752440
-#define L102A 3.0078125E-1
-#define L102B 2.48745663981195213739E-4
-#define L10EA 4.3359375E-1
-#define L10EB 7.00731903251827651129E-4
-
-#ifdef ANSIPROT
-extern double frexp ( double, int * );
-extern double ldexp ( double, int );
-extern double polevl ( double, void *, int );
-extern double p1evl ( double, void *, int );
-extern int isnan ( double );
-extern int isfinite ( double );
-#else
-double frexp(), ldexp(), polevl(), p1evl();
-int isnan(), isfinite();
-#endif
-extern double LOGE2, SQRT2, INFINITY, NAN;
-
-double log10(x)
-double x;
-{
-VOLATILE double z;
-double y;
-#ifdef DEC
-short *q;
-#endif
-int e;
-
-#ifdef NANS
-if( isnan(x) )
-	return(x);
-#endif
-#ifdef INFINITIES
-if( x == INFINITY )
-	return(x);
-#endif
-/* Test for domain */
-if( x <= 0.0 )
-	{
-	if( x == 0.0 )
-	        {
-		mtherr( fname, SING );
-		return( -INFINITY );
-	        }
-	else
-	        {
-		mtherr( fname, DOMAIN );
-		return( NAN );
-	        }
-	}
-
-/* separate mantissa from exponent */
-
-#ifdef DEC
-q = (short *)&x;
-e = *q;			/* short containing exponent */
-e = ((e >> 7) & 0377) - 0200;	/* the exponent */
-*q &= 0177;	/* strip exponent from x */
-*q |= 040000;	/* x now between 0.5 and 1 */
-#endif
-
-#ifdef IBMPC
-x = frexp( x, &e );
-/*
-q = (short *)&x;
-q += 3;
-e = *q;
-e = ((e >> 4) & 0x0fff) - 0x3fe;
-*q &= 0x0f;
-*q |= 0x3fe0;
-*/
-#endif
-
-/* Equivalent C language standard library function: */
-#ifdef UNK
-x = frexp( x, &e );
-#endif
-
-#ifdef MIEEE
-x = frexp( x, &e );
-#endif
-
-/* logarithm using log(1+x) = x - .5x**2 + x**3 P(x)/Q(x) */
-
-if( x < SQRTH )
-	{
-	e -= 1;
-	x = ldexp( x, 1 ) - 1.0; /*  2x - 1  */
-	}	
-else
-	{
-	x = x - 1.0;
-	}
-
-
-/* rational form */
-z = x*x;
-y = x * ( z * polevl( x, P, 6 ) / p1evl( x, Q, 6 ) );
-y = y - ldexp( z, -1 );   /*  y - 0.5 * x**2  */
-
-/* multiply log of fraction by log10(e)
- * and base 2 exponent by log10(2)
- */
-z = (x + y) * L10EB;  /* accumulate terms in order of size */
-z += y * L10EA;
-z += x * L10EA;
-z += e * L102B;
-z += e * L102A;
-
-
-return( z );
-}