log10.c: shortened code without changing functionality

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

1 files changed, 37 insertions, 90 deletions

diff --git a/src/log10.c b/src/log10.c
index a4ff9af..0d35c76 100644
--- a/src/log10.c
+++ b/src/log10.c
@@ -21,19 +21,17 @@ along with this program. If not, see http://www.gnu.org/licenses/ .
 #include <limits.h> /* for CHAR_BIT */
 #include "mpc-impl.h"
 
-/* Auxiliary functions which implement Ziv's strategy for special cases.
-   if flag = 0: compute only real part
-   if flag = 1: compute only imaginary
+/* Auxiliary function which implements Ziv's strategy for special cases.
    Exact cases should be dealt with separately. */
 static int
-mpc_log10_aux (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd, int flag, int nb)
+mpc_log10_aux (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd, int nb)
 {
   mp_prec_t prec = (MPFR_PREC_MIN > 4) ? MPFR_PREC_MIN : 4;
   mpc_t tmp;
   mpfr_t log10;
   int ok = 0, ret;
 
-  prec = mpfr_get_prec ((flag == 0) ? mpc_realref (rop) : mpc_imagref (rop));
+  prec = mpfr_get_prec (mpc_imagref (rop));
   prec += 10;
   mpc_init2 (tmp, prec);
   mpfr_init2 (log10, prec);
@@ -44,40 +42,17 @@ mpc_log10_aux (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd, int flag, int nb)
       /* In each case we have two roundings, thus the final value is
          x * (1+u)^2 where x is the exact value, and |u| <= 2^(-prec-1).
          Thus the error is always less than 3 ulps. */
-      switch (nb)
-        {
-        case 0: /* imag <- atan2(y/x) */
-          mpfr_atan2 (mpc_imagref (tmp), mpc_imagref (op), mpc_realref (op),
-                      MPC_RND_IM (rnd));
-          mpfr_div (mpc_imagref (tmp), mpc_imagref (tmp), log10, MPFR_RNDN);
-          ok = mpfr_can_round (mpc_imagref (tmp), prec - 2, MPFR_RNDN,
-                               MPFR_RNDZ, MPC_PREC_IM(rop) +
-                               (MPC_RND_IM (rnd) == MPFR_RNDN));
-          if (ok)
-            ret = mpfr_set (mpc_imagref (rop), mpc_imagref (tmp),
-                            MPC_RND_IM (rnd));
-          break;
-        case 1: /* real <- log(x) */
-          mpfr_log (mpc_realref (tmp), mpc_realref (op), MPC_RND_RE (rnd));
-          mpfr_div (mpc_realref (tmp), mpc_realref (tmp), log10, MPFR_RNDN);
-          ok = mpfr_can_round (mpc_realref (tmp), prec - 2, MPFR_RNDN,
-                               MPFR_RNDZ, MPC_PREC_RE(rop) +
-                               (MPC_RND_RE (rnd) == MPFR_RNDN));
-          if (ok)
-            ret = mpfr_set (mpc_realref (rop), mpc_realref (tmp),
-                            MPC_RND_RE (rnd));
-          break;
-        case 2: /* imag <- pi */
-          mpfr_const_pi (mpc_imagref (tmp), MPC_RND_IM (rnd));
-          mpfr_div (mpc_imagref (tmp), mpc_imagref (tmp), log10, MPFR_RNDN);
-          ok = mpfr_can_round (mpc_imagref (tmp), prec - 2, MPFR_RNDN,
-                               MPFR_RNDZ, MPC_PREC_IM(rop) +
-                               (MPC_RND_IM (rnd) == MPFR_RNDN));
-          if (ok)
-            ret = mpfr_set (mpc_imagref (rop), mpc_imagref (tmp),
-                            MPC_RND_IM (rnd));
-          break;
-        }
+      if (nb == 0)
+        mpfr_atan2 (mpc_imagref (tmp), mpc_imagref (op), mpc_realref (op),
+          MPC_RND_IM (rnd));
+      else
+        mpfr_const_pi (mpc_imagref (tmp), MPC_RND_IM (rnd));
+      mpfr_div (mpc_imagref (tmp), mpc_imagref (tmp), log10, MPFR_RNDN);
+      ok = mpfr_can_round (mpc_imagref (tmp), prec - 2, MPFR_RNDN,
+             MPFR_RNDZ, MPC_PREC_IM(rop) + (MPC_RND_IM (rnd) == MPFR_RNDN));
+      if (ok)
+        ret = mpfr_set (mpc_imagref (rop), mpc_imagref (tmp),
+                        MPC_RND_IM (rnd));
       prec += prec / 2;
       mpc_set_prec (tmp, prec);
       mpfr_set_prec (log10, prec);
@@ -97,46 +72,22 @@ mpc_log10 (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd)
   mpc_t ww;
   mpc_rnd_t invrnd;
 
-  /* special values: NaN and infinities: same as mpc_log */
-  if (!mpc_fin_p (op)) /* real or imaginary parts are NaN or Inf */
-    {
-      if (mpfr_nan_p (mpc_realref (op)))
-        {
-          if (mpfr_inf_p (mpc_imagref (op)))
-            /* (NaN, Inf) -> (+Inf, NaN) */
-            mpfr_set_inf (mpc_realref (rop), +1);
-          else
-            /* (NaN, xxx) -> (NaN, NaN) */
-            mpfr_set_nan (mpc_realref (rop));
-         mpfr_set_nan (mpc_imagref (rop));
-         inex_im = 0; /* Inf/NaN is exact */
-        }
-      else if (mpfr_nan_p (mpc_imagref (op)))
-        {
-          if (mpfr_inf_p (mpc_realref (op)))
-            /* (Inf, NaN) -> (+Inf, NaN) */
-            mpfr_set_inf (mpc_realref (rop), +1);
-          else
-            /* (xxx, NaN) -> (NaN, NaN) */
-            mpfr_set_nan (mpc_realref (rop));
-          mpfr_set_nan (mpc_imagref (rop));
-          inex_im = 0; /* Inf/NaN is exact */
-        }
-      else /* We have an infinity in at least one part. */
-        {
-          /* (+Inf, y) -> (+Inf, 0) for finite positive-signed y */
-          if (mpfr_inf_p (mpc_realref (op)) && mpfr_signbit (mpc_realref (op))
-              == 0 && mpfr_number_p (mpc_imagref (op)))
-            inex_im = mpfr_atan2 (mpc_imagref (rop), mpc_imagref (op),
-                                  mpc_realref (op), MPC_RND_IM (rnd));
-          else
-            /* (xxx, Inf) -> (+Inf, atan2(Inf/xxx))
-               (Inf, yyy) -> (+Inf, atan2(yyy/Inf)) */
-            inex_im = mpc_log10_aux (rop, op, rnd, 1, 0);
-          mpfr_set_inf (mpc_realref (rop), +1);
-        }
-      return MPC_INEX(0, inex_im);
-    }
+   /* special values: NaN and infinities: same as mpc_log */
+   if (!mpc_fin_p (op)) {
+      if (  mpfr_nan_p (mpc_realref (op))
+          || mpfr_nan_p (mpc_imagref (op))
+          || (   mpfr_inf_p (mpc_realref (op))
+              && mpfr_signbit (mpc_realref (op))== 0
+              && mpfr_number_p (mpc_imagref (op))))
+         return mpc_log (rop, op, rnd);
+      else {
+         /* (xxx, Inf) -> (+Inf, atan2(Inf/xxx))
+            (Inf, yyy) -> (+Inf, atan2(yyy/Inf)) */
+         inex_im = mpc_log10_aux (rop, op, rnd, 0);
+         mpfr_set_inf (mpc_realref (rop), +1);
+         return MPC_INEX(0, inex_im);
+      }
+   }
 
    /* special cases: real and purely imaginary numbers */
   re_cmp = mpfr_cmp_ui (mpc_realref (op), 0);
@@ -149,7 +100,7 @@ mpc_log10 (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd)
             inex_im = mpfr_atan2 (mpc_imagref (rop), mpc_imagref (op),
                                   mpc_realref (op), MPC_RND_IM (rnd));
           else
-            inex_im = mpc_log10_aux (rop, op, rnd, 1, 0);
+            inex_im = mpc_log10_aux (rop, op, rnd, 0);
           mpfr_set_inf (mpc_realref (rop), -1);
           inex_re = 0; /* -Inf is exact */
         }
