summaryrefslogtreecommitdiff
path: root/Build/source/libs/mpfr/mpfr-src/src/subnormal.c
diff options
context:
space:
mode:
Diffstat (limited to 'Build/source/libs/mpfr/mpfr-src/src/subnormal.c')
-rw-r--r--Build/source/libs/mpfr/mpfr-src/src/subnormal.c31
1 files changed, 24 insertions, 7 deletions
diff --git a/Build/source/libs/mpfr/mpfr-src/src/subnormal.c b/Build/source/libs/mpfr/mpfr-src/src/subnormal.c
index c9762505775..76fccf0576f 100644
--- a/Build/source/libs/mpfr/mpfr-src/src/subnormal.c
+++ b/Build/source/libs/mpfr/mpfr-src/src/subnormal.c
@@ -1,7 +1,7 @@
/* mpfr_subnormalize -- Subnormalize a floating point number
emulating sub-normal numbers.
-Copyright 2005-2019 Free Software Foundation, Inc.
+Copyright 2005-2020 Free Software Foundation, Inc.
Contributed by the AriC and Caramba projects, INRIA.
This file is part of the GNU MPFR Library.
@@ -70,7 +70,7 @@ mpfr_subnormalize (mpfr_ptr y, int old_inexact, mpfr_rnd_t rnd)
and since y is not a power of 2: 0.5*2^emin < Y < 1*2^emin
We also know the direction of the error thanks to ternary value. */
- if (rnd == MPFR_RNDN)
+ if (rnd == MPFR_RNDN || rnd == MPFR_RNDNA)
{
mp_limb_t *mant, rb ,sb;
mp_size_t s;
@@ -95,7 +95,7 @@ mpfr_subnormalize (mpfr_ptr y, int old_inexact, mpfr_rnd_t rnd)
Otherwise, rounding bit = 1, sticky bit = 0 and inexact = 0
So we have 0.1100000000000000000000000*2^emin exactly.
We return 0.1*2^(emin+1) according to the even-rounding
- rule on subnormals. */
+ rule on subnormals. Note the same holds for RNDNA. */
goto set_min_p1;
}
else if (MPFR_IS_LIKE_RNDZ (rnd, MPFR_IS_NEG (y)))
@@ -116,6 +116,7 @@ mpfr_subnormalize (mpfr_ptr y, int old_inexact, mpfr_rnd_t rnd)
{
mpfr_t dest;
mpfr_prec_t q;
+ mpfr_rnd_t rnd2;
int inexact, inex2;
MPFR_ASSERTD (MPFR_GET_EXP (y) > __gmpfr_emin);
@@ -129,17 +130,20 @@ mpfr_subnormalize (mpfr_ptr y, int old_inexact, mpfr_rnd_t rnd)
/* Round y in dest */
MPFR_SET_EXP (dest, MPFR_GET_EXP (y));
MPFR_SET_SIGN (dest, sign);
+ rnd2 = rnd == MPFR_RNDNA ? MPFR_RNDN : rnd;
MPFR_RNDRAW_EVEN (inexact, dest,
- MPFR_MANT (y), MPFR_PREC (y), rnd, sign,
+ MPFR_MANT (y), MPFR_PREC (y), rnd2, sign,
MPFR_SET_EXP (dest, MPFR_GET_EXP (dest) + 1));
if (MPFR_LIKELY (old_inexact != 0))
{
- if (MPFR_UNLIKELY (rnd == MPFR_RNDN &&
+ if (MPFR_UNLIKELY (rnd2 == MPFR_RNDN &&
(inexact == MPFR_EVEN_INEX ||
inexact == -MPFR_EVEN_INEX)))
{
- /* if both roundings are in the same direction, we have to go
- back in the other direction */
+ /* If both roundings are in the same direction,
+ we have to go back in the other direction.
+ For MPFR_RNDNA it is the same, since we are not
+ exactly in the middle case (old_inexact != 0). */
if (SAME_SIGN (inexact, old_inexact))
{
if (SAME_SIGN (inexact, MPFR_INT_SIGN (y)))
@@ -155,6 +159,19 @@ mpfr_subnormalize (mpfr_ptr y, int old_inexact, mpfr_rnd_t rnd)
else if (MPFR_UNLIKELY (inexact == 0))
inexact = old_inexact;
}
+ else if (rnd == MPFR_RNDNA &&
+ (inexact == MPFR_EVEN_INEX || inexact == -MPFR_EVEN_INEX))
+ {
+ /* We are in the middle case: since we used RNDN to round, we should
+ round in the opposite direction when inexact has the opposite
+ sign of y. */
+ if (!SAME_SIGN (inexact, MPFR_INT_SIGN (y)))
+ {
+ mpfr_nexttoinf (dest);
+ MPFR_ASSERTD(!MPFR_IS_INF (dest));
+ inexact = -inexact;
+ }
+ }
inex2 = mpfr_set (y, dest, rnd);
MPFR_ASSERTN (inex2 == 0);