[Gmp-commit] /var/hg/gmp-proj/mini-gmp: 6 new changesets

mercurial at gmplib.org mercurial at gmplib.org
Mon Jan 2 21:59:22 CET 2012


details:   /var/hg/gmp-proj/mini-gmp/rev/740aac61d72e
changeset: 32:740aac61d72e
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 10:57:53 2012 +0100
description:
Simplified mpz_gcd. New function mpz_make_odd.

details:   /var/hg/gmp-proj/mini-gmp/rev/85318940180e
changeset: 33:85318940180e
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 21:47:49 2012 +0100
description:
Use .PRECIOUS to keep object files.

details:   /var/hg/gmp-proj/mini-gmp/rev/fd2aaafed644
changeset: 34:fd2aaafed644
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 21:49:43 2012 +0100
description:
Define mpz_odd_p and mpz_even_p. Declare mpz_divexact and
mpz_divexact_ui as functions, deleting corresponding macros.

details:   /var/hg/gmp-proj/mini-gmp/rev/7f24526bad18
changeset: 35:7f24526bad18
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 21:54:55 2012 +0100
description:
Implemented mpz_gcdext.

Added MP_BITCNT_T_SWAP, fixed MPZ_PTR_SWAP and MPZ_SRCPTR_SWAP.

Simplified mpz_abs. Defined mpz_neg.

Implemented mpz_divexact and mpz_divexact_ui as functions.

Internal function mpz_divexact_2exp.

details:   /var/hg/gmp-proj/mini-gmp/rev/44e6a8fa316e
changeset: 36:44e6a8fa316e
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 21:55:37 2012 +0100
description:
More gcd testing, including gcdext.

details:   /var/hg/gmp-proj/mini-gmp/rev/4bdf488fa1a3
changeset: 37:4bdf488fa1a3
user:      Niels M?ller <nisse at lysator.liu.se>
date:      Mon Jan 02 21:58:47 2012 +0100
description:
Deleted gcdext debug code. Some whitespace cleanup.

diffstat:

 mini-gmp.c     |  465 +++++++++++++++++++++++++++++++++++++++-----------------
 mini-gmp.h     |    8 +-
 tests/Makefile |    3 +
 tests/t-gcd.c  |  135 +++++++++++++++-
 4 files changed, 456 insertions(+), 155 deletions(-)

diffs (truncated from 833 to 300 lines):

diff -r c2aa8afe583f -r 4bdf488fa1a3 mini-gmp.c
--- a/mini-gmp.c	Sun Jan 01 23:23:55 2012 +0100
+++ b/mini-gmp.c	Mon Jan 02 21:58:47 2012 +0100
@@ -34,16 +34,12 @@
      mpz_divexact_ui
      mpz_divisible_p
      mpz_divisible_ui_p
-     mpz_even_p
      mpz_export
-     mpz_gcd
-     mpz_gcd_ui
      mpz_getlimbn
      mpz_hamdist
      mpz_ior
      mpz_lcm_ui
      mpz_neg
-     mpz_odd_p
      mpz_popcount
      mpz_powm
      mpz_sqrtrem
@@ -72,6 +68,7 @@
 
 #define ABS(x) ((x) >= 0 ? (x) : -(x))
 
+#define MIN(a, b) ((a) < (b) ? (a) : (b))
 #define MAX(a, b) ((a) > (b) ? (a) : (b))
 
 #define assert_nocarry(x) do { \
@@ -204,7 +201,12 @@
     (x) = (y);								\
     (y) = __mp_size_t_swap__tmp;					\
   } while (0)
