[Gmp-commit] /home/hgfiles/gmp: Reorganized and generalized mpn_gcd_subdiv_step.
mercurial at gmplib.org
mercurial at gmplib.org
Thu Apr 29 16:59:55 CEST 2010
details: /home/hgfiles/gmp/rev/bae441964b30
changeset: 13578:bae441964b30
user: Niels M?ller <nisse at lysator.liu.se>
date: Thu Apr 29 16:53:38 2010 +0200
description:
Reorganized and generalized mpn_gcd_subdiv_step.
diffstat:
ChangeLog | 28 +++++
configure.in | 2 +-
gmp-impl.h | 29 ++++-
mpn/generic/gcd.c | 42 ++++++-
mpn/generic/gcd_subdiv_step.c | 133 +++++++++++++++++---------
mpn/generic/gcdext.c | 37 ++++--
mpn/generic/gcdext_lehmer.c | 123 +++++++++++++++++++++++-
mpn/generic/gcdext_subdiv_step.c | 197 ---------------------------------------
8 files changed, 313 insertions(+), 278 deletions(-)
diffs (truncated from 825 to 300 lines):
diff -r 2c730f91277b -r bae441964b30 ChangeLog
--- a/ChangeLog Wed Apr 28 16:06:07 2010 +0200
+++ b/ChangeLog Thu Apr 29 16:53:38 2010 +0200
@@ -1,3 +1,31 @@
+2010-04-29 Niels Möller <nisse at lysator.liu.se>
+
+ * configure.in (gmp_mpn_functions): Deleted gcdext_subdiv_step.
+
+ * mpn/generic/gcdext.c (mpn_gcdext): Use new generalized
+ mpn_gcd_subdiv_step.
+
+ * mpn/generic/gcdext_lehmer.c (gcdext_update): New function.
+ (gcdext_done): New function.
+ (gcdext_hook): New const hook struct.
+ (mpn_gcdext_lehmer_n): Use new generalized mpn_gcd_subdiv_step.
+
+ * mpn/generic/gcd.c (gcd_done): New function.
+ (gcd_hook): New const hook struct.
+ (mpn_gcd): Adapted to new mpn_gcd_subdiv_step interface.
+
+ * mpn/generic/gcd_subdiv_step.c (mpn_gcd_subdiv_step): Reorganized
+ function. Added hook function pointers to the argument list, so
+ the same function can be used for gcd, gcdext, and jacobi.
+
+ * gmp-impl.h (struct gcd_subdiv_step_hook): New struct.
+ (mpn_gcdext_subdiv_step): Deleted prototype.
+ (struct gcdext_ctx): New struct.
+ (gcdext_hook): Declare const struct.
+ (mpn_gcd_subdiv_step): Updated prototype.
+
+ * mpn/generic/gcdext_subdiv_step.c: Deleted file.
+
2010-04-28 Niels Möller <nisse at lysator.liu.se>
* tests/mpz/t-jac.c (check_data): Added some more test cases.
diff -r 2c730f91277b -r bae441964b30 configure.in
--- a/configure.in Wed Apr 28 16:06:07 2010 +0200
+++ b/configure.in Thu Apr 29 16:53:38 2010 +0200
@@ -2514,7 +2514,7 @@
rootrem sqrtrem get_str set_str scan0 scan1 popcount hamdist cmp \
perfsqr perfpow \
gcd_1 gcd gcdext_1 gcdext gcd_subdiv_step \
- gcdext_lehmer gcdext_subdiv_step \
+ gcdext_lehmer \
div_q tdiv_qr jacbase jacobi_lehmer get_d \
matrix22_mul matrix22_mul1_inverse_vector \
hgcd2 hgcd mullo_n mullo_basecase \
diff -r 2c730f91277b -r bae441964b30 gmp-impl.h
--- a/gmp-impl.h Wed Apr 28 16:06:07 2010 +0200
+++ b/gmp-impl.h Thu Apr 29 16:53:38 2010 +0200
@@ -3859,14 +3859,35 @@
#define mpn_hgcd __MPN (hgcd)
__GMP_DECLSPEC mp_size_t mpn_hgcd __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, struct hgcd_matrix *, mp_ptr));
+struct gcd_subdiv_step_hook
+{
+ /* Passes quotient. */
+ void (*update)(void *, mp_srcptr, mp_size_t, unsigned);
+ /* Passes final gcd (must be copied if it is to be retained). */
+ void (*done)(void *, mp_srcptr, mp_size_t, unsigned);
+};
/* Needs storage for the quotient */
#define MPN_GCD_SUBDIV_STEP_ITCH(n) (n)
#define mpn_gcd_subdiv_step __MPN(gcd_subdiv_step)
-__GMP_DECLSPEC mp_size_t mpn_gcd_subdiv_step __GMP_PROTO ((mp_ptr, mp_size_t *, mp_ptr, mp_ptr, mp_size_t, mp_ptr));
-
-#define mpn_gcdext_subdiv_step __MPN(gcdext_subdiv_step)
-__GMP_DECLSPEC mp_size_t mpn_gcdext_subdiv_step __GMP_PROTO ((mp_ptr, mp_size_t *, mp_ptr, mp_size_t *, mp_ptr, mp_ptr, mp_size_t, mp_ptr, mp_ptr, mp_size_t *, mp_ptr, mp_ptr));
+__GMP_DECLSPEC mp_size_t mpn_gcd_subdiv_step __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, const struct gcd_subdiv_step_hook *, void *, mp_ptr));
+
+struct gcdext_ctx
+{
+ /* Result parameters. */
+ mp_ptr gp;
+ mp_size_t gn;
+ mp_ptr up;
+ mp_size_t *usize;
+
+ /* Cofactors updated in each step. */
+ mp_size_t un;
+ mp_ptr u0, u1, tp;
+};
+
+#define gcdext_hook __gmp_gcdext_hook
+extern const struct gcd_subdiv_step_hook
+gcdext_hook;
#define MPN_GCDEXT_LEHMER_N_ITCH(n) (4*(n) + 3)
diff -r 2c730f91277b -r bae441964b30 mpn/generic/gcd.c
--- a/mpn/generic/gcd.c Wed Apr 28 16:06:07 2010 +0200
+++ b/mpn/generic/gcd.c Thu Apr 29 16:53:38 2010 +0200
@@ -18,6 +18,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/. */
+#include <stdlib.h> /* for NULL */
+
#include "gmp.h"
#include "gmp-impl.h"
#include "longlong.h"
@@ -51,6 +53,26 @@
#define CHOOSE_P(n) (2*(n) / 3)
#endif
+struct gcd_ctx
+{
+ mp_ptr gp;
+ mp_size_t gn;
+};
+
+static void
+gcd_done (void *p, mp_srcptr gp, mp_size_t gn, unsigned swapped)
+{
+ struct gcd_ctx *ctx = (struct gcd_ctx *) p;
+ MPN_COPY (ctx->gp, gp, gn);
+ ctx->gn = gn;
+}
+
+static const struct gcd_subdiv_step_hook
+gcd_hook = {
+ NULL,
+ gcd_done,
+};
+
#if GMP_NAIL_BITS > 0
/* Nail supports should be easy, replacing the sub_ddmmss with nails
* logic. */
@@ -113,7 +135,7 @@
mp_size_t scratch;
mp_size_t matrix_scratch;
- mp_size_t gn;
+ struct gcd_ctx ctx;
mp_ptr tp;
TMP_DECL;
@@ -162,11 +184,13 @@
if (mpn_zero_p (up, n))
{
MPN_COPY (gp, vp, n);
- gn = n;
+ ctx.gn = n;
goto done;
}
}
+ ctx.gp = gp;
+
#if TUNE_GCD_P
while (CHOOSE_P (n) > 0)
#else
@@ -189,7 +213,7 @@
else
{
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (gp, &gn, up, vp, n, tp);
+ n = mpn_gcd_subdiv_step (up, vp, n, &gcd_hook, &ctx, tp);
if (n == 0)
goto done;
}
@@ -233,16 +257,18 @@
subtraction followed by one division. */
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (gp, &gn, up, vp, n, tp);
+ n = mpn_gcd_subdiv_step (up, vp, n, &gcd_hook, &ctx, tp);
if (n == 0)
goto done;
}
}
+ ASSERT(up[n-1] | vp[n-1]);
+
if (n == 1)
{
*gp = mpn_gcd_1(up, 1, vp[0]);
- gn = 1;
+ ctx.gn = 1;
goto done;
}
@@ -257,7 +283,7 @@
if (vp[0] == 0)
{
*gp = mpn_gcd_1 (up, 2, vp[1]);
- gn = 1;
+ ctx.gn = 1;
goto done;
}
else if (! (vp[0] & 1))
@@ -268,11 +294,11 @@
vp[1] >>= r;
}
- gn = gcd_2(gp, up, vp);
+ ctx.gn = gcd_2(gp, up, vp);
done:
TMP_FREE;
- return gn;
+ return ctx.gn;
}
#ifdef TUNE_GCD_P
diff -r 2c730f91277b -r bae441964b30 mpn/generic/gcd_subdiv_step.c
--- a/mpn/generic/gcd_subdiv_step.c Wed Apr 28 16:06:07 2010 +0200
+++ b/mpn/generic/gcd_subdiv_step.c Thu Apr 29 16:53:38 2010 +0200
@@ -4,7 +4,7 @@
SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES. IN FACT, IT IS ALMOST
GUARANTEED THAT THEY'LL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE.
-Copyright 2003, 2004, 2005, 2008 Free Software Foundation, Inc.
+Copyright 2003, 2004, 2005, 2008, 2010 Free Software Foundation, Inc.
This file is part of the GNU MP Library.
@@ -27,18 +27,32 @@
/* Used when mpn_hgcd or mpn_hgcd2 has failed. Then either one of a or
b is small, or the difference is small. Perform one subtraction
- followed by one division. If the gcd is found, stores it in gp and
- *gn, and returns zero. Otherwise, compute the reduced a and b, and
- return the new size. */
+ followed by one division. Returns zero if the gcd is found.
+ Otherwise, compute the reduced a and b, and return the new size. */
-/* FIXME: Check when the smaller number is a single limb, and invoke
- * mpn_gcd_1. */
+/* About the hook functions.
+
+ hook->update (ctx, qp, qn, d)
+
+ is called when one of the numbers, or a multiple of it, is
+ subtracted from the other. d == 0 means q a is subtracted from b, d
+ == 1 means that q b is subtracted from a.
+
+ hook->done (ctx, gp, gn, d)
+
+ is called when the gcd is found. d == 0 if the last reduction
+ subtracted a from b, d == 1 if it subtracted b from a, and d == 2
+ if it is unknown which was the most recent reduction. */
+
mp_size_t
-mpn_gcd_subdiv_step (mp_ptr gp, mp_size_t *gn,
- mp_ptr ap, mp_ptr bp, mp_size_t n, mp_ptr tp)
+mpn_gcd_subdiv_step (mp_ptr ap, mp_ptr bp, mp_size_t n,
+ const struct gcd_subdiv_step_hook *hook, void *ctx,
+ mp_ptr tp)
{
mp_size_t an, bn;
+ int swapped;
+
ASSERT (n > 0);
ASSERT (ap[n-1] > 0 || bp[n-1] > 0);
@@ -46,59 +60,84 @@
MPN_NORMALIZE (ap, an);
MPN_NORMALIZE (bp, bn);
- if (UNLIKELY (an == 0))
- {
- return_b:
- MPN_COPY (gp, bp, bn);
- *gn = bn;
- return 0;
- }
- else if (UNLIKELY (bn == 0))
- {
- return_a:
- MPN_COPY (gp, ap, an);
- *gn = an;
- return 0;
- }
+ swapped = 0;
- /* Arrange so that a > b, subtract an -= bn, and maintain
+ /* Arrange so that a < b, subtract b -= a, and maintain
normalization. */
- if (an < bn)
- MPN_PTR_SWAP (ap, an, bp, bn);
- else if (an == bn)
+ if (an == bn)
{
int c;
MPN_CMP (c, ap, bp, an);
if (UNLIKELY (c == 0))
- goto return_a;
- else if (c < 0)
- MP_PTR_SWAP (ap, bp);
+ {
+ /* For gcdext, return the smallest of the two cofactors. */
+ hook->done (ctx, ap, an, 2);
+ return 0;
+ }
+ else if (c > 0)
+ {
+ MP_PTR_SWAP (ap, bp);
More information about the gmp-commit
mailing list