4.3BSD - /usr/src/usr.lib/libom/exp.c

   1: /*	@(#)exp.c	4.1	12/25/82	*/
   2: 
   3: /*
   4: 	exp returns the exponential function of its
   5: 	floating-point argument.
   6: 
   7: 	The coefficients are #1069 from Hart and Cheney. (22.35D)
   8: */
   9: 
  10: #include <errno.h>
  11: #include <math.h>
  12: 
  13: int errno;
  14: static double   p0  = .2080384346694663001443843411e7;
  15: static double   p1  = .3028697169744036299076048876e5;
  16: static double   p2  = .6061485330061080841615584556e2;
  17: static double   q0  = .6002720360238832528230907598e7;
  18: static double   q1  = .3277251518082914423057964422e6;
  19: static double   q2  = .1749287689093076403844945335e4;
  20: static double   log2e   = 1.4426950408889634073599247;
  21: static double   sqrt2   = 1.4142135623730950488016887;
  22: static double   maxf    = 10000;
  23: 
  24: double
  25: exp(arg)
  26: double arg;
  27: {
  28:     double fract;
  29:     double temp1, temp2, xsq;
  30:     int ent;
  31: 
  32:     if(arg == 0.)
  33:         return(1.);
  34:     if(arg < -maxf)
  35:         return(0.);
  36:     if(arg > maxf) {
  37:         errno = ERANGE;
  38:         return(HUGE);
  39:     }
  40:     arg *= log2e;
  41:     ent = floor(arg);
  42:     fract = (arg-ent) - 0.5;
  43:     xsq = fract*fract;
  44:     temp1 = ((p2*xsq+p1)*xsq+p0)*fract;
  45:     temp2 = ((1.0*xsq+q2)*xsq+q1)*xsq + q0;
  46:     return(ldexp(sqrt2*(temp2+temp1)/(temp2-temp1), ent));
  47: }

Defined functions

Defined variables