diff options
-rw-r--r-- | NEWS | 6 | ||||
-rw-r--r-- | TODO | 6 | ||||
-rw-r--r-- | doc/mpc.texi | 6 | ||||
-rw-r--r-- | src/div.c | 21 | ||||
-rw-r--r-- | tests/div.dat | 4 |
5 files changed, 18 insertions, 25 deletions
@@ -4,8 +4,10 @@ Changes in version 1.0: - Bug fixes: - mul and norm now return correct results even if there are over- or underflows during the computation - - asin, proj, sqr: wrong result when input variable has infinite part and - equals output variable + - div now returns a value indicating the effective rounding direction, + as the other functions. + - asin, proj, sqr: Wrong result when input variable has infinite part and + equals output variable is corrected. Changes in version 0.9: - New functions: @@ -1,9 +1,3 @@ -From Paul Zimmermann 1st September 2011: -return the exact ternary value in mpc_div: -- for directed rounding modes it is easy (assuming we are able to detect - exact cases) -- rounding to nearest to n bits is equivalent to directed rounding to n+1 bits - From Andreas Enge 31 August 2011: implement mul_karatsuba with three multiplications at precision around p, instead of two at precision 2*p and one at precision p diff --git a/doc/mpc.texi b/doc/mpc.texi index f8c6f42..8e400ef 100644 --- a/doc/mpc.texi +++ b/doc/mpc.texi @@ -881,11 +881,6 @@ Set @var{rop} to the square of @var{op} rounded according to @var{rnd}. @deftypefunx int mpc_ui_div (mpc_t @var{rop}, unsigned long int @var{op1}, mpc_t @var{op2}, mpc_rnd_t @var{rnd}) @deftypefunx int mpc_fr_div (mpc_t @var{rop}, mpfr_t @var{op1}, mpc_t @var{op2}, mpc_rnd_t @var{rnd}) Set @var{rop} to @var{op1}/@var{op2} rounded according to @var{rnd}. -For @code{mpc_div}, @code{mpc_ui_div} and @code{mpc_fr_div}, the usual -semantics of the return value is not yet implemented. We only guarantee -that a return value of 0 indicates an exact result. -But exact results may lead to a non-zero return value, and the sign of a -non-zero return value is not significant. @end deftypefun @deftypefun int mpc_neg (mpc_t @var{rop}, mpc_t @var{op}, mpc_rnd_t @var{rnd}) @@ -902,7 +897,6 @@ if @var{rop} and @var{op} are the same variable. @deftypefun int mpc_abs (mpfr_t @var{rop}, mpc_t @var{op}, mpfr_rnd_t @var{rnd}) Set the floating-point number @var{rop} to the absolute value of @var{op}, rounded in the direction @var{rnd}. -The returned value is zero iff the result is exact. @end deftypefun @deftypefun int mpc_norm (mpfr_t @var{rop}, mpc_t @var{op}, mpfr_rnd_t @var{rnd}) @@ -1,6 +1,6 @@ /* mpc_div -- Divide two complex numbers. -Copyright (C) 2002, 2003, 2004, 2005, 2008, 2009, 2010 INRIA +Copyright (C) 2002, 2003, 2004, 2005, 2008, 2009, 2010, 2011 INRIA This file is part of GNU MPC. @@ -227,9 +227,10 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) mpc_t res, c_conj; mpfr_t q; mpfr_prec_t prec; - int inexact_prod, inexact_norm, inexact_re, inexact_im, loops = 0; + int inex, inexact_prod, inexact_norm, inexact_re, inexact_im, loops = 0; int underflow_norm, overflow_norm; int underflow_re, overflow_re, underflow_im = 0, overflow_im = 0; + mpfr_rnd_t rnd_re = MPC_RND_RE (rnd), rnd_im = MPC_RND_IM (rnd); /* According to the C standard G.3, there are three types of numbers: */ /* finite (both parts are usual real numbers; contains 0), infinite */ @@ -299,7 +300,7 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) overflow_re = mpfr_overflow_p (); ok_re = !inexact_re || underflow_re || overflow_re || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDZ, - MPC_RND_RE(rnd), MPC_PREC_RE(a)); + GMP_RNDZ, MPC_PREC_RE(a) + (rnd_re == GMP_RNDN)); if (ok_re) /* compute imaginary part */ { mpfr_clear_underflow (); @@ -308,8 +309,8 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) 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, - MPC_RND_IM(rnd), MPC_PREC_IM(a)); + || mpfr_can_round (MPC_IM (res), prec - 4, GMP_RNDN, + GMP_RNDZ, MPC_PREC_IM(a) + (rnd_im == GMP_RNDN)); } } else { @@ -331,7 +332,7 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) overflow_re = mpfr_overflow_p (); ok_re = !inexact_re || underflow_re || overflow_re || mpfr_can_round (MPC_RE (res), prec - 4, GMP_RNDZ, - MPC_RND_RE(rnd), MPC_PREC_RE(a)); + GMP_RNDZ, MPC_PREC_RE(a) + (rnd_re == GMP_RNDN)); if (ok_re) /* compute imaginary part */ { mpfr_clear_underflow (); @@ -342,14 +343,16 @@ mpc_div (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd) overflow_im = mpfr_overflow_p (); ok_im = !inexact_im || underflow_im || overflow_im || mpfr_can_round (MPC_IM (res), prec - 4, GMP_RNDZ, - MPC_RND_IM(rnd), MPC_PREC_IM(a)); + GMP_RNDZ, MPC_PREC_IM(a) + (rnd_im == GMP_RNDN)); } } } while ((!ok_re || !ok_im) && !underflow_norm && !overflow_norm); - mpc_set (a, res, rnd); - + 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 */ if (overflow_re || underflow_norm) { mpfr_set_inf (MPC_RE (a), mpfr_sgn (MPC_RE (res))); diff --git a/tests/div.dat b/tests/div.dat index 4d7c1b5..ff2e1f3 100644 --- a/tests/div.dat +++ b/tests/div.dat @@ -2428,13 +2428,13 @@ - + 250 -inf 250 +inf 250 1 250 0 250 -1e-164895850 250 -1e-164895850 N N # small exact/inexact examples -#+ + 10 0b1.01010001 10 -0b1.1000101@-6 10 973 10 964 10 725 10 745 N N +- - 10 0b1.01010001 10 -0b1.1000101@-6 10 973 10 964 10 725 10 745 N N 0 0 10 -14 10 9 10 -837 10 637 10 63 10 -5 N N 0 0 2 2 2 -1 2 4 2 3 2 1 2 2 N N + - 4 1.375 4 1.25 4 15 4 14 4 11 4 0 N N # Bug 20080923 -#+ - 4 0b11@527 4 -0b111@-489 4 -0b11@-206 4 0 4 -0b1@-733 4 -0b101@-1750 N Z ++ + 4 0b11@527 4 -0b111@-489 4 -0b11@-206 4 0 4 -0b1@-733 4 -0b101@-1750 N Z # potential intermediate over- or underflow 0 0 10 1 10 0 10 0 10 0b1@536870912 10 0 10 0b1@536870912 N N |