#define __normdf norm	/* for now */
/*
 * double frexp(value, eptr)
 * double value;
 * int *eptr;
 *
 * returns significand (|significand| < 1)
 *	   in *eptr returns n such that value = significand * 2**n
 *
 * ++jrb	bammi@dsrgsun.ces.cwru.edu
 *
 */
#include "flonum.h"

#define BIAS 1023
#define B1   1022

double frexp(value, eptr)
double value;
#ifdef SHORTLIB
short *eptr;
#else
int *eptr;
#endif
{
    struct bitdouble *res = (struct bitdouble *) &value;
    unsigned int expo, sign;
#ifdef __STDC__
    double __normdf( double, int, int, int );
#else
    extern double __normdf();
#endif
    
    expo = res->exp;
    sign = res->sign;
    res->exp = 0;
    res->sign = 0;
    *eptr = expo - B1;
    if(expo != 0)
	res->exp = 1;	/* put in hidden bit */
    else
	if((res->mant1 == 0) && (res->mant2 == 0))
	{
	    *eptr = 0;
	    return 0.0;	/* no point in normalizing (exp will be wrong) */
	}

    return __normdf(value, B1, sign, 0);
}

#ifdef TEST /* combined test for frexp and ldexp */
	/* 32 bit ints for this test please */

#ifdef __MSHORT__
#error "please run this test in 32 bit int mode"
#endif

#include <stdio.h>

#define NTEST	100000L
#define MAXRAND 0x7fffffffL	/* assuming 32 bit ints */
#define ABS(X)	( (X < 0) ? (-X) : X )
extern long rand();

int main()
{
    double sig, e, expected, result, max_abs_err;
    int twoexp;
    register long i;
    register long errs = 0;
    register int s;
#ifdef __STDC__
    double ldexp(double, int);
    double frexp(double, int *);
#else
    extern double ldexp();
    extern double frexp();
#endif

    max_abs_err = 0.0;
    for(i = 0; i < NTEST; i++)
    {
	expected = ((double)(s = rand()) + (double)rand())/(double)(MAXRAND);
	if(s > (MAXRAND >> 1)) expected = -expected;

	sig = frexp(expected, &twoexp);
	if(ABS(sig) >= 1.0)
	{
	    printf("sig > 1, %.4e: %.4e  %d\n", expected, sig, twoexp);
	    errs++;
	    continue;
	}
	
	result = ldexp(sig, twoexp);
	
	e = (expected == 0.0) ? result : (result - expected)/expected;
	if(e < 0) e = (-e);
	if(e > 1.0e-6)
	{
	    printf("%.8e: %.8e  E %.8e\n",
		   expected, result, e);
	    errs++;
	}
	if (e > max_abs_err) max_abs_err = e;
    }
    
    printf("%ld Error(s), Max abs err %.14e\n", errs, max_abs_err);
    return errs;
}
#endif /* TEST */