libm/float/exp2f.c

   1 /*                                                      exp2f.c
   2  *
   3  *      Base 2 exponential function
   4  *
   5  *
   6  *
   7  * SYNOPSIS:
   8  *
   9  * float x, y, exp2f();
  10  *
  11  * y = exp2f( x );
  12  *
  13  *
  14  *
  15  * DESCRIPTION:
  16  *
  17  * Returns 2 raised to the x power.
  18  *
  19  * Range reduction is accomplished by separating the argument
  20  * into an integer k and fraction f such that
  21  *     x    k  f
  22  *    2  = 2  2.
  23  *
  24  * A polynomial approximates 2**x in the basic range [-0.5, 0.5].
  25  *
  26  *
  27  * ACCURACY:
  28  *
  29  *                      Relative error:
  30  * arithmetic   domain     # trials      peak         rms
  31  *    IEEE     -127,+127    100000      1.7e-7      2.8e-8
  32  *
  33  *
  34  * See exp.c for comments on error amplification.
  35  *
  36  *
  37  * ERROR MESSAGES:
  38  *
  39  *   message         condition      value returned
  40  * exp underflow    x < -MAXL2        0.0
  41  * exp overflow     x > MAXL2         MAXNUMF
  42  *
  43  * For IEEE arithmetic, MAXL2 = 127.
  44  */
  45 \f
  46
  47 /*
  48 Cephes Math Library Release 2.2:  June, 1992
  49 Copyright 1984, 1987, 1988, 1992 by Stephen L. Moshier
  50 Direct inquiries to 30 Frost Street, Cambridge, MA 02140
  51 */
  52
  53
  54
  55 #include <math.h>
  56 static char fname[] = {"exp2f"};
  57
  58 static float P[] = {
  59  1.535336188319500E-004,
  60  1.339887440266574E-003,
  61  9.618437357674640E-003,
  62  5.550332471162809E-002,
  63  2.402264791363012E-001,
  64  6.931472028550421E-001
  65 };
  66 #define MAXL2 127.0
  67 #define MINL2 -127.0
  68
  69
  70
  71 extern float MAXNUMF;
  72
  73 float polevlf(float, float *, int), floorf(float), ldexpf(float, int);
  74
  75 float exp2f( float xx )
  76 {
  77 float x, px;
  78 int i0;
  79
  80 x = xx;
  81 if( x > MAXL2)
  82         {
  83         mtherr( fname, OVERFLOW );
  84         return( MAXNUMF );
  85         }
  86
  87 if( x < MINL2 )
  88         {
  89         mtherr( fname, UNDERFLOW );
  90         return(0.0);
  91         }
  92
  93 /* The following is necessary because range reduction blows up: */
  94 if( x == 0 )
  95         return(1.0);
  96
  97 /* separate into integer and fractional parts */
  98 px = floorf(x);
  99 i0 = px;
 100 x = x - px;
 101
 102 if( x > 0.5 )
 103         {
 104         i0 += 1;
 105         x -= 1.0;
 106         }
 107
 108 /* rational approximation
 109  * exp2(x) = 1.0 +  xP(x)
 110  */
 111 px = 1.0 + x * polevlf( x, P, 5 );
 112
 113 /* scale by power of 2 */
 114 px = ldexpf( px, i0 );
 115 return(px);
 116 }