[Gmp-commit] /home/hgfiles/gmp: 2 new changesets

mercurial at gmplib.org mercurial at gmplib.org
Mon Dec 21 15:17:29 CET 2009


details:   /home/hgfiles/gmp/rev/fa4b6ce53583
changeset: 13156:fa4b6ce53583
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Dec 21 14:55:55 2009 +0100
description:
New function and test for Toom-6.5

details:   /home/hgfiles/gmp/rev/62c33581e638
changeset: 13157:62c33581e638
user:      Marco Bodrato <bodrato at mail.dm.unipi.it>
date:      Mon Dec 21 15:17:18 2009 +0100
description:
New shared helper function for high degree Toom.

diffstat:

 ChangeLog                            |   14 +
 configure.in                         |    4 +-
 gmp-impl.h                           |   12 +
 mpn/generic/toom63_mul.c             |   35 --
 mpn/generic/toom6h_mul.c             |  501 +++++++++++++++++++++++++++++++++++
 mpn/generic/toom_couple_handling.c   |   70 ++++
 mpn/generic/toom_interpolate_12pts.c |  350 ++++++++++++++++++++++++
 tests/mpn/Makefile.am                |    2 +-
 8 files changed, 951 insertions(+), 37 deletions(-)

diffs (truncated from 1061 to 300 lines):

diff -r 4036bcf04357 -r 62c33581e638 ChangeLog
--- a/ChangeLog	Mon Dec 21 14:23:46 2009 +0100
+++ b/ChangeLog	Mon Dec 21 15:17:18 2009 +0100
@@ -21,6 +21,20 @@
 
 2009-12-21  Marco Bodrato <bodrato at mail.dm.unipi.it>
 
+	* mpn/generic/toom_couple_handling.c: New file, helper function
+	for high degree Toom.
+	* gmp-impl.h: Provide corresponding declaration.
+	* configure.in (gmp_mpn_functions): List toom_couple_handling.
+	* mpn/generic/toom6h_mul.c: Use shared toom_couple_handling.
+	* mpn/generic/toom63_mul.c: Likewise.
+
+	* mpn/generic/toom6h_mul.c: New file.
+	* mpn/generic/toom_interpolate_12pts.c: New file.
+	* gmp-impl.h: Provide corresponding declarations.
+	* configure.in (gmp_mpn_functions): List toom_interpolate_12pts,
+	toom6h_mul.
+	* tests/mpn/t-toom6h.c: New test program.
+
 	* tests/mpn/t-mulmod_bnm1.c (ref_mulmod_bnm1): Use ref_mul.
 	* tests/mpn/t-sqrmod_bnm1.c (ref_sqrmod_bnm1): Likewise.
 
diff -r 4036bcf04357 -r 62c33581e638 configure.in
--- a/configure.in	Mon Dec 21 14:23:46 2009 +0100
+++ b/configure.in	Mon Dec 21 15:17:18 2009 +0100
@@ -2506,11 +2506,13 @@
   toom22_mul toom32_mul toom42_mul toom52_mul toom62_mul		   \
   toom33_mul toom43_mul toom53_mul toom63_mul				   \
   toom44_mul								   \
+  toom6h_mul								   \
+  toom_couple_handling							   \
   toom2_sqr toom3_sqr toom4_sqr						   \
   toom_eval_dgr3_pm1 toom_eval_dgr3_pm2 				   \
   toom_eval_pm1 toom_eval_pm2 toom_eval_pm2exp	   			   \
   toom_interpolate_5pts toom_interpolate_6pts toom_interpolate_7pts	   \
-  toom_interpolate_8pts							   \
+  toom_interpolate_8pts toom_interpolate_12pts				   \
   invertappr invert binvert mulmod_bnm1 sqrmod_bnm1			   \
   sbpi1_div_q sbpi1_div_qr sbpi1_divappr_q				   \
   dcpi1_div_q dcpi1_div_qr dcpi1_divappr_q				   \