-
+#define MP_BITCNT_T_SWAP(x,y)			\
+  do {						\
+    mp_bitcnt_t __mp_bitcnt_t_swap__tmp = (x);	\
+    (x) = (y);					\
+    (y) = __mp_bitcnt_t_swap__tmp;		\
+  } while (0)
 #define MP_PTR_SWAP(x, y)						\
   do {									\
     mp_ptr __mp_ptr_swap__tmp = (x);					\
@@ -231,13 +233,13 @@
 
 #define MPZ_PTR_SWAP(x, y)						\
   do {									\
-    mpz_ptr __mpz_ptr_swap__tmp = (x);					\
+    __mpz_struct *__mpz_ptr_swap__tmp = (x);				\
     (x) = (y);								\
     (y) = __mpz_ptr_swap__tmp;						\
   } while (0)
 #define MPZ_SRCPTR_SWAP(x, y)						\
   do {									\
-    mpz_srcptr __mpz_srcptr_swap__tmp = (x);				\
+    const __mpz_struct *__mpz_srcptr_swap__tmp = (x);			\
     (x) = (y);								\
     (y) = __mpz_srcptr_swap__tmp;					\
   } while (0)
@@ -886,7 +888,7 @@
 	  assert_nocarry (mpn_lshift (tp, dp, dn, shift));
 	  nh = mpn_lshift (np, np, nn, shift);
 	  dp = tp;
-	  d1 = dp[dn - 1];	  
+	  d1 = dp[dn - 1];
 	}
       else
 	nh = 0;
@@ -1304,18 +1306,22 @@
 		   v->_mp_d, ABS(v->_mp_size));
 }
 
-void mpz_abs (mpz_t r, const mpz_t u)
+void
+mpz_abs (mpz_t r, const mpz_t u)
 {
-  mp_size_t n = ABS (u->_mp_size);
-
   if (r != u)
-    {
-      mp_ptr rp;
-      rp = MPZ_REALLOC (r, n);
-      mpn_copyi (rp, u->_mp_d, n);
-    }
-
-  r->_mp_size = n;
+    mpz_set (r, u);
+
+  r->_mp_size = ABS (r->_mp_size);
+}
+
+void
+mpz_neg (mpz_t r, const mpz_t u)
+{
+  if (r != u)
+    mpz_set (r, u);
+
+  r->_mp_size = -r->_mp_size;
 }
 
 void
@@ -1568,7 +1574,12 @@
   mp_ptr rp;
 
   un = ABS (u->_mp_size);
-  
+  if (un == 0)
+    {
+      r->_mp_size = 0;
+      return;
+    }
+
   limbs = bits / GMP_LIMB_BITS;
   shift = bits % GMP_LIMB_BITS;
 
@@ -1589,7 +1600,7 @@
   r->_mp_size = (u->_mp_size < 0) ? - rn : rn;
 }
 
-enum div_round_mode { DIV_FLOOR, DIV_CEIL, DIV_TRUNC };
+enum div_round_mode { DIV_FLOOR, DIV_CEIL, DIV_TRUNC, DIV_EXACT };
 
 /* Allows q or r to be zero. */
 static void
@@ -1619,6 +1630,8 @@
 
   if (nn < dn)
     {
+      assert (mode != DIV_EXACT);
+
       if (mode == DIV_CEIL && qs >= 0)
 	{
 	  /* q = 1, r = n - d */
@@ -1673,6 +1686,9 @@
 	  tq->_mp_size = qs < 0 ? -qn : qn;
 	}
       rn = mpn_normalized_size (np, dn);
+
+      assert (mode != DIV_EXACT || rn == 0);
+
       if (r)
 	{
 	  mpn_copyi (tr->_mp_d, np, rn);
@@ -1759,6 +1775,12 @@
   mpz_div_qr (NULL, r, n, d, DIV_TRUNC);
 }
 
+void
+mpz_divexact (mpz_t q, const mpz_t n, const mpz_t d)
+{
+  mpz_div_qr (q, NULL, n, d, DIV_EXACT);
+}
+
 unsigned long
 mpz_div_qr_ui (mpz_t q, mpz_t r,
 	       const mpz_t n, unsigned long d, enum div_round_mode mode)
@@ -1866,6 +1888,12 @@
   return mpz_div_qr_ui (NULL, NULL, n, d, DIV_TRUNC);
 }
 
+void
+mpz_divexact_ui (mpz_t q, const mpz_t n, unsigned long d)
+{
+  assert_nocarry (mpz_tdiv_q_ui (q, n, d));
+}
+
 static mp_limb_t
 limb_gcd (mp_limb_t u, mp_limb_t v)
 {
@@ -1927,7 +1955,7 @@
       return v;
     }
 
