diff options
author | enge <enge@211d60ee-9f03-0410-a15a-8952a2c7a4e4> | 2011-09-02 18:48:08 +0000 |
---|---|---|
committer | enge <enge@211d60ee-9f03-0410-a15a-8952a2c7a4e4> | 2011-09-02 18:48:08 +0000 |
commit | 62e34d750708e9121ac790841cbd715630fb111f (patch) | |
tree | 81d39cb1e719d853d981bcf0dbe0c0a5ac932d52 | |
parent | 80d7ee0c37146ab7acfb22766d3d4770a76ced8c (diff) | |
download | mpc-62e34d750708e9121ac790841cbd715630fb111f.tar.gz |
div.c: fixed bad application of mpfr_can_round trick - do not put in the same
directed rounding mode twice
div.dat: added example exposing the bug
git-svn-id: svn://scm.gforge.inria.fr/svn/mpc/trunk@1078 211d60ee-9f03-0410-a15a-8952a2c7a4e4
-rw-r--r-- | src/div.c | 25 | ||||
-rw-r--r-- | tests/div.dat | 3 |
2 files changed, 13 insertions, 15 deletions
@@ -18,7 +18,6 @@ You should have received a copy of the GNU Lesser General Public License along with this program. If not, see http://www.gnu.org/licenses/ . */ -#include <stdio.h> /* for MPC_ASSERT */ #include "mpc-impl.h" static int @@ -304,7 +303,6 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) if (!mpfr_zero_p (MPC_IM (res))) mpfr_set_inf (MPC_IM (res), mpfr_sgn (MPC_IM (res))); } - /* divide the product by the norm */ if (inexact_norm == 0 && (inexact_re == 0 || inexact_im == 0)) { @@ -317,7 +315,7 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) underflow_re = mpfr_underflow_p (); overflow_re = mpfr_overflow_p (); ok_re = !inexact_re || underflow_re || overflow_re - || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDZ, + || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDN, GMP_RNDZ, MPC_PREC_RE(a) + (rnd_re == GMP_RNDN)); if (ok_re) /* compute imaginary part */ { @@ -338,39 +336,36 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) /* if 1/q is inexact, the approximations of the real and imaginary part below will be inexact, unless RE(res) or IM(res) is zero */ - inexact_re = inexact_re || !mpfr_zero_p (MPC_RE (res)); - inexact_im = inexact_im || !mpfr_zero_p (MPC_IM (res)); + inexact_re |= ~mpfr_zero_p (MPC_RE (res)); + inexact_im |= ~mpfr_zero_p (MPC_IM (res)); } mpfr_clear_underflow (); mpfr_clear_overflow (); - inexact_re = mpfr_mul (MPC_RE (res), MPC_RE (res), q, GMP_RNDZ) - || inexact_re; + inexact_re |= mpfr_mul (MPC_RE (res), MPC_RE (res), q, GMP_RNDZ); underflow_re = mpfr_underflow_p (); overflow_re = mpfr_overflow_p (); ok_re = !inexact_re || underflow_re || overflow_re - || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDZ, + || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDN, GMP_RNDZ, MPC_PREC_RE(a) + (rnd_re == GMP_RNDN)); if (ok_re) /* compute imaginary part */ { mpfr_clear_underflow (); mpfr_clear_overflow (); - inexact_im = mpfr_mul (MPC_IM (res), MPC_IM (res), q, GMP_RNDZ) - || inexact_im; + inexact_im != mpfr_mul (MPC_IM (res), MPC_IM (res), q, GMP_RNDZ); underflow_im = mpfr_underflow_p (); overflow_im = mpfr_overflow_p (); ok_im = !inexact_im || underflow_im || overflow_im - || mpfr_can_round (MPC_IM (res), prec - 4, GMP_RNDZ, + || mpfr_can_round (MPC_IM (res), prec - 4, GMP_RNDN, GMP_RNDZ, MPC_PREC_IM(a) + (rnd_im == GMP_RNDN)); } } - } - while ((!ok_re || !ok_im) && !underflow_norm && !overflow_norm - && !underflow_prod && !overflow_prod); + } while ((!ok_re || !ok_im) && !underflow_norm && !overflow_norm + && !underflow_prod && !overflow_prod); inex = mpc_set (a, res, rnd); inexact_re = MPC_INEX_RE (inex); inexact_im = MPC_INEX_IM (inex); - + /* fix values and inexact flags in case of overflow/underflow */ /* FIXME: heuristic, certainly does not cover all cases */ if (overflow_re || (underflow_norm && !underflow_prod)) { diff --git a/tests/div.dat b/tests/div.dat index 2459b59..8bbebf7 100644 --- a/tests/div.dat +++ b/tests/div.dat @@ -2440,6 +2440,9 @@ # overflow (reported by Emmanuel Thome) - + 250 -inf 250 +inf 250 1 250 0 250 -1e-164895850 250 -1e-164895850 N N +# bug found by tgeneric of ui_div ++ + 2 0b1.1@256 2 0b1.1@-2758 34 52349199244 2 0 2 0b1.1@-221 2 -0b1@-3234 U N + # cases that should yield 1, but cannot be handled due to intermediate # over- or underflows # current result: (@NaN@ 0) |