diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/add1.c | 2 | ||||
| -rw-r--r-- | src/mpfr-impl.h | 21 | ||||
| -rw-r--r-- | src/sub1.c | 32 |
3 files changed, 49 insertions, 6 deletions
diff --git a/src/add1.c b/src/add1.c index 9f2859690..1e61eb905 100644 --- a/src/add1.c +++ b/src/add1.c @@ -41,7 +41,7 @@ mpfr_add1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (MPFR_UNLIKELY (MPFR_IS_UBF (b))) { - exp = mpfr_ubf_zexp2exp (MPFR_ZEXP (b)); + exp = MPFR_UBF_GET_EXP (b); if (exp > __gmpfr_emax) return mpfr_overflow (a, rnd_mode, MPFR_SIGN (b));; } diff --git a/src/mpfr-impl.h b/src/mpfr-impl.h index a083c7eab..c1d0b9306 100644 --- a/src/mpfr-impl.h +++ b/src/mpfr-impl.h @@ -2420,11 +2420,28 @@ __MPFR_DECLSPEC mpfr_exp_t mpfr_ubf_diff_exp (mpfr_srcptr, mpfr_srcptr); } #endif -#define MPFR_ZEXP(x) \ - ((void) (x)->_mpfr_exp /* to check that x has a correct type */, \ +/* Get the _mpfr_zexp field (pointer to a mpz_t) of a UBF object. + For practical reasons, the type of the argument x can be either + mpfr_ubf_ptr or mpfr_ptr, since the latter is used in functions + that accept both MPFR numbers and UBF's; this is checked by the + code "(x)->_mpfr_exp" (the "sizeof" prevents an access, which + could be invalid when MPFR_ZEXP(x) is used for an assignment, + and also avoids breaking the aliasing rules if they are dealt + with in the future). + This macro can be used when building a UBF. So we do not check + that the _mpfr_exp field has the value MPFR_EXP_UBF. */ +#define MPFR_ZEXP(x) \ + ((void) sizeof ((x)->_mpfr_exp), \ ((mpfr_ubf_ptr) (x))->_mpfr_zexp) +/* If x is a UBF, clear its mpz_t exponent. */ #define MPFR_UBF_CLEAR_EXP(x) \ ((void) (MPFR_IS_UBF (x) && (mpz_clear (MPFR_ZEXP (x)), 0))) +/* Like MPFR_GET_EXP, but accepts UBF (with exponent saturated to + the interval [MPFR_EXP_MIN,MPFR_EXP_MAX]). */ +#define MPFR_UBF_GET_EXP(x) \ + (MPFR_IS_UBF (x) ? mpfr_ubf_zexp2exp (MPFR_ZEXP (x)) : \ + MPFR_GET_EXP ((mpfr_ptr) (x))) + #endif /* __MPFR_IMPL_H__ */ diff --git a/src/sub1.c b/src/sub1.c index 5b475d9d3..6ee6ad544 100644 --- a/src/sub1.c +++ b/src/sub1.c @@ -91,9 +91,20 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (MPFR_UNLIKELY (MPFR_IS_UBF (b) || MPFR_IS_UBF (c))) { - exp_b = MPFR_IS_UBF (b) ? - mpfr_ubf_zexp2exp (MPFR_ZEXP (b)) : MPFR_GET_EXP (b); + exp_b = MPFR_UBF_GET_EXP (b); + /* Early underflow detection. Rare, but a test is needed anyway + since in the "MAX (aq, bq) + 2 <= diff_exp" branch, the exponent + may decrease and MPFR_EXP_MIN would yield an integer overflow. */ + if (MPFR_UNLIKELY (exp_b < __gmpfr_emin - 1)) + { + if (rnd_mode == MPFR_RNDN) + rnd_mode = MPFR_RNDZ; + return mpfr_underflow (a, rnd_mode, MPFR_SIGN(a)); + } diff_exp = mpfr_ubf_diff_exp (b, c); + /* mpfr_set4 below used with MPFR_RNDF does not support UBF. */ + if (rnd_mode == MPFR_RNDF) + rnd_mode = MPFR_RNDN; } else { @@ -185,7 +196,13 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (MPFR_UNLIKELY (exp_a > __gmpfr_emax)) return mpfr_overflow (a, rnd_mode, MPFR_SIGN (a)); if (MPFR_UNLIKELY (exp_a < __gmpfr_emin)) - goto underflow; + { + if (rnd_mode == MPFR_RNDN && + (exp_a < __gmpfr_emin - 1 || + (inexact * MPFR_INT_SIGN (a) >= 0 && mpfr_powerof2_raw (a)))) + rnd_mode = MPFR_RNDZ; + return mpfr_underflow (a, rnd_mode, MPFR_SIGN(a)); + } MPFR_SET_EXP (a, exp_a); MPFR_RET (inexact); } @@ -660,6 +677,15 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (MPFR_LIKELY(cancel)) { cancel -= add_exp; /* OK: add_exp is an int equal to 0 or 1 */ + MPFR_ASSERTD (cancel >= 0); + /* Detect an underflow case to avoid a possible integer overflow + with UBF in the computation of exp_a. */ + if (MPFR_UNLIKELY (exp_b < __gmpfr_emin - 1)) + { + if (rnd_mode == MPFR_RNDN) + rnd_mode = MPFR_RNDZ; + return mpfr_underflow (a, rnd_mode, MPFR_SIGN(a)); + } exp_a = exp_b - cancel; /* The following assertion corresponds to a limitation of the MPFR implementation. It may fail with a 32-bit ABI and huge precisions, |