
/* log base 10 x */

#include <math.h>

#define N_COEFF 5

static double _log_coeff[] = 
	{
	0.868591718,
	0.289335524,
	0.177522071,
	0.094375476,
	0.191337714
	};

#define SQRT_10		3.16227766
#define INV_SQRT_10	0.316227766

double log(x)
double x;
{
  double x0;			/* x folded between 1/sqr(10) and sqr(10) */
  double x1;			/* (x0 - 1) / (x0 + 1) */
  double x2;			/* x1 ^ 2 */
  double xn;			/* x0 ^ 2n+1 */

  double result, offset;
  int i;

  if (x == 1.0)
	return(0.0);
  if (x < 0.0)
  	/* error */
	return(0);		/* need NaN here */
	  
  offset = 0.0;
  x0 = x;
  while (x0 > SQRT_10)
  	{
	offset += 1.0;
	x0 /= 10.0;
	}
  while (x0 < INV_SQRT_10)
  	{
	offset -= 1.0;
	x0 *= 10.0;
	}
  xn = x1 = (x0 - 1.0) / (x0 + 1.0);
  x2 = x1 * x1;
  result = 0.0;
  for (i = 0 ; i < N_COEFF ; i++)
	{
	result += (xn * _log_coeff[i]);
	xn = xn * x2;
	}
  return(result + offset);
}
