diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/math/modf.c | 89 | ||||
| -rw-r--r-- | src/math/modff.c | 28 | 
2 files changed, 45 insertions, 72 deletions
| diff --git a/src/math/modf.c b/src/math/modf.c index cca3b652..de45069f 100644 --- a/src/math/modf.c +++ b/src/math/modf.c @@ -1,68 +1,37 @@ -/* origin: FreeBSD /usr/src/lib/msun/src/s_modf.c */ -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ -/* - * modf(double x, double *iptr) - * return fraction part of x, and return x's integral part in *iptr. - * Method: - *      Bit twiddling. - * - * Exception: - *      No exception. - */ - -#include "libm.h" +#include <math.h> +#include <stdint.h>  double modf(double x, double *iptr)  { -	int32_t i0,i1,j0; -	uint32_t i; +	union {double x; uint64_t n;} u = {x}; +	uint64_t mask; +	int e; -	EXTRACT_WORDS(i0, i1, x); -	j0 = ((i0>>20) & 0x7ff) - 0x3ff; /* exponent of x */ -	if (j0 < 20) {  /* integer part in high x */ -		if (j0 < 0) {  /* |x| < 1 */ -			INSERT_WORDS(*iptr, i0 & 0x80000000, 0); /* *iptr = +-0 */ -			return x; -		} -		i = 0x000fffff >> j0; -		if (((i0&i)|i1) == 0) {  /* x is integral */ -			uint32_t high; -			*iptr = x; -			GET_HIGH_WORD(high, x); -			INSERT_WORDS(x, high & 0x80000000, 0);  /* return +-0 */ -			return x; -		} -		INSERT_WORDS(*iptr, i0&~i, 0); -		return x - *iptr; -	} else if (j0 > 51) {  /* no fraction part */ -		uint32_t high; -		if (j0 == 0x400) {  /* inf/NaN */ -			*iptr = x; -			return 0.0 / x; -		} +	e = (int)(u.n>>52 & 0x7ff) - 0x3ff; + +	/* no fractional part */ +	if (e >= 52) {  		*iptr = x; -		GET_HIGH_WORD(high, x); -		INSERT_WORDS(x, high & 0x80000000, 0);  /* return +-0 */ -		return x; -	} else {               /* fraction part in low x */ -		i = (uint32_t)0xffffffff >> (j0 - 20); -		if ((i1&i) == 0) {  /* x is integral */ -			uint32_t high; -			*iptr = x; -			GET_HIGH_WORD(high, x); -			INSERT_WORDS(x, high & 0x80000000, 0);  /* return +-0 */ +		if (e == 0x400 && u.n<<12 != 0) /* nan */  			return x; -		} -		INSERT_WORDS(*iptr, i0, i1&~i); -		return x - *iptr; +		u.n &= (uint64_t)1<<63; +		return u.x; +	} + +	/* no integral part*/ +	if (e < 0) { +		u.n &= (uint64_t)1<<63; +		*iptr = u.x; +		return x; +	} + +	mask = (uint64_t)-1>>12 >> e; +	if ((u.n & mask) == 0) { +		*iptr = x; +		u.n &= (uint64_t)1<<63; +		return u.x;  	} +	u.n &= ~mask; +	*iptr = u.x; +	return x - *iptr;  } diff --git a/src/math/modff.c b/src/math/modff.c index bf6e4ced..84d0b82a 100644 --- a/src/math/modff.c +++ b/src/math/modff.c @@ -1,33 +1,37 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h>  float modff(float x, float *iptr)  { -	uint32_t u, mask; +	union {float x; uint32_t n;} u = {x}; +	uint32_t mask;  	int e; -	GET_FLOAT_WORD(u, x); -	e = (int)(u>>23 & 0xff) - 0x7f; +	e = (int)(u.n>>23 & 0xff) - 0x7f;  	/* no fractional part */  	if (e >= 23) {  		*iptr = x; -		if (e == 0x80 && u<<9 != 0) /* nan */ +		if (e == 0x80 && u.n<<9 != 0) { /* nan */  			return x; -		SET_FLOAT_WORD(x, u & 0x80000000); -		return x; +		} +		u.n &= 0x80000000; +		return u.x;  	}  	/* no integral part */  	if (e < 0) { -		SET_FLOAT_WORD(*iptr, u & 0x80000000); +		u.n &= 0x80000000; +		*iptr = u.x;  		return x;  	}  	mask = 0x007fffff>>e; -	if ((u & mask) == 0) { +	if ((u.n & mask) == 0) {  		*iptr = x; -		SET_FLOAT_WORD(x, u & 0x80000000); -		return x; +		u.n &= 0x80000000; +		return u.x;  	} -	SET_FLOAT_WORD(*iptr, u & ~mask); +	u.n &= ~mask; +	*iptr = u.x;  	return x - *iptr;  } | 
