[Gmp-commit] /var/hg/gmp-5.1: Compute inverse as floor(B^2/(dh+1)), per Niels...

[mercurial at gmplib.org]

details:   /var/hg/gmp-5.1/rev/292ef78028cd
changeset: 15436:292ef78028cd
user:      Torbjorn Granlund <tege at gmplib.org>
date:      Wed Sep 25 19:54:41 2013 +0200
description:
Compute inverse as floor(B^2/(dh+1)), per Niels' suggestion.

diffstat:

 .bootstrap                  |   4 ++++
 ChangeLog                   |   6 ++++++
 NEWS                        |  18 ++++++++++++++++++
 mpn/generic/sb_div_sec.c    |  15 +++------------
 mpn/generic/sbpi1_div_sec.c |   7 ++-----
 5 files changed, 33 insertions(+), 17 deletions(-)

diffs (114 lines):

diff -r a06cfe94884b -r 292ef78028cd .bootstrap
--- a/.bootstrap	Fri Sep 20 20:40:15 2013 +0200
+++ b/.bootstrap	Wed Sep 25 19:54:41 2013 +0200
@@ -10,6 +10,10 @@
 # script.
 aclocal && libtoolize && autoconf && autoheader && automake
 
+cp -L ltmain.sh foo
+rm ltmain.sh
+mv foo ltmain.sh
+
 cat >doc/version.texi <<EOF
 @set UPDATED 19 January 2038
 @set UPDATED-MONTH January 2038
diff -r a06cfe94884b -r 292ef78028cd ChangeLog
--- a/ChangeLog	Fri Sep 20 20:40:15 2013 +0200
+++ b/ChangeLog	Wed Sep 25 19:54:41 2013 +0200
@@ -2,6 +2,12 @@
 
 	* doc/gmp.texi: Declare countless of function arguments as 'const'.
 
+2013-07-15  Torbjorn Granlund  <tege at gmplib.org>
+
+	* mpn/generic/sb_div_sec.c: Compute inverse as floor(B^2/(dh+1)), per
+	Niels' suggestion.
+	* mpn/generic/sbpi1_div_sec.c: Remove inverse rounding-up code.
+
 2013-07-12  Torbjorn Granlund  <tege at gmplib.org>
 
 	* mpn/generic/sbpi1_div_sec.c: Partial rewrite.
diff -r a06cfe94884b -r 292ef78028cd NEWS
--- a/NEWS	Fri Sep 20 20:40:15 2013 +0200
+++ b/NEWS	Wed Sep 25 19:54:41 2013 +0200
@@ -8,6 +8,24 @@
 Changes between GMP version 5.1.1 and 5.1.2
 
   BUGS FIXED
+  * The internal functions mpn_sbpi1_div_qr_sec mpn_sbpi1_div_r_sec could
+    compute garbage with a low probability.  They are now rewritten, and the
+    test code has been improved.
+
+  * The documentation now correctly says 'const' for input arguments.
+
+  SPEEDUPS
+  * None.
+
+  FEATURES
+  * None.
+
+  MISC
+  * ?
+
+Changes between GMP version 5.1.1 and 5.1.2
+
+  BUGS FIXED
   * A bug in mpz_powm_ui triggered by base arguments of at least 15000 decimal
     digits or mod arguments of at least 7500 decimal digits has been fixed.
 
diff -r a06cfe94884b -r 292ef78028cd mpn/generic/sb_div_sec.c
--- a/mpn/generic/sb_div_sec.c	Fri Sep 20 20:40:15 2013 +0200
+++ b/mpn/generic/sb_div_sec.c	Wed Sep 25 19:54:41 2013 +0200
@@ -81,18 +81,9 @@
       np2 = np;
     }
 
-  if (dn == 1)
-    {
-      d0 = dp2[dn - 1];
-      invert_limb (inv32, d0);
-    }
-  else
-    {
-      d1 = dp2[dn - 1];
-      d0 = dp2[dn - 2];
-      invert_pi1 (dinv, d1, d0);
-      inv32 = dinv.inv32;
-    }
+  d0 = dp2[dn - 1];
+  d0 += (~d0 != 0);
+  invert_limb (inv32, d0);
 
   /* We add nn + dn to tp here, not nn + 1 + dn, as expected.  This is since nn
      here will have been incremented.  */
diff -r a06cfe94884b -r 292ef78028cd mpn/generic/sbpi1_div_sec.c
--- a/mpn/generic/sbpi1_div_sec.c	Fri Sep 20 20:40:15 2013 +0200
+++ b/mpn/generic/sbpi1_div_sec.c	Wed Sep 25 19:54:41 2013 +0200
@@ -36,7 +36,7 @@
    too large (B is the limb base, D is the divisor, and i is the induction
    variable); the subsequent step will handle the extra partial remainder bits.
 
-   WIth that partial remainder reduction, each step generates a quotient "half
+   With that partial remainder reduction, each step generates a quotient "half
    limb".  The outer loop generates two quotient half limbs, an upper (q1h) and
    a lower (q0h) which are stored sparsely in separate limb arrays.  These
    arrays are added at the end; using separate arrays avoids data-dependent
@@ -48,7 +48,7 @@
    remainders, which we reduce later, as described above.
 
    In order to keep quotients from getting too big, corresponding to a negative
-   partial remainder, we use an inverse which is sligtly smaller than usually.
+   partial remainder, we use an inverse which is slightly smaller than usually.
 */
 
 #if OPERATION_sbpi1_div_qr_sec
@@ -94,9 +94,6 @@
 #endif
     }
 
-  /* Decremenet inverse to keep quotient half limbs from being too large.  */
-  dinv -= dinv != 0;				/* FIXME: cmp-to-int */
-
   /* Create a divisor copy shifted half a limb.  */
   hp = tp;					/* (dn + 1) limbs */
   hp[dn] = mpn_lshift (hp, dp, dn, GMP_NUMB_BITS / 2);