/* log - returns natural log and log base 10 */

/*
Based on Yale CS Dept C library.

Coefficients are #2705 from Hart & Cheney. (19.38D)

It calls frexp and itof.
*/

extern	double	frexp();
extern	double	itof();

static	double	log2	= 0.693147180559945309e0;
static	double	ln10	= 2.302585092994045684;
static	double	sqrto2	= 0.707106781186547524e0;
static	double	p0	= -.240139179559210510e2;
static	double	p1	= 0.309572928215376501e2;
static	double	p2	= -.963769093368686593e1;
static	double	p3	= 0.421087371217979714e0;
static	double	q0	= -.120069589779605255e2;
static	double	q1	= 0.194809660700889731e2;
static	double	q2	= -.891110902798312337e1;

double	log(arg)

double	arg;
{
	double	x,z,zsq,temp;
	int	exp;

	if (arg <= 0.0)
	    return(-.99e30);
	x = frexp(arg,&exp);
	while (x < sqrto2)
	{
	    x = x*2.0;
	    exp--;
	}
	z = (x-1.0)/(x+1.0);
	zsq = z*z;
	temp = ((p3*zsq+p2)*zsq+p1)*zsq+p0;
	temp = temp/(((zsq+q2)*zsq+q1)*zsq+q0);
	temp = temp*z+itof(exp)*log2;
	return(temp);
}

double	log10(arg)

double	arg;
{
	return(log(arg)/ln10);
}