-  return limb_gcd (mpn_div_qr_1 (NULL, u->_mp_d, un, v), v); 
+  return limb_gcd (mpn_div_qr_1 (NULL, u->_mp_d, un, v), v);
 }
 
 #define MAKE_ODD(xp, n) do {					\
@@ -1942,156 +1970,309 @@
       }								\
   } while (0)
 
+static mp_bitcnt_t
+mpz_make_odd (mpz_t r, const mpz_t u)
+{
+  mp_size_t un, rn, i;
+  mp_ptr rp;
+  unsigned shift;
+
+  un = ABS (u->_mp_size);
+  assert (un > 0);
+
+  for (i = 0; u->_mp_d[i] == 0; i++)
+    ;
+
+  count_trailing_zeros (shift, u->_mp_d[i]);
+
+  rn = un - i;
+  rp = MPZ_REALLOC (r, rn);
+  if (shift > 0)
+    {
+      mpn_rshift (rp, u->_mp_d + i, rn, shift);
+      rn -= (rp[rn-1] == 0);
+    }
+  else
+    mpn_copyi (rp, u->_mp_d + i, rn);
+
+  r->_mp_size = rn;
+  return i * GMP_LIMB_BITS + shift;
+}
+
 void
 mpz_gcd (mpz_t g, const mpz_t u, const mpz_t v)
 {
-  mp_srcptr orig_up, orig_vp;
-  mp_ptr tp, up, vp;  
-  mp_size_t un, vn, zero_limbs, alloc;
-  mp_limb_t ul, vl;
-  unsigned vshift, zero_bits;
-
-  un = ABS (u->_mp_size);
-  vn = ABS (v->_mp_size);
-
-  if (un == 0)
+  mpz_t tu, tv;
+  mp_bitcnt_t uz, vz, gz;
+
+  if (u->_mp_size == 0)
     {
       mpz_abs (g, v);
       return;
     }
-  if (vn == 0)
+  if (v->_mp_size == 0)
     {
       mpz_abs (g, u);
       return;
     }
 
-  orig_up = u->_mp_d;
-  orig_vp = v->_mp_d;
-
-  for (zero_limbs = 0; (orig_up[0] | orig_vp[0]) == 0; zero_limbs++)
+  mpz_init (tu);
+  mpz_init (tv);
+
+  uz = mpz_make_odd (tu, u);
+  vz = mpz_make_odd (tv, v);
+  gz = MIN(uz, vz);
+
+  if (tu->_mp_size < tv->_mp_size)
+    mpz_swap (tu, tv);
+
+  mpz_tdiv_r (tu, tu, tv);
+  if (tu->_mp_size == 0)
     {
-      orig_up++; un--;
-      orig_vp++; vn--;
+      mpz_swap (g, tv);
     }
-
-  count_trailing_zeros (zero_bits, orig_up[0] | orig_vp[0]);
-
-  while (orig_up[0] == 0)
-    {
-      orig_up++; un--;
-    }
-
-  while (orig_vp[0] == 0)
-    {
-      orig_vp++; vn--;
-    }
-      
-  if (un < vn)
-    MPN_SRCPTR_SWAP (orig_up, un, orig_vp, vn);
-
-  count_trailing_zeros (vshift, orig_vp[0]);
-  
-  if (vn == 1)
-    {
-      vl = orig_vp[0] >> vshift;
-      ul = mpn_div_qr_1 (NULL, orig_up, un, vl);
-      goto gcd_1;
-    }
-  else if (vn == 2 && ( (orig_vp[1] >> vshift) == 0))
-    {
-      vl = (orig_vp[0] >> vshift) | (orig_vp[1] << (GMP_LIMB_BITS - vshift));
-      ul = mpn_div_qr_1 (NULL, orig_up, un, vl);
-      goto gcd_1;
-    }
-
-  alloc = MAX (un, 2*vn);
-      
-  tp = xalloc_limbs (alloc);
-  up = tp;
-
-  mpn_copyi (up, orig_up, un);
-  
-  mpn_div_qr (NULL, up, un, orig_vp, vn);
-  un = mpn_normalized_size (up, vn);
+  else


More information about the gmp-commit mailing list