log10.c

OMAP1030 处理器的ARM 侧硬件测试代码 OMAP1030 是TI的双核处理器

/****************************************************************************/
/*  log10  v2.24                                                           */
/*  Copyright (c) 1995-2002 Texas Instruments Incorporated                  */
/****************************************************************************/
#include <math.h>
#include <values.h>
#include <errno.h>

/****************************************************************************/
/*  LOG10() - Logarithm							    */
/*									    */
/*  result = log10(e) * log(x)						    */
/****************************************************************************/
double log10(double x)
{
    double a, b, f, r, w, z, znum;
    int n;

    /*********************************************************************/
    /* check for domain and range errors                                 */
    /*********************************************************************/
    if (x <= 0) { errno = (x == 0) ? ERANGE : EDOM; return (-HUGE_VAL); }

    /*********************************************************************/
    /* f = mantissa(x), n = exponent(x)                                  */
    /*********************************************************************/
    f = frexp(x, &n);

    /*********************************************************************/
    /* for numbers < sqrt(0.5)                                           */
    /*********************************************************************/
    if (f <= SQRTHALF)
    {
	--n;
	znum = f - 0.5;
	z    = znum / (znum * 0.5 + 0.5);
    }

    /*********************************************************************/
    /* for numbers > sqrt(0.5)                                           */
    /*********************************************************************/
    else
    {
	znum = (f - 0.5) - 0.5;
	z    = znum / (f * 0.5 + 0.5);
    }

    /*********************************************************************/
    /* determine polynomial expression                                   */
    /*********************************************************************/
    w = z * z;

#if BITS<=24
    a = A0;
    b = w + B0;
#elif BITS>=25 && BITS<=32
    a = A1 * w + A0;
    b = w + B0;
#elif BITS>=33 && BITS<=48
    a = (A2 * w + A1) * w + A0;
    b = (w + B1) * w + B0;
#else
    a = (A2 * w + A1) * w + A0;
    b = ((w + B2) * w + B1) * w + B0;
#endif

    /*********************************************************************/
    /* calculate the natural log of (mant x) / 2		         */
    /*********************************************************************/
    r = z + z * w * (a / b);

    /*********************************************************************/
    /* log10(x) = (ln (mant x) + 2 * (exp x)) * log10(e) (but more       */
    /* mathematically stable)	                                         */
    /*********************************************************************/
    return (((n * C4 + r) + n * C3) * LOG10e);
}