diff options
| -rw-r--r-- | ChangeLog | 7 | ||||
| -rw-r--r-- | math/math_private.h | 19 | ||||
| -rw-r--r-- | sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c | 67 | ||||
| -rw-r--r-- | sysdeps/x86_64/fpu/math_private.h | 17 |
4 files changed, 110 insertions, 0 deletions
diff --git a/ChangeLog b/ChangeLog index c29660d652..b2d98e4ee7 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,5 +1,12 @@ 2009-08-24 Ulrich Drepper <drepper@redhat.com> + * math/math_private.h (ieee_double_shape_type): Add uint64_t word to + union. + (EXTRACT_WORDS64, INSERT_WORDS64): Define. + * sysdeps/x86_64/fpu/math_private.h (EXTRACT_WORDS64, INSERT_WORDS64): + Redefine. + * sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c: New file. + * sysdeps/x86_64/fpu/math_private.h: Add specialized GET_FLOAT_WORD and SET_FLOAT_WORD definitions. diff --git a/math/math_private.h b/math/math_private.h index 129646f8c5..fade7e1181 100644 --- a/math/math_private.h +++ b/math/math_private.h @@ -17,6 +17,7 @@ #define _MATH_PRIVATE_H_ #include <endian.h> +#include <stdint.h> #include <sys/types.h> /* The original fdlibm code used statements like: @@ -43,6 +44,7 @@ typedef union u_int32_t msw; u_int32_t lsw; } parts; + uint64_t word; } ieee_double_shape_type; #endif @@ -57,6 +59,7 @@ typedef union u_int32_t lsw; u_int32_t msw; } parts; + uint64_t word; } ieee_double_shape_type; #endif @@ -89,6 +92,14 @@ do { \ (i) = gl_u.parts.lsw; \ } while (0) +/* Get all in one, efficient on 64-bit machines. */ +#define EXTRACT_WORDS64(i,d) \ +do { \ + ieee_double_shape_type gh_u; \ + gh_u.value = (d); \ + (i) = gh_u.word; \ +} while (0) + /* Set a double from two 32 bit ints. */ #define INSERT_WORDS(d,ix0,ix1) \ @@ -99,6 +110,14 @@ do { \ (d) = iw_u.value; \ } while (0) +/* Get all in one, efficient on 64-bit machines. */ +#define INSERT_WORDS64(i,d) \ +do { \ + ieee_double_shape_type iw_u; \ + iw_u.word = (i); \ + (d) = iw_u.value; \ +} while (0) + /* Set the more significant 32 bits of a double from an int. */ #define SET_HIGH_WORD(d,v) \ diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c new file mode 100644 index 0000000000..9123fdc7bd --- /dev/null +++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c @@ -0,0 +1,67 @@ +/* @(#)s_ceil.c 5.1 93/09/24 */ +/* + * ==================================================== + * 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. + * ==================================================== + */ + +/* + * ceil(x) + * Return x rounded toward -inf to integral value + * Method: + * Bit twiddling. + * Exception: + * Inexact flag raised if x not equal to ceil(x). + */ + +#include "math.h" +#include "math_private.h" + +#ifdef __STDC__ +static const double huge = 1.0e300; +#else +static double huge = 1.0e300; +#endif + +#ifdef __STDC__ + double __ceil(double x) +#else + double __ceil(x) + double x; +#endif +{ + int64_t i0,i; + int32_t j0; + EXTRACT_WORDS64(i0,x); + j0 = ((i0>>52)&0x7ff)-0x3ff; + if(j0<=51) { + if(j0<0) { /* raise inexact if x != 0 */ + if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */ + if(i0<0) {i0=INT64_C(0x8000000000000000);} + else if(i0!=0) { i0=INT64_C(0x3ff0000000000000);} + } + } else { + i = INT64_C(0x000fffffffffffff)>>j0; + if((i0&i)==0) return x; /* x is integral */ + if(huge+x>0.0) { /* raise inexact flag */ + if(i0>0) i0 += UINT64_C(0x0010000000000000)>>j0; + i0 &= (~i); + } + } + } else { + if(j0==0x400) return x+x; /* inf or NaN */ + else return x; /* x is integral */ + } + INSERT_WORDS64(x,i0); + return x; +} +weak_alias (__ceil, ceil) +#ifdef NO_LONG_DOUBLE +strong_alias (__ceil, __ceill) +weak_alias (__ceil, ceill) +#endif diff --git a/sysdeps/x86_64/fpu/math_private.h b/sysdeps/x86_64/fpu/math_private.h index 8f4b792e77..4be753654a 100644 --- a/sysdeps/x86_64/fpu/math_private.h +++ b/sysdeps/x86_64/fpu/math_private.h @@ -22,6 +22,23 @@ while (0) /* We can do a few things better on x86-64. */ /* Direct movement of float into integer register. */ +#undef EXTRACT_WORDS64 +#define EXTRACT_WORDS64(i,d) \ +do { \ + long int i_; \ + asm ("movd %1, %0" : "=rm" (i_) : "x" (d)); \ + (i) = i_; \ +} while (0) + +/* And the reverse. */ +#undef INSERT_WORDS64 +#define INSERT_WORDS64(d,i) \ +do { \ + long int i_ = i; \ + asm ("movd %1, %0" : "=x" (d) : "rm" (i_)); \ +} while (0) + +/* Direct movement of float into integer register. */ #undef GET_FLOAT_WORD #define GET_FLOAT_WORD(i,d) \ do { \ |