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