diff options
-rw-r--r-- | math/exp.c | 7 | ||||
-rw-r--r-- | math/exp2.c | 5 | ||||
-rw-r--r-- | math/log.c | 1 | ||||
-rw-r--r-- | math/log2.c | 1 | ||||
-rw-r--r-- | math/math_config.h | 26 | ||||
-rw-r--r-- | math/pow.c | 11 |
6 files changed, 50 insertions, 1 deletions
@@ -33,6 +33,10 @@ #define C5 __exp_data.poly[8 - EXP_POLY_ORDER] #define C6 __exp_data.poly[9 - EXP_POLY_ORDER] +/* Handle inputs that may overflow or underflow when computing the result + that is scale*(1+tmp), the exponent bits of scale might have overflown + into the sign bit so that needs correction before sbits is used as a + double, ki is only used to determine the sign of the exponent. */ static inline double specialcase (double_t tmp, uint64_t sbits, uint64_t ki) { @@ -71,12 +75,15 @@ specialcase (double_t tmp, uint64_t sbits, uint64_t ki) return check_uflow (y); } +/* Top 12 bits of a double (sign and exponent bits). */ static inline uint32_t top12 (double x) { return asuint64 (x) >> 52; } +/* Computes exp(x+xtail) where |xtail| < 2^-8/N and |xtail| <= |x|. + If hastail is 0 then xtail is assumed to be 0 too. */ static inline double exp_inline (double x, double xtail, int hastail) { diff --git a/math/exp2.c b/math/exp2.c index 3d8aadf..a64a26c 100644 --- a/math/exp2.c +++ b/math/exp2.c @@ -31,6 +31,10 @@ #define C5 __exp_data.exp2_poly[4] #define C6 __exp_data.exp2_poly[5] +/* Handle inputs that may overflow or underflow when computing the result + that is scale*(1+tmp), the exponent bits of scale might have overflown + into the sign bit so that needs correction before sbits is used as a + double, ki is only used to determine the sign of the exponent. */ static inline double specialcase (double_t tmp, uint64_t sbits, uint64_t ki) { @@ -69,6 +73,7 @@ specialcase (double_t tmp, uint64_t sbits, uint64_t ki) return check_uflow (y); } +/* Top 12 bits of a double (sign and exponent bits). */ static inline uint32_t top12 (double x) { @@ -30,6 +30,7 @@ #define N (1 << LOG_TABLE_BITS) #define OFF 0x3fe6000000000000 +/* Top 16 bits of a double. */ static inline uint32_t top16 (double x) { diff --git a/math/log2.c b/math/log2.c index 2b3960d..60c2e44 100644 --- a/math/log2.c +++ b/math/log2.c @@ -30,6 +30,7 @@ #define N (1 << LOG2_TABLE_BITS) #define OFF 0x3fe6000000000000 +/* Top 16 bits of a double. */ static inline uint32_t top16 (double x) { diff --git a/math/math_config.h b/math/math_config.h index 3986d48..ee44468 100644 --- a/math/math_config.h +++ b/math/math_config.h @@ -66,12 +66,16 @@ #if HAVE_FAST_ROUND # define TOINT_INTRINSICS 1 +/* Round x to nearest int, ties have to be rounded consistently with + converttoint. The result is in [-2^31+1, 2^31-1]. */ static inline double_t roundtoint (double_t x) { return round (x); } +/* Convert x to nearest int, ties have to be rounded consistently with + roundtoint. The result is in [-2^31+1, 2^31-1]. */ static inline uint64_t converttoint (double_t x) { @@ -240,28 +244,48 @@ eval_as_double (double x) # define unlikely(x) (x) #endif -/* Error handling tail calls for special cases, with sign argument. */ +/* Error handling tail calls for special cases, with a sign argument. + The sign of the return value is set if the argument is non-zero. */ + +/* The result overflows. */ HIDDEN float __math_oflowf (uint32_t); +/* The result underflows to 0 in nearest rounding mode. */ HIDDEN float __math_uflowf (uint32_t); +/* The result underflows to 0 in some directed rounding mode only. */ HIDDEN float __math_may_uflowf (uint32_t); +/* Division by zero. */ HIDDEN float __math_divzerof (uint32_t); +/* The result overflows. */ HIDDEN double __math_oflow (uint32_t); +/* The result underflows to 0 in nearest rounding mode. */ HIDDEN double __math_uflow (uint32_t); +/* The result underflows to 0 in some directed rounding mode only. */ HIDDEN double __math_may_uflow (uint32_t); +/* Division by zero. */ HIDDEN double __math_divzero (uint32_t); + /* Error handling using input checking. */ + +/* Invalid input unless it is a quiet NaN. */ HIDDEN float __math_invalidf (float); +/* Invalid input unless it is a quiet NaN. */ HIDDEN double __math_invalid (double); + /* Error handling using output checking, only for errno setting. */ + +/* Check if the result overflowed to infinity. */ HIDDEN double __math_check_oflow (double); +/* Check if the result underflowed to 0. */ HIDDEN double __math_check_uflow (double); +/* Check if the result overflowed to infinity. */ static inline double check_oflow (double x) { return WANT_ERRNO ? __math_check_oflow (x) : x; } +/* Check if the result underflowed to 0. */ static inline double check_uflow (double x) { @@ -34,12 +34,16 @@ ulperr_exp: 0.509 ULP (ULP error of exp, 0.511 ULP without fma) #define N (1 << POW_LOG_TABLE_BITS) #define OFF 0x3fe6955500000000 +/* Top 12 bits of a double (sign and exponent bits). */ static inline uint32_t top12 (double x) { return asuint64 (x) >> 52; } +/* Compute y+tail = log(x) where the rounded result is y and tail has about + additional 15 bits precision. The bit representation of x if in ix, but + normalized in the subnormal range using sign bit too for the exponent. */ static inline double_t log_inline (uint64_t ix, double_t *tail) { @@ -122,6 +126,10 @@ log_inline (uint64_t ix, double_t *tail) #define C5 __exp_data.poly[8 - EXP_POLY_ORDER] #define C6 __exp_data.poly[9 - EXP_POLY_ORDER] +/* Handle inputs that may overflow or underflow when computing the result + that is scale*(1+tmp), the exponent bits of scale might have overflown + into the sign bit so that needs correction before sbits is used as a + double, ki is only used to determine the sign of the exponent. */ static inline double specialcase (double_t tmp, uint64_t sbits, uint64_t ki) { @@ -165,6 +173,8 @@ specialcase (double_t tmp, uint64_t sbits, uint64_t ki) #define SIGN_BIAS (0x800 << EXP_TABLE_BITS) +/* Computes sign*exp(x+xtail) where |xtail| < 2^-8/N and |xtail| <= |x|. + The sign_bias argument is SIGN_BIAS or 0 and sets the sign to -1 or 1. */ static inline double exp_inline (double x, double xtail, uint32_t sign_bias) { @@ -257,6 +267,7 @@ checkint (uint64_t iy) return 2; } +/* Returns 1 if input is the bit representation of 0, infinity or nan. */ static inline int zeroinfnan (uint64_t i) { |