diff -r 4036bcf04357 -r 62c33581e638 gmp-impl.h
--- a/gmp-impl.h	Mon Dec 21 14:23:46 2009 +0100
+++ b/gmp-impl.h	Mon Dec 21 15:17:18 2009 +0100
@@ -1071,6 +1071,12 @@
 #define mpn_toom_interpolate_8pts __MPN(toom_interpolate_8pts)
 __GMP_DECLSPEC void      mpn_toom_interpolate_8pts __GMP_PROTO ((mp_ptr, mp_size_t, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
 
+#define mpn_toom_interpolate_12pts __MPN(toom_interpolate_12pts)
+__GMP_DECLSPEC void      mpn_toom_interpolate_12pts __GMP_PROTO ((mp_ptr, mp_ptr, mp_ptr, mp_ptr, mp_size_t, mp_size_t, int, mp_ptr));
+
+#define   mpn_toom_couple_handling __MPN(toom_couple_handling)
+__GMP_DECLSPEC void      toom_couple_handling __GMP_PROTO ((mp_ptr, mp_size_t, mp_ptr, int, mp_size_t, int, int));
+
 #define   mpn_toom_eval_dgr3_pm1 __MPN(toom_eval_dgr3_pm1)
 __GMP_DECLSPEC int mpn_toom_eval_dgr3_pm1 __GMP_PROTO ((mp_ptr, mp_ptr, mp_srcptr, mp_size_t, mp_size_t, mp_ptr));
 
@@ -1125,6 +1131,12 @@
 #define   mpn_toom4_sqr __MPN(toom4_sqr)
 __GMP_DECLSPEC void      mpn_toom4_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
 
+#define   mpn_toom6h_mul __MPN(toom6h_mul)
+__GMP_DECLSPEC void      mpn_toom6h_mul __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
+
+#define   mpn_toom6h_sqr __MPN(toom6h_sqr)
+__GMP_DECLSPEC void      mpn_toom6h_sqr __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
+
 #define   mpn_fft_best_k __MPN(fft_best_k)
 __GMP_DECLSPEC int       mpn_fft_best_k __GMP_PROTO ((mp_size_t, int)) ATTRIBUTE_CONST;
 
diff -r 4036bcf04357 -r 62c33581e638 mpn/generic/toom63_mul.c
--- a/mpn/generic/toom63_mul.c	Mon Dec 21 14:23:46 2009 +0100
+++ b/mpn/generic/toom63_mul.c	Mon Dec 21 15:17:18 2009 +0100
@@ -59,41 +59,6 @@
   return 0;
 }
 
-static void
-toom_couple_handling (mp_ptr pp, mp_size_t n, mp_ptr np, int nsign, mp_size_t off, int ps, int ns)
-{
-  if (nsign) {
-#ifdef HAVE_NATIVE_mpn_rsh1sub_n
-    mpn_rsh1sub_n (np, pp, np, n);
-#else
-    mpn_sub_n (np, pp, np, n);
-    mpn_rshift (np, np, n, 1);
-#endif
-  } else {
-#ifdef HAVE_NATIVE_mpn_rsh1add_n
-    mpn_rsh1add_n (np, pp, np, n);
-#else
-    mpn_add_n (np, pp, np, n);
-    mpn_rshift (np, np, n, 1);
-#endif
-  }
-
-#ifdef HAVE_NATIVE_mpn_rsh1sub_n
-  if (ps == 1)
-    mpn_rsh1sub_n (pp, pp, np, n);
-  else
-#endif
-  {
-    mpn_sub_n (pp, pp, np, n);
-    if (ps > 0)
-      mpn_rshift (pp, pp, n, ps);
-  }
-  if (ns > 0)
-    mpn_rshift (np, np, n, ns);
-  pp[n-1] += mpn_add_n (pp+off, pp+off, np, n>>1);
-  ASSERT_NOCARRY( mpn_add_1(pp+n-1, np+off, n-off, pp[n-1]) );
-}
-
 static int
 abs_sub_add_n (mp_ptr rm, mp_ptr rp, mp_srcptr rs, mp_size_t n) {
   int result;
diff -r 4036bcf04357 -r 62c33581e638 mpn/generic/toom6h_mul.c
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mpn/generic/toom6h_mul.c	Mon Dec 21 15:17:18 2009 +0100
@@ -0,0 +1,501 @@
+/* Implementation of the algorithm for Toom-Cook 6.5-way.
+
+   Contributed to the GNU project by Marco Bodrato.
+
+   THE FUNCTION IN THIS FILE IS INTERNAL WITH A MUTABLE INTERFACE.  IT IS ONLY
+   SAFE TO REACH IT THROUGH DOCUMENTED INTERFACES.  IN FACT, IT IS ALMOST
+   GUARANTEED THAT IT WILL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE.
+
+Copyright 2009 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation; either version 3 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+License for more details.
+
+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/.  */
+
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+
+#if HAVE_NATIVE_mpn_sublsh_n
+#define DO_mpn_sublsh_n(dst,src,n,s,ws) mpn_sublsh_n(dst,dst,src,n,s)
+#else
+static mp_limb_t
+DO_mpn_sublsh_n(mp_ptr dst, mp_srcptr src, mp_size_t n, unsigned int s, mp_ptr ws)
+{
+#if USE_MUL_1 && 0
+  return mpn_submul_1(dst,src,n,CNST_LIMB(1) <<(s));
+#else
+  mp_limb_t __cy;
+  __cy = mpn_lshift(ws,src,n,s);
+  return    __cy + mpn_sub_n(dst,dst,ws,n);
+#endif
+}
+#endif
+
+#if HAVE_NATIVE_mpn_addlsh_n
+#define DO_mpn_addlsh_n(dst,src,n,s,ws) mpn_addlsh_n(dst,dst,src,n,s)
+#else
+static mp_limb_t
+DO_mpn_addlsh_n(mp_ptr dst, mp_srcptr src, mp_size_t n, unsigned int s, mp_ptr ws)
+{
+#if USE_MUL_1 && 0
+  return mpn_addmul_1(dst,src,n,CNST_LIMB(1) <<(s));
+#else
+  mp_limb_t __cy;
+  __cy = mpn_lshift(ws,src,n,s);
+  return    __cy + mpn_add_n(dst,dst,ws,n);
+#endif
+}
+#endif
+
+#if HAVE_NATIVE_mpn_subrsh
+#define DO_mpn_subrsh(dst,nd,src,ns,s,ws) mpn_subrsh(dst,nd,src,ns,s)
+#else
+/* FIXME: This is not a correct definition, it assumes no carry */
+#define DO_mpn_subrsh(dst,nd,src,ns,s,ws)				\
+do {									\
+  mp_limb_t __cy;							\
+  MPN_DECR_U (dst, nd, src[0] >> s);					\
+  __cy = DO_mpn_sublsh_n (dst, src + 1, ns - 1, GMP_NUMB_BITS - s, ws);	\
+  MPN_DECR_U (dst + ns - 1, nd - ns + 1, __cy);				\
+} while (0)
+#endif
+
+
+/* Stores |{ap,n}-{bp,n}| in {rp,n}, returns the sign. */
+static int
+abs_sub_n (mp_ptr rp, mp_srcptr ap, mp_srcptr bp, mp_size_t n)
+{
+  mp_limb_t  x, y;
+  while (--n >= 0)
+    {
+      x = ap[n];
+      y = bp[n];
+      if (x != y)
+	{
+	  n++;
+	  if (x > y)
+	    {
+	      mpn_sub_n (rp, ap, bp, n);
+	      return 0;
+	    }
+	  else
+	    {
+	      mpn_sub_n (rp, bp, ap, n);
+	      return 1;
+	    }
+	}
+      rp[n] = 0;
+    }
+  return 0;
+}
+
+static int
+abs_sub_add_n (mp_ptr rm, mp_ptr rp, mp_srcptr rs, mp_size_t n) {
+  int result;
+  result = abs_sub_n (rm, rp, rs, n);
+  ASSERT_NOCARRY(mpn_add_n (rp, rp, rs, n));
+  return result;
+}
+
+static int
+mpn_toom_ev_rsh(mp_ptr rp, mp_ptr rm,
+		mp_srcptr ap, unsigned int q, mp_size_t n, mp_size_t t,
+		unsigned int s, mp_ptr ws, mp_ptr scratch)
+{
+  unsigned int i;
+  /* {ap,q*n+t} -> {rp,n+1} {rm,n+1} , with {ws, n+1}*/
+  ASSERT( n >= t );
+  ASSERT( s != 0 ); /* or _ev_pm1 should be used */
+  ASSERT( q > 1 );
+  ASSERT( s*q < GMP_NUMB_BITS );
+  rp[n] = mpn_lshift(rp, ap, n, s*q);
+  ws[n] = mpn_lshift(ws, ap+n, n, s*(q-1));
+  if( (q & 1) != 0) {
+    ASSERT_NOCARRY(mpn_add(ws,ws,n+1,ap+n*q,t));
+    rp[n] += DO_mpn_addlsh_n(rp, ap+n*(q-1), n, s, scratch);
+  } else {
+    ASSERT_NOCARRY(mpn_add(rp,rp,n+1,ap+n*q,t));
+  }
+  for(i=2; i<q-1; i++)
+  {
+    rp[n] += DO_mpn_addlsh_n(rp, ap+n*i, n, s*(q-i), scratch);
+    i++;
+    ws[n] += DO_mpn_addlsh_n(ws, ap+n*i, n, s*(q-i), scratch);
+  };
+  return abs_sub_add_n (rm, rp, ws, n + 1);
+}
+
+
+#if GMP_NUMB_BITS < 21
+#error Not implemented.
+#endif
+
+
+/* FIXME: tuneup should decide the threshold */
+#ifndef MUL_TOOM6H_THRESHOLD
+#define MUL_TOOM6H_THRESHOLD 256
+#endif
+
+#if TUNE_PROGRAM_BUILD
+#define MAYBE_mul_basecase 1
+#define MAYBE_mul_toom22   1
+#define MAYBE_mul_toom33   1
+#define MAYBE_mul_toom6h   1
+#else
+#define MAYBE_mul_basecase						\
+  (MUL_TOOM6H_THRESHOLD < 6 * MUL_TOOM22_THRESHOLD)
+#define MAYBE_mul_toom22						\
+  (MUL_TOOM6H_THRESHOLD < 6 * MUL_TOOM33_THRESHOLD)
+#define MAYBE_mul_toom33						\
+  (MUL_TOOM6H_THRESHOLD < 6 * MUL_TOOM44_THRESHOLD)
+#define MAYBE_mul_toom6h						\
+  (MUL_FFT_THRESHOLD >= 6 * MUL_TOOM6H_THRESHOLD)
+#endif
+
+#define TOOM6H_MUL_N_REC(p, a, b, n, ws)				\
+  do {									\
+    if (MAYBE_mul_basecase						\
+	&& BELOW_THRESHOLD (n, MUL_TOOM22_THRESHOLD))			\
+      mpn_mul_basecase (p, a, n, b, n);					\
+    else if (MAYBE_mul_toom22						\
+	     && BELOW_THRESHOLD (n, MUL_TOOM33_THRESHOLD))		\
+      mpn_toom22_mul (p, a, n, b, n, ws);				\
+    else if (MAYBE_mul_toom33						\
+	     && BELOW_THRESHOLD (n, MUL_TOOM44_THRESHOLD))		\
+      mpn_toom33_mul (p, a, n, b, n, ws);				\
+    else if (! MAYBE_mul_toom6h						\
+	     || BELOW_THRESHOLD (n, MUL_TOOM6H_THRESHOLD))		\


More information about the gmp-commit mailing list