Annotation of 42BSD/usr.lib/libm/exp.c, revision 1.1.1.1

1.1       root        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: }

unix.superglobalmegacorp.com

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.