summaryrefslogtreecommitdiff
path: root/libm/double/exp2.c
diff options
context:
space:
mode:
authorEric Andersen <andersen@codepoet.org>2001-11-22 14:04:29 +0000
committerEric Andersen <andersen@codepoet.org>2001-11-22 14:04:29 +0000
commit7ce331c01ce6eb7b3f5c715a38a24359da9c6ee2 (patch)
tree3a7e8476e868ae15f4da1b7ce26b2db6f434468c /libm/double/exp2.c
parentc117dd5fb183afb1a4790a6f6110d88704be6bf8 (diff)
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
Diffstat (limited to 'libm/double/exp2.c')
-rw-r--r--libm/double/exp2.c183
1 files changed, 0 insertions, 183 deletions
diff --git a/libm/double/exp2.c b/libm/double/exp2.c
deleted file mode 100644
index be5bdfd0c..000000000
--- a/libm/double/exp2.c
+++ /dev/null
@@ -1,183 +0,0 @@
-/* exp2.c
- *
- * Base 2 exponential function
- *
- *
- *
- * SYNOPSIS:
- *
- * double x, y, exp2();
- *
- * y = exp2( x );
- *
- *
- *
- * DESCRIPTION:
- *
- * Returns 2 raised to the x power.
- *
- * Range reduction is accomplished by separating the argument
- * into an integer k and fraction f such that
- * x k f
- * 2 = 2 2.
- *
- * A Pade' form
- *
- * 1 + 2x P(x**2) / (Q(x**2) - x P(x**2) )
- *
- * approximates 2**x in the basic range [-0.5, 0.5].
- *
- *
- * ACCURACY:
- *
- * Relative error:
- * arithmetic domain # trials peak rms
- * IEEE -1022,+1024 30000 1.8e-16 5.4e-17
- *
- *
- * See exp.c for comments on error amplification.
- *
- *
- * ERROR MESSAGES:
- *
- * message condition value returned
- * exp underflow x < -MAXL2 0.0
- * exp overflow x > MAXL2 MAXNUM
- *
- * For DEC arithmetic, MAXL2 = 127.
- * For IEEE arithmetic, MAXL2 = 1024.
- */
-
-
-/*
-Cephes Math Library Release 2.8: June, 2000
-Copyright 1984, 1995, 2000 by Stephen L. Moshier
-*/
-
-
-
-#include <math.h>
-
-#ifdef UNK
-static double P[] = {
- 2.30933477057345225087E-2,
- 2.02020656693165307700E1,
- 1.51390680115615096133E3,
-};
-static double Q[] = {
-/* 1.00000000000000000000E0,*/
- 2.33184211722314911771E2,
- 4.36821166879210612817E3,
-};
-#define MAXL2 1024.0
-#define MINL2 -1024.0
-#endif
-
-#ifdef DEC
-static unsigned short P[] = {
-0036675,0027102,0122327,0053227,
-0041241,0116724,0115412,0157355,
-0042675,0036404,0101733,0132226,
-};
-static unsigned short Q[] = {
-/*0040200,0000000,0000000,0000000,*/
-0042151,0027450,0077732,0160744,
-0043210,0100661,0077550,0056560,
-};
-#define MAXL2 127.0
-#define MINL2 -127.0
-#endif
-
-#ifdef IBMPC
-static unsigned short P[] = {
-0xead3,0x549a,0xa5c8,0x3f97,
-0x5bde,0x9361,0x33ba,0x4034,
-0x7693,0x907b,0xa7a0,0x4097,
-};
-static unsigned short Q[] = {
-/*0x0000,0x0000,0x0000,0x3ff0,*/
-0x5c3c,0x0ffb,0x25e5,0x406d,
-0x0bae,0x2fed,0x1036,0x40b1,
-};
-#define MAXL2 1024.0
-#define MINL2 -1022.0
-#endif
-
-#ifdef MIEEE
-static unsigned short P[] = {
-0x3f97,0xa5c8,0x549a,0xead3,
-0x4034,0x33ba,0x9361,0x5bde,
-0x4097,0xa7a0,0x907b,0x7693,
-};
-static unsigned short Q[] = {
-/*0x3ff0,0x0000,0x0000,0x0000,*/
-0x406d,0x25e5,0x0ffb,0x5c3c,
-0x40b1,0x1036,0x2fed,0x0bae,
-};
-#define MAXL2 1024.0
-#define MINL2 -1022.0
-#endif
-
-#ifdef ANSIPROT
-extern double polevl ( double, void *, int );
-extern double p1evl ( double, void *, int );
-extern double floor ( double );
-extern double ldexp ( double, int );
-extern int isnan ( double );
-extern int isfinite ( double );
-#else
-double polevl(), p1evl(), floor(), ldexp();
-int isnan(), isfinite();
-#endif
-#ifdef INFINITIES
-extern double INFINITY;
-#endif
-extern double MAXNUM;
-
-double exp2(x)
-double x;
-{
-double px, xx;
-short n;
-
-#ifdef NANS
-if( isnan(x) )
- return(x);
-#endif
-if( x > MAXL2)
- {
-#ifdef INFINITIES
- return( INFINITY );
-#else
- mtherr( "exp2", OVERFLOW );
- return( MAXNUM );
-#endif
- }
-
-if( x < MINL2 )
- {
-#ifndef INFINITIES
- mtherr( "exp2", UNDERFLOW );
-#endif
- return(0.0);
- }
-
-xx = x; /* save x */
-/* separate into integer and fractional parts */
-px = floor(x+0.5);
-n = px;
-x = x - px;
-
-/* rational approximation
- * exp2(x) = 1 + 2xP(xx)/(Q(xx) - P(xx))
- * where xx = x**2
- */
-xx = x * x;
-px = x * polevl( xx, P, 2 );
-x = px / ( p1evl( xx, Q, 2 ) - px );
-x = 1.0 + ldexp( x, 1 );
-
-/* scale by power of 2 */
-x = ldexp( x, n );
-return(x);
-}