Annotation of 43BSD/usr.lib/libom/log.c, revision 1.1.1.1

1.1       root        1: /*     @(#)log.c       4.2     1/22/85 */
                      2: 
                      3: /*
                      4:        log returns the natural logarithm of its floating
                      5:        point argument.
                      6: 
                      7:        The coefficients are #2705 from Hart & Cheney. (19.38D)
                      8: 
                      9:        It calls frexp.
                     10: */
                     11: 
                     12: #include <errno.h>
                     13: #include <math.h>
                     14: 
                     15: int    errno;
                     16: double frexp();
                     17: static double  log2    = 0.693147180559945309e0;
                     18: static double  ln10    = 2.302585092994045684;
                     19: static double  sqrto2  = 0.707106781186547524e0;
                     20: static double  p0      = -.240139179559210510e2;
                     21: static double  p1      = 0.309572928215376501e2;
                     22: static double  p2      = -.963769093377840513e1;
                     23: static double  p3      = 0.421087371217979714e0;
                     24: static double  q0      = -.120069589779605255e2;
                     25: static double  q1      = 0.194809660700889731e2;
                     26: static double  q2      = -.891110902798312337e1;
                     27: 
                     28: double
                     29: log(arg)
                     30: double arg;
                     31: {
                     32:        double x,z, zsq, temp;
                     33:        int exp;
                     34: 
                     35:        if(arg <= 0.) {
                     36:                errno = EDOM;
                     37:                return(-HUGE);
                     38:        }
                     39:        x = frexp(arg,&exp);
                     40:        while(x<0.5) {
                     41:                x = x*2;
                     42:                exp = exp-1;
                     43:        }
                     44:        if(x<sqrto2) {
                     45:                x = 2*x;
                     46:                exp = exp-1;
                     47:        }
                     48: 
                     49:        z = (x-1)/(x+1);
                     50:        zsq = z*z;
                     51: 
                     52:        temp = ((p3*zsq + p2)*zsq + p1)*zsq + p0;
                     53:        temp = temp/(((1.0*zsq + q2)*zsq + q1)*zsq + q0);
                     54:        temp = temp*z + exp*log2;
                     55:        return(temp);
                     56: }
                     57: 
                     58: double
                     59: log10(arg)
                     60: double arg;
                     61: {
                     62: 
                     63:        return(log(arg)/ln10);
                     64: }

unix.superglobalmegacorp.com

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