diff options
| -rw-r--r-- | src/internal/libm.h | 122 | ||||
| -rw-r--r-- | src/math/__fpclassify.c | 11 | ||||
| -rw-r--r-- | src/math/__fpclassifyf.c | 11 | ||||
| -rw-r--r-- | src/math/copysign.c | 11 | ||||
| -rw-r--r-- | src/math/copysignf.c | 17 | ||||
| -rw-r--r-- | src/math/fabs.c | 11 | ||||
| -rw-r--r-- | src/math/fabsf.c | 11 | ||||
| -rw-r--r-- | src/math/nextafter.c | 30 | ||||
| -rw-r--r-- | src/math/nextafterf.c | 30 | ||||
| -rw-r--r-- | src/math/nexttoward.c | 27 | ||||
| -rw-r--r-- | src/math/nexttowardf.c | 25 | 
11 files changed, 140 insertions, 166 deletions
| diff --git a/src/internal/libm.h b/src/internal/libm.h index 99448a08..946c310d 100644 --- a/src/internal/libm.h +++ b/src/internal/libm.h @@ -21,16 +21,6 @@  #include "libc.h" -union fshape { -	float value; -	uint32_t bits; -}; - -union dshape { -	double value; -	uint64_t bits; -}; -  #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024  #elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384 && __BYTE_ORDER == __LITTLE_ENDIAN  union ldshape { @@ -58,86 +48,86 @@ union ldshape {  #error Unsupported long double representation  #endif -#define FORCE_EVAL(x) do {                          \ -	if (sizeof(x) == sizeof(float)) {           \ -		volatile float __x;                 \ -		__x = (x);                          \ -	} else if (sizeof(x) == sizeof(double)) {   \ -		volatile double __x;                \ -		__x = (x);                          \ -	} else {                                    \ -		volatile long double __x;           \ -		__x = (x);                          \ -	}                                           \ +#define FORCE_EVAL(x) do {                        \ +	if (sizeof(x) == sizeof(float)) {         \ +		volatile float __x;               \ +		__x = (x);                        \ +	} else if (sizeof(x) == sizeof(double)) { \ +		volatile double __x;              \ +		__x = (x);                        \ +	} else {                                  \ +		volatile long double __x;         \ +		__x = (x);                        \ +	}                                         \  } while(0)  /* Get two 32 bit ints from a double.  */ -#define EXTRACT_WORDS(hi,lo,d)                                  \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.value = (d);                                              \ -  (hi) = __u.bits >> 32;                                        \ -  (lo) = (uint32_t)__u.bits;                                    \ +#define EXTRACT_WORDS(hi,lo,d)                    \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.f = (d);                                    \ +  (hi) = __u.i >> 32;                             \ +  (lo) = (uint32_t)__u.i;                         \  } while (0)  /* Get the more significant 32 bit int from a double.  */ -#define GET_HIGH_WORD(i,d)                                      \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.value = (d);                                              \ -  (i) = __u.bits >> 32;                                         \ +#define GET_HIGH_WORD(hi,d)                       \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.f = (d);                                    \ +  (hi) = __u.i >> 32;                             \  } while (0)  /* Get the less significant 32 bit int from a double.  */ -#define GET_LOW_WORD(i,d)                                       \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.value = (d);                                              \ -  (i) = (uint32_t)__u.bits;                                     \ +#define GET_LOW_WORD(lo,d)                        \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.f = (d);                                    \ +  (lo) = (uint32_t)__u.i;                         \  } while (0)  /* Set a double from two 32 bit ints.  */ -#define INSERT_WORDS(d,hi,lo)                                   \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.bits = ((uint64_t)(hi) << 32) | (uint32_t)(lo);           \ -  (d) = __u.value;                                              \ +#define INSERT_WORDS(d,hi,lo)                     \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.i = ((uint64_t)(hi)<<32) | (uint32_t)(lo);  \ +  (d) = __u.f;                                    \  } while (0)  /* Set the more significant 32 bits of a double from an int.  */ -#define SET_HIGH_WORD(d,hi)                                     \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.value = (d);                                              \ -  __u.bits &= 0xffffffff;                                       \ -  __u.bits |= (uint64_t)(hi) << 32;                             \ -  (d) = __u.value;                                              \ +#define SET_HIGH_WORD(d,hi)                       \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.f = (d);                                    \ +  __u.i &= 0xffffffff;                            \ +  __u.i |= (uint64_t)(hi) << 32;                  \ +  (d) = __u.f;                                    \  } while (0)  /* Set the less significant 32 bits of a double from an int.  */ -#define SET_LOW_WORD(d,lo)                                      \ -do {                                                            \ -  union dshape __u;                                             \ -  __u.value = (d);                                              \ -  __u.bits &= 0xffffffff00000000ull;                            \ -  __u.bits |= (uint32_t)(lo);                                   \ -  (d) = __u.value;                                              \ +#define SET_LOW_WORD(d,lo)                        \ +do {                                              \ +  union {double f; uint64_t i;} __u;              \ +  __u.f = (d);                                    \ +  __u.i &= 0xffffffff00000000ull;                 \ +  __u.i |= (uint32_t)(lo);                        \ +  (d) = __u.f;                                    \  } while (0)  /* Get a 32 bit int from a float.  */ -#define GET_FLOAT_WORD(i,d)                                     \ -do {                                                            \ -  union fshape __u;                                             \ -  __u.value = (d);                                              \ -  (i) = __u.bits;                                               \ +#define GET_FLOAT_WORD(w,d)                       \ +do {                                              \ +  union {float f; uint32_t i;} __u;               \ +  __u.f = (d);                                    \ +  (w) = __u.i;                                    \  } while (0)  /* Set a float from a 32 bit int.  */ -#define SET_FLOAT_WORD(d,i)                                     \ -do {                                                            \ -  union fshape __u;                                             \ -  __u.bits = (i);                                               \ -  (d) = __u.value;                                              \ +#define SET_FLOAT_WORD(d,w)                       \ +do {                                              \ +  union {float f; uint32_t i;} __u;               \ +  __u.i = (w);                                    \ +  (d) = __u.f;                                    \  } while (0)  /* fdlibm kernel functions */ diff --git a/src/math/__fpclassify.c b/src/math/__fpclassify.c index c9dd0275..f7c0e2df 100644 --- a/src/math/__fpclassify.c +++ b/src/math/__fpclassify.c @@ -1,10 +1,11 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h>  int __fpclassify(double x)  { -	union dshape u = { x }; -	int e = u.bits>>52 & 0x7ff; -	if (!e) return u.bits<<1 ? FP_SUBNORMAL : FP_ZERO; -	if (e==0x7ff) return u.bits<<12 ? FP_NAN : FP_INFINITE; +	union {double f; uint64_t i;} u = {x}; +	int e = u.i>>52 & 0x7ff; +	if (!e) return u.i<<1 ? FP_SUBNORMAL : FP_ZERO; +	if (e==0x7ff) return u.i<<12 ? FP_NAN : FP_INFINITE;  	return FP_NORMAL;  } diff --git a/src/math/__fpclassifyf.c b/src/math/__fpclassifyf.c index 8149087b..fd00eb1b 100644 --- a/src/math/__fpclassifyf.c +++ b/src/math/__fpclassifyf.c @@ -1,10 +1,11 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h>  int __fpclassifyf(float x)  { -	union fshape u = { x }; -	int e = u.bits>>23 & 0xff; -	if (!e) return u.bits<<1 ? FP_SUBNORMAL : FP_ZERO; -	if (e==0xff) return u.bits<<9 ? FP_NAN : FP_INFINITE; +	union {float f; uint32_t i;} u = {x}; +	int e = u.i>>23 & 0xff; +	if (!e) return u.i<<1 ? FP_SUBNORMAL : FP_ZERO; +	if (e==0xff) return u.i<<9 ? FP_NAN : FP_INFINITE;  	return FP_NORMAL;  } diff --git a/src/math/copysign.c b/src/math/copysign.c index 038b8b4c..b09331b6 100644 --- a/src/math/copysign.c +++ b/src/math/copysign.c @@ -1,11 +1,8 @@  #include "libm.h"  double copysign(double x, double y) { -	union dshape ux, uy; - -	ux.value = x; -	uy.value = y; -	ux.bits &= (uint64_t)-1>>1; -	ux.bits |= uy.bits & (uint64_t)1<<63; -	return ux.value; +	union {double f; uint64_t i;} ux={x}, uy={y}; +	ux.i &= -1ULL/2; +	ux.i |= uy.i & 1ULL<<63; +	return ux.f;  } diff --git a/src/math/copysignf.c b/src/math/copysignf.c index 47ab37e4..0af6ae9b 100644 --- a/src/math/copysignf.c +++ b/src/math/copysignf.c @@ -1,11 +1,10 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h> -float copysignf(float x, float y) { -	union fshape ux, uy; - -	ux.value = x; -	uy.value = y; -	ux.bits &= (uint32_t)-1>>1; -	ux.bits |= uy.bits & (uint32_t)1<<31; -	return ux.value; +float copysignf(float x, float y) +{ +	union {float f; uint32_t i;} ux={x}, uy={y}; +	ux.i &= 0x7fffffff; +	ux.i |= uy.i & 0x80000000; +	return ux.f;  } diff --git a/src/math/fabs.c b/src/math/fabs.c index 6e28f1e5..e8258cfd 100644 --- a/src/math/fabs.c +++ b/src/math/fabs.c @@ -1,10 +1,9 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h>  double fabs(double x)  { -	union dshape u; - -	u.value = x; -	u.bits &= (uint64_t)-1 / 2; -	return u.value; +	union {double f; uint64_t i;} u = {x}; +	u.i &= -1ULL/2; +	return u.f;  } diff --git a/src/math/fabsf.c b/src/math/fabsf.c index 516f1104..4efc8d68 100644 --- a/src/math/fabsf.c +++ b/src/math/fabsf.c @@ -1,10 +1,9 @@ -#include "libm.h" +#include <math.h> +#include <stdint.h>  float fabsf(float x)  { -	union fshape u; - -	u.value = x; -	u.bits &= (uint32_t)-1 / 2; -	return u.value; +	union {float f; uint32_t i;} u = {x}; +	u.i &= 0x7fffffff; +	return u.f;  } diff --git a/src/math/nextafter.c b/src/math/nextafter.c index 9ee82518..ab5795a4 100644 --- a/src/math/nextafter.c +++ b/src/math/nextafter.c @@ -1,35 +1,31 @@  #include "libm.h" -#define SIGN ((uint64_t)1<<63) -  double nextafter(double x, double y)  { -	union dshape ux, uy; +	union {double f; uint64_t i;} ux={x}, uy={y};  	uint64_t ax, ay;  	int e;  	if (isnan(x) || isnan(y))  		return x + y; -	ux.value = x; -	uy.value = y; -	if (ux.bits == uy.bits) +	if (ux.i == uy.i)  		return y; -	ax = ux.bits & ~SIGN; -	ay = uy.bits & ~SIGN; +	ax = ux.i & -1ULL/2; +	ay = uy.i & -1ULL/2;  	if (ax == 0) {  		if (ay == 0)  			return y; -		ux.bits = (uy.bits & SIGN) | 1; -	} else if (ax > ay || ((ux.bits ^ uy.bits) & SIGN)) -		ux.bits--; +		ux.i = (uy.i & 1ULL<<63) | 1; +	} else if (ax > ay || ((ux.i ^ uy.i) & 1ULL<<63)) +		ux.i--;  	else -		ux.bits++; -	e = ux.bits >> 52 & 0x7ff; -	/* raise overflow if ux.value is infinite and x is finite */ +		ux.i++; +	e = ux.i >> 52 & 0x7ff; +	/* raise overflow if ux.f is infinite and x is finite */  	if (e == 0x7ff)  		FORCE_EVAL(x+x); -	/* raise underflow if ux.value is subnormal or zero */ +	/* raise underflow if ux.f is subnormal or zero */  	if (e == 0) -		FORCE_EVAL(x*x + ux.value*ux.value); -	return ux.value; +		FORCE_EVAL(x*x + ux.f*ux.f); +	return ux.f;  } diff --git a/src/math/nextafterf.c b/src/math/nextafterf.c index 22b61dce..75a09f7d 100644 --- a/src/math/nextafterf.c +++ b/src/math/nextafterf.c @@ -1,34 +1,30 @@  #include "libm.h" -#define SIGN 0x80000000 -  float nextafterf(float x, float y)  { -	union fshape ux, uy; +	union {float f; uint32_t i;} ux={x}, uy={y};  	uint32_t ax, ay, e;  	if (isnan(x) || isnan(y))  		return x + y; -	ux.value = x; -	uy.value = y; -	if (ux.bits == uy.bits) +	if (ux.i == uy.i)  		return y; -	ax = ux.bits & ~SIGN; -	ay = uy.bits & ~SIGN; +	ax = ux.i & 0x7fffffff; +	ay = uy.i & 0x7fffffff;  	if (ax == 0) {  		if (ay == 0)  			return y; -		ux.bits = (uy.bits & SIGN) | 1; -	} else if (ax > ay || ((ux.bits ^ uy.bits) & SIGN)) -		ux.bits--; +		ux.i = (uy.i & 0x80000000) | 1; +	} else if (ax > ay || ((ux.i ^ uy.i) & 0x80000000)) +		ux.i--;  	else -		ux.bits++; -	e = ux.bits & 0x7f800000; -	/* raise overflow if ux.value is infinite and x is finite */ +		ux.i++; +	e = ux.i & 0x7f800000; +	/* raise overflow if ux.f is infinite and x is finite */  	if (e == 0x7f800000)  		FORCE_EVAL(x+x); -	/* raise underflow if ux.value is subnormal or zero */ +	/* raise underflow if ux.f is subnormal or zero */  	if (e == 0) -		FORCE_EVAL(x*x + ux.value*ux.value); -	return ux.value; +		FORCE_EVAL(x*x + ux.f*ux.f); +	return ux.f;  } diff --git a/src/math/nexttoward.c b/src/math/nexttoward.c index 6f32eca4..827ee5c3 100644 --- a/src/math/nexttoward.c +++ b/src/math/nexttoward.c @@ -6,40 +6,37 @@ double nexttoward(double x, long double y)  	return nextafter(x, y);  }  #else -#define SIGN ((uint64_t)1<<63) -  double nexttoward(double x, long double y)  { -	union dshape ux; +	union {double f; uint64_t i;} ux = {x};  	int e;  	if (isnan(x) || isnan(y))  		return x + y;  	if (x == y)  		return y; -	ux.value = x;  	if (x == 0) { -		ux.bits = 1; +		ux.i = 1;  		if (signbit(y)) -			ux.bits |= SIGN; +			ux.i |= 1ULL<<63;  	} else if (x < y) {  		if (signbit(x)) -			ux.bits--; +			ux.i--;  		else -			ux.bits++; +			ux.i++;  	} else {  		if (signbit(x)) -			ux.bits++; +			ux.i++;  		else -			ux.bits--; +			ux.i--;  	} -	e = ux.bits>>52 & 0x7ff; -	/* raise overflow if ux.value is infinite and x is finite */ +	e = ux.i>>52 & 0x7ff; +	/* raise overflow if ux.f is infinite and x is finite */  	if (e == 0x7ff)  		FORCE_EVAL(x+x); -	/* raise underflow if ux.value is subnormal or zero */ +	/* raise underflow if ux.f is subnormal or zero */  	if (e == 0) -		FORCE_EVAL(x*x + ux.value*ux.value); -	return ux.value; +		FORCE_EVAL(x*x + ux.f*ux.f); +	return ux.f;  }  #endif diff --git a/src/math/nexttowardf.c b/src/math/nexttowardf.c index 9a693b1a..bbf172f9 100644 --- a/src/math/nexttowardf.c +++ b/src/math/nexttowardf.c @@ -2,35 +2,34 @@  float nexttowardf(float x, long double y)  { -	union fshape ux; +	union {float f; uint32_t i;} ux = {x};  	uint32_t e;  	if (isnan(x) || isnan(y))  		return x + y;  	if (x == y)  		return y; -	ux.value = x;  	if (x == 0) { -		ux.bits = 1; +		ux.i = 1;  		if (signbit(y)) -			ux.bits |= 0x80000000; +			ux.i |= 0x80000000;  	} else if (x < y) {  		if (signbit(x)) -			ux.bits--; +			ux.i--;  		else -			ux.bits++; +			ux.i++;  	} else {  		if (signbit(x)) -			ux.bits++; +			ux.i++;  		else -			ux.bits--; +			ux.i--;  	} -	e = ux.bits & 0x7f800000; -	/* raise overflow if ux.value is infinite and x is finite */ +	e = ux.i & 0x7f800000; +	/* raise overflow if ux.f is infinite and x is finite */  	if (e == 0x7f800000)  		FORCE_EVAL(x+x); -	/* raise underflow if ux.value is subnormal or zero */ +	/* raise underflow if ux.f is subnormal or zero */  	if (e == 0) -		FORCE_EVAL(x*x + ux.value*ux.value); -	return ux.value; +		FORCE_EVAL(x*x + ux.f*ux.f); +	return ux.f;  } | 
