1 /* @(#)log.c 4.2 01/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 } 65