@@ -167,14 +118,10 @@ mpc_log10 (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd)
             rnd_im = INV_RND (MPC_RND_IM (rnd));
           else
             rnd_im = MPC_RND_IM (rnd);
-          ww->re[0] = *mpc_realref (op);
-          MPFR_CHANGE_SIGN (ww->re);
-          ww->im[0] = *mpc_imagref (op);
-          if (mpfr_cmp_ui (ww->re, 1) == 0)
-            inex_re = mpfr_set_ui (mpc_realref (rop), 0, MPC_RND_RE (rnd));
-          else
-            inex_re = mpc_log10_aux (rop, ww, rnd, 0, 1);
-          inex_im = mpc_log10_aux (rop, op, MPC_RND (0,rnd_im), 1, 2);
+          w [0] = *mpc_realref (op);
+          MPFR_CHANGE_SIGN (w);
+          inex_re = mpfr_log10 (mpc_realref (rop), w, MPC_RND_RE (rnd));
+          inex_im = mpc_log10_aux (rop, op, MPC_RND (0,rnd_im), 2);
           if (negative_zero)
             {
               mpc_conj (rop, rop, MPC_RNDNN);
@@ -188,7 +135,7 @@ mpc_log10 (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd)
        if (im_cmp > 0)
          {
            inex_re = mpfr_log10 (mpc_realref (rop), mpc_imagref (op), MPC_RND_RE (rnd));
-           inex_im = mpc_log10_aux (rop, op, rnd, 1, 2);
+           inex_im = mpc_log10_aux (rop, op, rnd, 2);
            /* division by 2 does not change the ternary flag */
            mpfr_div_2ui (mpc_imagref (rop), mpc_imagref (rop), 1, MPFR_RNDN);
          }
@@ -198,7 +145,7 @@ mpc_log10 (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd)
            MPFR_CHANGE_SIGN (w);
            inex_re = mpfr_log10 (mpc_realref (rop), w, MPC_RND_RE (rnd));
            invrnd = MPC_RND (0, INV_RND (MPC_RND_IM (rnd)));
-           inex_im = mpc_log10_aux (rop, op, invrnd, 1, 2);
+           inex_im = mpc_log10_aux (rop, op, invrnd, 2);
            /* division by 2 does not change the ternary flag */
            mpfr_div_2ui (mpc_imagref (rop), mpc_imagref (rop), 1, MPFR_RNDN);
            mpfr_neg (mpc_imagref (rop), mpc_imagref (rop), MPFR_RNDN);