From adca15f4321606f21c681579829f6929aa37561a Mon Sep 17 00:00:00 2001 From: vlefevre Date: Tue, 11 Feb 2020 12:22:29 +0000 Subject: [tests/tfmma.c] Merged the latest tests from the trunk (r12752,12759,13688,13696). [src/sub1.c] Bug fix: the underflow case (possible with UBF, e.g. via mpfr_fmma or mpfr_fmms) was not tested in the case c small. (merged changeset r13694 from the trunk) git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/branches/4.0@13697 280ebfd0-de03-0410-8827-d642c229c3f4 --- src/sub1.c | 24 ++++++++------ tests/tfmma.c | 100 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 114 insertions(+), 10 deletions(-) diff --git a/src/sub1.c b/src/sub1.c index 9337cc235..5b475d9d3 100644 --- a/src/sub1.c +++ b/src/sub1.c @@ -123,11 +123,11 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (rnd_mode == MPFR_RNDF) return mpfr_set4 (a, b, MPFR_RNDZ, MPFR_SIGN (a)); - MPFR_EXP (a) = exp_b; /* may be up to MPFR_EXP_MAX */ + exp_a = exp_b; /* may be any out-of-range value due to UBF */ MPFR_RNDRAW_EVEN (inexact, a, MPFR_MANT (b), bq, rnd_mode, MPFR_SIGN (a), - if (MPFR_EXP (a) != MPFR_EXP_MAX) - ++ MPFR_EXP (a)); + if (exp_a != MPFR_EXP_MAX) + exp_a ++); MPFR_LOG_MSG (("inexact=%d\n", inexact)); if (inexact == 0) { @@ -139,7 +139,7 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) if (! MPFR_IS_LIKE_RNDZ (rnd_mode, MPFR_IS_NEG (a))) { inexact = MPFR_INT_SIGN (a); - goto check_overflow; + goto end_of_c_small; } } else /* inexact != 0 */ @@ -164,7 +164,7 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) which means we get a wrong rounded result if x == 1, i.e. inexact == MPFR_EVEN_INEX (for positive numbers). */ if (MPFR_LIKELY (inexact != MPFR_EVEN_INEX * MPFR_INT_SIGN (a))) - goto check_overflow; + goto end_of_c_small; } /* We need to take the value preceding |a|. We can't use mpfr_nexttozero due to a possible out-of-range exponent. @@ -174,16 +174,20 @@ mpfr_sub1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) mpn_sub_1 (ap, ap, an, MPFR_LIMB_ONE << sh); if (MPFR_UNLIKELY (MPFR_LIMB_MSB (ap[an-1]) == 0)) { - MPFR_EXP (a) --; + exp_a --; /* The following is valid whether an = 1 or an > 1. */ ap[an-1] |= MPFR_LIMB_HIGHBIT; } inexact = - MPFR_INT_SIGN (a); - check_overflow: - if (MPFR_UNLIKELY (MPFR_EXP (a) > __gmpfr_emax)) + end_of_c_small: + /* The underflow case is possibly only with UBF. The overflow case + is also possible with normal FP due to rounding. */ + if (MPFR_UNLIKELY (exp_a > __gmpfr_emax)) return mpfr_overflow (a, rnd_mode, MPFR_SIGN (a)); - else - MPFR_RET (inexact); + if (MPFR_UNLIKELY (exp_a < __gmpfr_emin)) + goto underflow; + MPFR_SET_EXP (a, exp_a); + MPFR_RET (inexact); } /* reserve a space to store b aligned with the result, i.e. shifted by diff --git a/tests/tfmma.c b/tests/tfmma.c index 602603377..a4243cca2 100644 --- a/tests/tfmma.c +++ b/tests/tfmma.c @@ -480,6 +480,103 @@ bug20170405 (void) mpfr_clears (x, y, z, (mpfr_ptr) 0); } +static void +underflow_tests (void) +{ + mpfr_t x, y, z; + mpfr_prec_t p; + mpfr_exp_t emin; + + emin = mpfr_get_emin (); + mpfr_set_emin (-17); + for (p = MPFR_PREC_MIN; p <= 1024; p++) + { + mpfr_inits2 (p, x, y, (mpfr_ptr) 0); + mpfr_init2 (z, p + 1); + mpfr_set_str_binary (x, "1e-10"); + mpfr_set_str_binary (y, "1e-11"); + mpfr_clear_underflow (); + mpfr_fmms (z, x, x, y, y, MPFR_RNDN); + /* the exact result is 2^-20-2^-22, and should underflow */ + MPFR_ASSERTN(mpfr_underflow_p ()); + MPFR_ASSERTN(mpfr_zero_p (z)); + MPFR_ASSERTN(mpfr_signbit (z) == 0); + mpfr_clears (x, y, z, (mpfr_ptr) 0); + } + mpfr_set_emin (emin); +} + +static void +bug20180604 (void) +{ + mpfr_t x, y, yneg, z; + mpfr_exp_t emin; + int ret; + + emin = mpfr_get_emin (); + mpfr_set_emin (-1073741821); + mpfr_inits2 (564, x, y, yneg, (mpfr_ptr) 0); + mpfr_init2 (z, 2256); + mpfr_set_str_binary (x, "1.10010000111100110011001101111111101000111001011000110100110010000101000100010001000000111100010000101001011011111001111000110101111100101111001100001100011101100100011110000000011000010110111100111000100101010001011111010111011001110010001011101111001011001110110000010000011100010001010001011100100110111110101001001111001011101111110011101110101010110100010010111011111100010101111100011110111001011111101110101101101110100101111010000101011110100000000110111101000001100001000100010110100111010011011010110011100111010000101110010101111001011100110101100001100E-737194993"); + mpfr_set_str_binary (y, "-1.00101000100001001101011110100010110011101010011011010111100110101011111100000100101000111010111101100100110010001110011011100100110110000001011001000111101111101111110101100110111000000011000001101001010100010010001110001000011010000100111001001100101111111100010101110101001101101101111010100011011110001000010000010100011000011000010110101100000111111110111001100100100001101011111011100101110111000100101010110100010011101010110010100110100111000000100111101101101000000011110000100110100100011000010011110010001010000110100011111101101101110001110001101101010E-737194903"); + + mpfr_clear_underflow (); + ret = mpfr_fmms (z, x, x, y, y, MPFR_RNDN); + MPFR_ASSERTN(mpfr_underflow_p ()); + MPFR_ASSERTN(mpfr_zero_p (z)); + MPFR_ASSERTN(mpfr_signbit (z) == 1); + MPFR_ASSERTN(ret > 0); + + mpfr_clear_underflow (); + ret = mpfr_fmms (z, y, y, x, x, MPFR_RNDN); + MPFR_ASSERTN(mpfr_underflow_p ()); + MPFR_ASSERTN(mpfr_zero_p (z)); + MPFR_ASSERTN(mpfr_signbit (z) == 0); + MPFR_ASSERTN(ret < 0); + + mpfr_neg (yneg, y, MPFR_RNDN); + mpfr_clear_underflow (); + ret = mpfr_fmms (z, x, x, y, yneg, MPFR_RNDN); + MPFR_ASSERTN(mpfr_underflow_p ()); + MPFR_ASSERTN(mpfr_zero_p (z)); + MPFR_ASSERTN(mpfr_signbit (z) == 0); + MPFR_ASSERTN(ret < 0); + + mpfr_clear_underflow (); + ret = mpfr_fmms (z, y, yneg, x, x, MPFR_RNDN); + MPFR_ASSERTN(mpfr_underflow_p ()); + MPFR_ASSERTN(mpfr_zero_p (z)); + MPFR_ASSERTN(mpfr_signbit (z) == 1); + MPFR_ASSERTN(ret > 0); + + mpfr_clears (x, y, yneg, z, (mpfr_ptr) 0); + mpfr_set_emin (emin); +} + +/* this bug was discovered from mpc_mul */ +static void +bug20200206 (void) +{ + mpfr_exp_t emin = mpfr_get_emin (); + mpfr_t xre, xim, yre, yim, zre; + + mpfr_set_emin (-1073); + mpfr_inits2 (53, xre, xim, yre, yim, zre, (mpfr_ptr) 0); + mpfr_set_str (xre, "-0x8.294611b331c8p-904", 16, MPFR_RNDN); + mpfr_set_str (xim, "-0x1.859278c2992fap-676", 16, MPFR_RNDN); + mpfr_set_str (yre, "0x9.ac54802a95f8p-820", 16, MPFR_RNDN); + mpfr_set_str (yim, "0x3.17e59e7612aap-412", 16, MPFR_RNDN); + mpfr_fmms (zre, xre, yre, xim, yim, MPFR_RNDN); + if (mpfr_regular_p (zre) && mpfr_get_exp (zre) < -1073) + { + printf ("Error, mpfr_fmms returns an out-of-range exponent:\n"); + mpfr_dump (zre); + exit (1); + } + mpfr_clears (xre, xim, yre, yim, zre, (mpfr_ptr) 0); + mpfr_set_emin (emin); +} + /* Test double-rounding cases in mpfr_set_1_2, which is called when all the precisions are the same. With set.c before r13347, this triggers errors for neg=0. */ @@ -548,6 +645,9 @@ main (int argc, char *argv[]) { tests_start_mpfr (); + bug20200206 (); + bug20180604 (); + underflow_tests (); random_tests (); zero_tests (); max_tests (); -- cgit v1.2.1