fixed problems with mpfr-2.1.1 (can't use MPFR_EXP for x=0)

git-svn-id: svn://scm.gforge.inria.fr/svn/mpc/trunk@35 211d60ee-9f03-0410-a15a-8952a2c7a4e4

3 files changed, 64 insertions, 10 deletions

diff --git a/mul.c b/mul.c
index c756852..90b35a6 100644
--- a/mul.c
+++ b/mul.c
@@ -290,9 +290,16 @@ mpc_mul_karatsuba (mpc_ptr rop, mpc_srcptr op1, mpc_srcptr op2, mpc_rnd_t rnd)
 	 if (inexact == 0)
 	   {
 	     mp_prec_t prec_x;
-	     prec_x = (MPFR_EXP(v) > MPFR_EXP(w)) ? MPFR_EXP(v) - MPFR_EXP(w)
-	       : MPFR_EXP(w) - MPFR_EXP(v);
-	     prec_x += MAX(prec_v, prec_w) + 1;
+             if (MPFR_IS_ZERO(v))
+               prec_x = prec_w;
+             else if (MPFR_IS_ZERO(w))
+               prec_x = prec_v;
+             else
+               {
+                 prec_x = (MPFR_EXP(v) > MPFR_EXP(w)) ? MPFR_EXP(v) - MPFR_EXP(w)
+                   : MPFR_EXP(w) - MPFR_EXP(v);
+                 prec_x += MAX(prec_v, prec_w) + 1;
+               }
            /* +1 is necessary for a potential carry */
 	     /* ensure we do not use a too large precision */
 	     if (prec_x > prec_u)
diff --git a/norm.c b/norm.c
index ab7b303..bdb3f7a 100644
--- a/norm.c
+++ b/norm.c
@@ -36,7 +36,8 @@ mpc_norm (mpfr_ptr a, mpc_srcptr b, mp_rnd_t rnd)
   mpfr_init (u);
   mpfr_init (v);
 
-  if (2 * SAFE_ABS (mp_exp_t, MPFR_EXP (MPC_RE (b)) - MPFR_EXP (MPC_IM (b)))
+  if (!MPFR_IS_ZERO(MPC_RE(b)) && !MPFR_IS_ZERO(MPC_IM(b)) &&
+      2 * SAFE_ABS (mp_exp_t, MPFR_EXP (MPC_RE (b)) - MPFR_EXP (MPC_IM (b)))
        > (unsigned int) prec)
     /* If real and imaginary part have very different magnitudes, then the */
     /* generic code increases the precision too much. Instead, compute the */
diff --git a/tdiv.c b/tdiv.c
index f3ba637..d93f894 100644
--- a/tdiv.c
+++ b/tdiv.c
@@ -78,7 +78,20 @@ mpc_div_ref (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd)
       inex_re = mpfr_mul (u, MPC_RE(b), MPC_RE(c), rnd_re); /* 1 */
       inex_re |= mpfr_mul (v, MPC_IM(b), MPC_IM(c), rnd_re); /* 1 */
       inex_re |= mpfr_add (t, u, v, rnd_re);
-      cancel = MAX(MPFR_EXP(u), MPFR_EXP(v)) - MPFR_EXP(t);
+      if (MPFR_IS_ZERO(u))
+        {
+          if (MPFR_IS_ZERO(v))
+            cancel = 0;
+          else
+            cancel = MPFR_EXP(v) - MPFR_EXP(t);
+        }
+      else if (MPFR_IS_ZERO(v))
+        cancel = MPFR_EXP(u) - MPFR_EXP(t);
+      else
+        {
+          cancel = (MPFR_IS_ZERO(t)) ? prec
+            : MAX(MPFR_EXP(u), MPFR_EXP(v)) - MPFR_EXP(t);
+        }
       /* err(t) <= [1 + 2*2^cancel] ulp(t) */
       inex_re |= mpfr_mul (t, t, q, rnd_re) || inexact_q;
       /* err(t) <= [1 + 3*(1 + 2*2^cancel) + 14] ulp(t)
@@ -98,9 +111,20 @@ mpc_div_ref (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd)
         {
           inex_im = mpfr_mul (u, MPC_RE(b), MPC_IM(c), INV_RND(rnd_im));
           inex_im |= mpfr_mul (v, MPC_IM(b), MPC_RE(c), rnd_im);
-          cancel = MAX(MPFR_EXP(u), MPFR_EXP(v));
+          if (MPFR_IS_ZERO(u))
+            {
+              if (MPFR_IS_ZERO(v))
+                cancel = 0;
+              else
+                cancel = MPFR_EXP(v);
+            }
+          else if (MPFR_IS_ZERO(v))
+            cancel = MPFR_EXP(u);
+          else
+            cancel = MAX(MPFR_EXP(u), MPFR_EXP(v));
           inex_im |= mpfr_sub (u, v, u, rnd_im);
-          cancel = cancel - MPFR_EXP(u);
+          if (!MPFR_IS_ZERO(u))
+            cancel = cancel - MPFR_EXP(u);
           inex_im |= mpfr_mul (u, u, q, rnd_im) || inexact_q;
           err = (cancel <= 1) ? 5 : ((cancel == 2) ? 6 :
                                      ((cancel <= 4) ? 7 : cancel + 3));
@@ -135,7 +159,18 @@ mpc_div_ref (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd)
       inex_re = mpfr_mul (u, MPC_RE(b), MPC_RE(c), rnd_re); /* 1 */
       inex_re |= mpfr_mul (v, MPC_IM(b), MPC_IM(c), rnd_re); /* 1 */
       inex_re |= mpfr_add (t, u, v, rnd_re);
-      cancel = MAX(MPFR_EXP(u), MPFR_EXP(v)) - MPFR_EXP(t);
+      if (MPFR_IS_ZERO(u))
+        {
+          if (MPFR_IS_ZERO(v))
+            cancel = 0;
+          else
+            cancel = MPFR_EXP(v) - MPFR_EXP(t);
+        }
+      else if (MPFR_IS_ZERO(v))
+        cancel = MPFR_EXP(u) - MPFR_EXP(t);
+      else
+        cancel = (MPFR_IS_ZERO(t)) ? prec
+          : MAX(MPFR_EXP(u), MPFR_EXP(v)) - MPFR_EXP(t);
       /* err(t) <= [1 + 2*2^cancel] ulp(t) */
       inex_re |= mpfr_div (t, t, q, rnd_re) || inexact_q;
       /* err(t) <= [1 + 2*(1 + 2*2^cancel) + 6] ulp(t)
@@ -154,9 +189,20 @@ mpc_div_ref (mpc_ptr a, mpc_srcptr b, mpc_srcptr c, mpc_rnd_t rnd)
         {
           inex_im = mpfr_mul (u, MPC_RE(b), MPC_IM(c), INV_RND(rnd_im));
           inex_im |= mpfr_mul (v, MPC_IM(b), MPC_RE(c), rnd_im);
-          cancel = MAX(MPFR_EXP(u), MPFR_EXP(v));
+          if (MPFR_IS_ZERO(u))
+            {
+              if (MPFR_IS_ZERO(v))
+                cancel = 0;
+              else
+                cancel = MPFR_EXP(v);
+            }
+          else if (MPFR_IS_ZERO(v))
+            cancel = MPFR_EXP(u);
+          else
+            cancel = MAX(MPFR_EXP(u), MPFR_EXP(v));
           inex_im |= mpfr_sub (u, v, u, rnd_im);
-          cancel = cancel - MPFR_EXP(u);
+          if (!MPFR_IS_ZERO(u))
+            cancel = cancel - MPFR_EXP(u);
           inex_im |= mpfr_div (u, u, q, rnd_im) || inexact_q;
           err = (cancel <= 0) ? 4 : cancel + 3;
           ok = (inex_im == 0) || ((err < prec) && mpfr_can_round (u,