Clean up

Add ZivLoop


git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/trunk@3301 280ebfd0-de03-0410-8827-d642c229c3f4

1 files changed, 42 insertions, 43 deletions

diff --git a/expm1.c b/expm1.c
index 795287d5d..be33600b0 100644
--- a/expm1.c
+++ b/expm1.c
@@ -1,6 +1,6 @@
 /* mpfr_expm1 -- Compute exp(x)-1
 
-Copyright 2001, 2002, 2003, 2004 Free Software Foundation.
+Copyright 2001, 2002, 2003, 2004, 2005 Free Software Foundation.
 
 This file is part of the MPFR Library.
 
@@ -29,86 +29,85 @@ MA 02111-1307, USA. */
 int
 mpfr_expm1 (mpfr_ptr y, mpfr_srcptr x , mp_rnd_t rnd_mode) 
 {
-  int inexact = 0;
+  int inexact;
+  MPFR_SAVE_EXPO_DECL (expo);
 
-  if (MPFR_UNLIKELY( MPFR_IS_SINGULAR(x) ))
+  if (MPFR_UNLIKELY (MPFR_IS_SINGULAR (x)))
     {
-      if (MPFR_IS_NAN(x))
+      if (MPFR_IS_NAN (x))
 	{
-	  MPFR_SET_NAN(y);
+	  MPFR_SET_NAN (y);
 	  MPFR_RET_NAN;
 	}
       /* check for inf or -inf (expm1(-inf)=-1) */
-      else if (MPFR_IS_INF(x))
+      else if (MPFR_IS_INF (x))
 	{ 
-	  if (MPFR_IS_POS(x))
+	  if (MPFR_IS_POS (x))
 	    {
-	      MPFR_SET_INF(y);
-	      MPFR_SET_POS(y);
-	      return 0;
+	      MPFR_SET_INF (y);
+	      MPFR_SET_POS (y);
+	      MPFR_RET (0);
 	    }
 	  else
-	    return mpfr_set_si(y, -1, rnd_mode);
+	    return mpfr_set_si (y, -1, rnd_mode);
 	}
       else
 	{
-          MPFR_ASSERTD(MPFR_IS_ZERO(x));
-	  MPFR_SET_ZERO(y);   /* expm1(+/- 0) = +/- 0 */
-	  MPFR_SET_SAME_SIGN(y,x);
-	  MPFR_RET(0);
+          MPFR_ASSERTD (MPFR_IS_ZERO (x));
+	  MPFR_SET_ZERO (y);   /* expm1(+/- 0) = +/- 0 */
+	  MPFR_SET_SAME_SIGN (y, x);
+	  MPFR_RET (0);
 	}
     }
-  /* Useless due to mpfr_set
-     MPFR_CLEAR_FLAGS(y);*/
 
+  MPFR_SAVE_EXPO_MARK (expo);
   /* General case */
   {
     /* Declaration of the intermediary variable */
-    mpfr_t te,t;
-
+    mpfr_t t;
     /* Declaration of the size variable */
     mp_prec_t Nx = MPFR_PREC(x);   /* Precision of input variable */
-    mp_prec_t Ny = MPFR_PREC(y);   /* Precision of input variable */
-    
-    mp_prec_t Nt;   /* Precision of the intermediary variable */
-    long int err;  /* Precision of error */
-                
+    mp_prec_t Ny = MPFR_PREC(y);   /* Precision of input variable */    
+    mp_prec_t Nt;       /* Precision of the intermediary variable */
+    mp_exp_t err, exp_te;                   /* Precision of error */
+    MPFR_ZIV_DECL (loop);
+
     /* compute the precision of intermediary variable */
-    Nt=MAX(Nx,Ny);
+    Nt = MAX (Nx, Ny);
     /* the optimal number of bits : see algorithms.ps */
-    Nt=Nt+5+MPFR_INT_CEIL_LOG2(Nt);
+    Nt = Nt + 5 + MPFR_INT_CEIL_LOG2 (Nt);
 
     /* initialise of intermediary	variable */
-    mpfr_init(t);             
-    mpfr_init(te);             
+    mpfr_init2 (t, Nt);
 
     /* First computation of cosh */
-    do
+    MPFR_ZIV_INIT (loop, Nt);
+    for (;;)
       {
-
-        /* reactualisation of the precision */
-        mpfr_set_prec(t,Nt);
-        mpfr_set_prec(te,Nt);
-
         /* compute expm1 */
-        mpfr_exp(te,x,GMP_RNDN);        /* exp(x)*/
-        mpfr_sub_ui(t,te,1,GMP_RNDN);   /* exp(x)-1 */
+        mpfr_exp (t, x, GMP_RNDN);         /* exp(x)*/
+	exp_te = MPFR_GET_EXP (t);         /* FIXME: exp(x) may overflow! */
+        mpfr_sub_ui (t, t, 1, GMP_RNDN);   /* exp(x)-1 */
 
         /* error estimate */
         /*err=Nt-(__gmpfr_ceil_log2(1+pow(2,MPFR_EXP(te)-MPFR_EXP(t))));*/
-        err = Nt - (MAX (MPFR_GET_EXP (te) - MPFR_GET_EXP (t), 0) + 1);
+        err = Nt - (MAX (exp_te - MPFR_GET_EXP (t), 0) + 1);
+
+	if (MPFR_LIKELY (mpfr_can_round (t, err, GMP_RNDN, GMP_RNDZ,
+                                         Ny + (rnd_mode == GMP_RNDN))))
+	  break;
 
         /* actualisation of the precision */
-        Nt += 10;
+	MPFR_ZIV_NEXT (loop, Nt);
+        mpfr_set_prec (t, Nt);
       }
-    while ((err < 0) || !mpfr_can_round (t, err, GMP_RNDN, GMP_RNDZ,
-                                         Ny + (rnd_mode == GMP_RNDN)));
- 
+    MPFR_ZIV_FREE (loop);
+
     inexact = mpfr_set (y, t, rnd_mode);
 
     mpfr_clear (t);
-    mpfr_clear (te);
   }
 
-  return inexact;
+  MPFR_SAVE_EXPO_FREE (expo);
+  return mpfr_check_range (y, inexact, rnd_mode);
 }