[Gmp-commit] /home/hgfiles/gmp: Some comments and threshold handling.

[mercurial at gmplib.org]

details:   /home/hgfiles/gmp/rev/535b98be60b3
changeset: 13101:535b98be60b3
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Thu Dec 17 08:03:49 2009 +0100
description:
Some comments and threshold handling.

diffstat:

 ChangeLog                |   5 +++++
 mpn/generic/invert.c     |  15 ++++++++-------
 mpn/generic/invertappr.c |  14 ++++++++++++--
 3 files changed, 25 insertions(+), 9 deletions(-)

diffs (93 lines):

diff -r 3135d56c2945 -r 535b98be60b3 ChangeLog
--- a/ChangeLog	Thu Dec 17 00:45:43 2009 +0100
+++ b/ChangeLog	Thu Dec 17 08:03:49 2009 +0100
@@ -1,3 +1,8 @@
+2009-12-15  Marco Bodrato <bodrato at mail.dm.unipi.it>
+
+	* mpn/generic/invert.c: Added some comment.
+	* mpn/generic/invertappr.c: Slightly better threshold handling.
+
 2009-12-17  Torbjorn Granlund  <tege at gmplib.org>
 
 	* mpn/generic/bdiv_q.c (mpn_bdiv_q_itch): Rewrite.
diff -r 3135d56c2945 -r 535b98be60b3 mpn/generic/invert.c
--- a/mpn/generic/invert.c	Thu Dec 17 00:45:43 2009 +0100
+++ b/mpn/generic/invert.c	Thu Dec 17 08:03:49 2009 +0100
@@ -23,6 +23,8 @@
 You should have received a copy of the GNU Lesser General Public License
 along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
 
+/* FIXME: Remove NULL and TMP_*, as soon as all the callers properly
+   allocate and pass the scratch to the function. */
 #include <stdlib.h>		/* for NULL */
 
 #include "gmp.h"
@@ -33,7 +35,6 @@
 #define INV_APPR_THRESHOLD (INV_NEWTON_THRESHOLD)
 #endif
 
-
 void
 mpn_invert (mp_ptr ip, mp_srcptr dp, mp_size_t n, mp_ptr scratch)
 {
@@ -64,17 +65,17 @@
       mpn_tdiv_qr (scratch, ip, 0, xp, 2 * n, dp, n);
       MPN_COPY (ip, scratch, n);
     } else { /* Use approximated inverse; correct the result if needed. */
-      mp_limb_t cy;
+      mp_limb_t e; /* The possible error in the approximate inverse */
 
-      cy = mpn_ni_invertappr (ip, dp, n, scratch);
+      ASSERT ( mpn_invert_itch (n) >= mpn_invertappr_itch (n) )
+      e = mpn_ni_invertappr (ip, dp, n, scratch);
 
-      if (cy) {
+      if (e) { /* Assume the error can only be "0" (no error) or "1". */
 	/* Code to detect and correct the "off by one" approximation. */
 	mpn_mul_n (scratch, ip, dp, n);
 	ASSERT_NOCARRY (mpn_add_n (scratch + n, scratch + n, dp, n));
-	if ( ! mpn_add (scratch, scratch, 2*n, dp, n)) {
-	  MPN_INCR_U (ip, n, 1);
-	}
+	if (! mpn_add (scratch, scratch, 2*n, dp, n))
+	  MPN_INCR_U (ip, n, 1); /* The value was wrong, correct it.  */
       }
     }
     TMP_FREE;
diff -r 3135d56c2945 -r 535b98be60b3 mpn/generic/invertappr.c
--- a/mpn/generic/invertappr.c	Thu Dec 17 00:45:43 2009 +0100
+++ b/mpn/generic/invertappr.c	Thu Dec 17 08:03:49 2009 +0100
@@ -57,7 +57,12 @@
 #define NPOWS \
  ((sizeof(mp_size_t) > 6 ? 48 : 8*sizeof(mp_size_t)) - LOG2C (INV_NEWTON_THRESHOLD))
 #define MAYBE_dcpi1_divappr \
-  (INV_NEWTON_THRESHOLD < 2 * DC_DIVAPPR_Q_THRESHOLD)
+  (INV_NEWTON_THRESHOLD < DC_DIVAPPR_Q_THRESHOLD)
+#if (INV_NEWTON_THRESHOLD > INV_MULMOD_BNM1_THRESHOLD) && \
+    (INV_APPR_THRESHOLD > INV_MULMOD_BNM1_THRESHOLD)
+#undef  INV_MULMOD_BNM1_THRESHOLD
+#define INV_MULMOD_BNM1_THRESHOLD 0 /* always when Newton */
+#endif
 #endif
 
 /* All the three functions mpn{,_bc,_ni}_invertappr (ip, dp, n, scratch), take
@@ -93,13 +98,18 @@
 
   /* Compute a base value of r limbs. */
   if (n == 1)
-    invert_limb (*(ip),*(dp));
+    invert_limb (*ip, *dp);
   else {
     mp_size_t i;
     xp = tp + n + 2;				/* 2 * n limbs */
+
     for (i = n - 1; i >= 0; i--)
       xp[i] = ~CNST_LIMB(0);
     mpn_com_n (xp + n, dp, n);
+
+    /* Now xp contains B^2n - {dp,n}*B^n - 1 */
+
+    /* FIXME: if mpn_*pi1_divappr_q handles n==2, use it! */
     if (n == 2) {
       mpn_tdiv_qr (tp, ip, 0, xp, 2 * n, dp, n);
       MPN_COPY (ip, tp, n);