[gcalctool] Simplify modifications of MP.t



commit 91c392759f25a96a240e349df38972d7eb5c4334
Author: Robert Ancell <robert ancell gmail com>
Date:   Sun May 10 18:37:48 2009 +1000

    Simplify modifications of MP.t
---
 src/mp-trigonometric.c |   24 +++++++++++++-----------
 src/mp.c               |   34 +++++++++++++++++-----------------
 2 files changed, 30 insertions(+), 28 deletions(-)

diff --git a/src/mp-trigonometric.c b/src/mp-trigonometric.c
index 05a9d9a..b103b17 100644
--- a/src/mp-trigonometric.c
+++ b/src/mp-trigonometric.c
@@ -60,7 +60,7 @@ mp_compare_mp_to_int(const MPNumber *x, int i)
 static void
 mpsin1(const MPNumber *x, MPNumber *z, int do_sin)
 {
-    int i, b2, ts;
+    int i, b2;
     MPNumber t1, t2;
 
     mpchk();
@@ -86,19 +86,19 @@ mpsin1(const MPNumber *x, MPNumber *z, int do_sin)
 
     z->sign = 0;
     i = 1;
-    ts = MP.t;
     if (do_sin != 0) {
         mp_set_from_mp(&t1, z);
         i = 2;
     }
 
     /* POWER SERIES LOOP.  REDUCE T IF POSSIBLE */
-    do {
-        MP.t = t1.exponent + ts + 2;
-        if (MP.t <= 2)
+    for (; ; i+= 2) {
+        int t, ts;
+        
+        t = MP.t + t1.exponent + 2;
+        if (t <= 2)
             break;
-
-        MP.t = min(MP.t,ts);
+        t = min(t, MP.t);
 
         /* PUT R(I3) FIRST IN CASE ITS DIGITS ARE MAINLY ZERO */
         mp_multiply(&t2, &t1, &t1);
@@ -106,19 +106,21 @@ mpsin1(const MPNumber *x, MPNumber *z, int do_sin)
         /*  IF I*(I+1) IS NOT REPRESENTABLE AS AN INTEGER, THE FOLLOWING
          *  DIVISION BY I*(I+1) HAS TO BE SPLIT UP.
          */
+        ts = MP.t;
+        MP.t = t;
         if (i > b2) {
             mp_divide_integer(&t1, -i, &t1);
             mp_divide_integer(&t1, i + 1, &t1);
         } else {
             mp_divide_integer(&t1, -i * (i + 1), &t1);
         }
-
-        i += 2;
         MP.t = ts;
+
         mp_add(&t1, z, z);
-    } while(t1.sign != 0);
+        if (t1.sign == 0)
+            break;
+    }
 
-    MP.t = ts;
     if (do_sin == 0)
         mp_add_integer(z, 1, z);
 }
diff --git a/src/mp.c b/src/mp.c
index f85ca80..905420d 100644
--- a/src/mp.c
+++ b/src/mp.c
@@ -939,7 +939,7 @@ mp_epowy(const MPNumber *x, MPNumber *z)
 void
 mpexp1(const MPNumber *x, MPNumber *y)
 {
-    int i, q, ib, ic, ts;
+    int i, q, ib, ic;
     float rlb;
     MPNumber t1, t2;
 
@@ -984,24 +984,27 @@ mpexp1(const MPNumber *x, MPNumber *y)
     }
     mp_set_from_mp(&t1, y);
     mp_set_from_mp(&t1, &t2);
-    i = 1;
-    ts = MP.t;
 
     /* SUM SERIES, REDUCING T WHERE POSSIBLE */
-    do {
-        MP.t = ts + 2 + t2.exponent - y->exponent;
-        if (MP.t <= 2)
+    for (i = 2; ; i++)  {
+        int t, ts;
+        
+        t = MP.t + 2 + t2.exponent - y->exponent;
+        if (t <= 2)
             break;
+        t = min(t, MP.t);
 
-        MP.t = min(MP.t,ts);
+        ts = MP.t;
+        MP.t = t;
         mp_multiply(&t1, &t2, &t2);
-        ++i;
         mp_divide_integer(&t2, i, &t2);
         MP.t = ts;
+
         mp_add(&t2, y, y);
-    } while (t2.sign != 0);
+        if (t2.sign == 0)
+            break;
+    }
 
-    MP.t = ts;
     if (q <= 0)
         return;
 
@@ -1241,8 +1244,7 @@ mplns(const MPNumber *x, MPNumber *y)
         return;
     }
 
-    /* SAVE T AND GET STARTING APPROXIMATION TO -LN(1+X) */
-    ts = MP.t;
+    /* GET STARTING APPROXIMATION TO -LN(1+X) */
     mp_set_from_mp(x, &t2);
     mp_divide_integer(x, 4, &t1);
     mp_add_fraction(&t1, -1, 3, &t1);
@@ -1255,6 +1257,7 @@ mplns(const MPNumber *x, MPNumber *y)
     /* START NEWTON ITERATION USING SMALL T, LATER INCREASE */
 
     /* Computing MAX */
+    ts = MP.t;
     MP.t = max(5, 13 - (MP.b << 1));
     if (MP.t <= ts)
     {
@@ -1276,12 +1279,10 @@ mplns(const MPNumber *x, MPNumber *y)
              */
             ts3 = MP.t;
             MP.t = ts;
-
             do {
                 ts2 = MP.t;
                 MP.t = (MP.t + it0) / 2;
             } while (MP.t > ts3);
-
             MP.t = ts2;
         }
         
@@ -1290,10 +1291,10 @@ mplns(const MPNumber *x, MPNumber *y)
             mperr("*** ERROR OCCURRED IN MPLNS, NEWTON ITERATION NOT CONVERGING PROPERLY ***");
         }
     }
+    MP.t = ts;
 
     /* REVERSE SIGN OF Y AND RETURN */
     y->sign = -y->sign;
-    MP.t = ts;
 }
 
 
@@ -2082,9 +2083,8 @@ mp_root(const MPNumber *x, int n, MPNumber *z)
             mperr("*** ERROR OCCURRED IN MP_ROOT, NEWTON ITERATION NOT CONVERGING PROPERLY ***");
         }
     }
-
-    /* RESTORE T */
     MP.t = ts;
+
     if (n >= 0) {
         mp_pwr_integer(&t1, n - 1, &t1);
         mp_multiply(x, &t1, z);



[Date Prev][Date Next]   [Thread Prev][Thread Next]   [Thread Index] [Date Index] [Author Index]