[Gmp-commit] /var/hg/gmp: 2 new changesets
mercurial at gmplib.org
mercurial at gmplib.org
Thu May 19 20:56:47 CEST 2011
details: /var/hg/gmp/rev/dc7e918260d8
changeset: 14186:dc7e918260d8
user: Niels M?ller <nisse at lysator.liu.se>
date: Thu May 19 20:47:26 2011 +0200
description:
New argument to mpn_gcd_subdiv_step.
details: /var/hg/gmp/rev/f901b9fcfdd3
changeset: 14187:f901b9fcfdd3
user: Niels M?ller <nisse at lysator.liu.se>
date: Thu May 19 20:55:23 2011 +0200
description:
Additional fixes for mpn_gcd_subdiv_step change.
diffstat:
ChangeLog | 12 ++++++
gmp-impl.h | 2 +-
mpn/generic/gcd.c | 4 +-
mpn/generic/gcd_subdiv_step.c | 84 +++++++++++++++++++++++++++++++-----------
mpn/generic/gcdext.c | 6 ++-
mpn/generic/gcdext_lehmer.c | 5 ++-
mpn/generic/jacobi_lehmer.c | 2 +-
7 files changed, 86 insertions(+), 29 deletions(-)
diffs (263 lines):
diff -r 71efb0367192 -r f901b9fcfdd3 ChangeLog
--- a/ChangeLog Tue May 17 23:05:59 2011 +0200
+++ b/ChangeLog Thu May 19 20:55:23 2011 +0200
@@ -1,3 +1,15 @@
+2011-05-19 Niels Möller <nisse at lysator.liu.se>
+
+ * mpn/generic/gcd_subdiv_step.c (mpn_gcd_subdiv_step): New
+ argument s, for use by hgcd.
+ * gmp-impl.h (mpn_gcd_subdiv_step): Update declaration.
+
+ * mpn/generic/gcd.c (mpn_gcd): Pass s = 0 to mpn_gcd_subdiv_step.
+ * mpn/generic/gcdext.c (mpn_gcdext): Likewise. Also added an ASSERT.
+ * mpn/generic/gcdext_lehmer.c (mpn_gcdext_lehmer_n): Likewise.
+ (mpn_gcdext_hook): Added some ASSERTs.
+ * mpn/generic/jacobi_lehmer.c (mpn_jacobi_lehmer): Likewise.
+
2011-05-17 Niels Möller <nisse at lysator.liu.se>
* doc/gmp.texi (mpn_gcd, mpn_gcdext): Document input requirements:
diff -r 71efb0367192 -r f901b9fcfdd3 gmp-impl.h
--- a/gmp-impl.h Tue May 17 23:05:59 2011 +0200
+++ b/gmp-impl.h Thu May 19 20:55:23 2011 +0200
@@ -3944,7 +3944,7 @@
#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_ptr, mp_size_t, gcd_subdiv_step_hook *, void *, mp_ptr));
+__GMP_DECLSPEC mp_size_t mpn_gcd_subdiv_step __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_size_t, gcd_subdiv_step_hook *, void *, mp_ptr));
struct gcdext_ctx
{
diff -r 71efb0367192 -r f901b9fcfdd3 mpn/generic/gcd.c
--- a/mpn/generic/gcd.c Tue May 17 23:05:59 2011 +0200
+++ b/mpn/generic/gcd.c Thu May 19 20:55:23 2011 +0200
@@ -210,7 +210,7 @@
else
{
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (up, vp, n, gcd_hook, &ctx, tp);
+ n = mpn_gcd_subdiv_step (up, vp, n, 0, gcd_hook, &ctx, tp);
if (n == 0)
goto done;
}
@@ -254,7 +254,7 @@
subtraction followed by one division. */
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (up, vp, n, &gcd_hook, &ctx, tp);
+ n = mpn_gcd_subdiv_step (up, vp, n, 0, &gcd_hook, &ctx, tp);
if (n == 0)
goto done;
}
diff -r 71efb0367192 -r f901b9fcfdd3 mpn/generic/gcd_subdiv_step.c
--- a/mpn/generic/gcd_subdiv_step.c Tue May 17 23:05:59 2011 +0200
+++ b/mpn/generic/gcd_subdiv_step.c Thu May 19 20:55:23 2011 +0200
@@ -29,8 +29,14 @@
/* 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. Returns zero if the gcd is found.
- Otherwise, compute the reduced a and b, and return the new size. */
+ followed by one division. The normal case is to compute the reduced
+ a and b, and return the new size.
+
+ If s == 0 (used for gcd and gcdext), returns zero if the gcd is
+ found.
+
+ If s > 0, don't reduce to size <= s, and return zero if no
+ reduction is possible (if either a, b or |a-b| is of size <= s). */
/* The hook function is called as
@@ -38,15 +44,15 @@
in the following cases:
- + If A == B at the start, G is the gcd, Q is NULL, d = -1.
+ + If A = B at the start, G is the gcd, Q is NULL, d = -1.
+ If one input is zero at the start, G is the gcd, Q is NULL,
- d = 0 if A == G and d = 1 if B == G.
+ d = 0 if A = G and d = 1 if B = G.
Otherwise, if d = 0 we have just subtracted a multiple of A from B,
and if d = 1 we have subtracted a multiple of B from A.
- + If A == B after subtraction, G is the gcd, Q is NULL.
+ + If A = B after subtraction, G is the gcd, Q is NULL.
+ If we get a zero remainder after division, G is the gcd, Q is the
quotient.
@@ -56,7 +62,7 @@
*/
mp_size_t
-mpn_gcd_subdiv_step (mp_ptr ap, mp_ptr bp, mp_size_t n,
+mpn_gcd_subdiv_step (mp_ptr ap, mp_ptr bp, mp_size_t n, mp_size_t s,
gcd_subdiv_step_hook *hook, void *ctx,
mp_ptr tp)
{
@@ -82,8 +88,10 @@
MPN_CMP (c, ap, bp, an);
if (UNLIKELY (c == 0))
{
- /* For gcdext, return the smallest of the two cofactors. */
- hook (ctx, ap, an, NULL, 0, -1);
+ /* For gcdext, return the smallest of the two cofactors, so
+ pass d = -1. */
+ if (s == 0)
+ hook (ctx, ap, an, NULL, 0, -1);
return 0;
}
else if (c > 0)
@@ -99,16 +107,27 @@
MPN_PTR_SWAP (ap, an, bp, bn);
swapped ^= 1;
}
- if (an == 0)
- {
- hook (ctx, bp, bn, NULL, 0, swapped ^ 1);
- return 0;
- }
}
+ if (an <= s)
+ {
+ if (s == 0)
+ hook (ctx, bp, bn, NULL, 0, swapped ^ 1);
+ return 0;
+ }
+
ASSERT_NOCARRY (mpn_sub (bp, bp, bn, ap, an));
MPN_NORMALIZE (bp, bn);
ASSERT (bn > 0);
+ if (bn <= s)
+ {
+ /* Undo subtraction. */
+ mp_limb_t cy = mpn_add (bp, ap, an, bp, bn);
+ if (cy > 0)
+ bp[an] = cy;
+ return 0;
+ }
+
/* Arrange so that a < b */
if (an == bn)
{
@@ -116,7 +135,12 @@
MPN_CMP (c, ap, bp, an);
if (UNLIKELY (c == 0))
{
- hook (ctx, bp, bn, NULL, 0, swapped);
+ if (s > 0)
+ /* Just record subtraction and return */
+ hook (ctx, NULL, 0, &one, 1, swapped);
+ else
+ /* Found gcd. */
+ hook (ctx, bp, bn, NULL, 0, swapped);
return 0;
}
@@ -141,14 +165,30 @@
mpn_tdiv_qr (tp, bp, 0, bp, bn, ap, an);
qn = bn - an + 1;
- if (mpn_zero_p (bp, an))
+ bn = an;
+ MPN_NORMALIZE (bp, bn);
+
+ if (UNLIKELY (bn <= s))
{
- hook (ctx, ap, an, tp, qn, swapped);
- return 0;
+ if (s == 0)
+ {
+ hook (ctx, ap, an, tp, qn, swapped);
+ return 0;
+ }
+
+ /* Quotient is one too large, so decrement it and add back A. */
+ if (bn > 0)
+ {
+ mp_limb_t cy = mpn_add (bp, ap, an, bp, bn);
+ if (cy)
+ bp[an++] = cy;
+ }
+ else
+ MPN_COPY (bp, ap, an);
+
+ MPN_DECR_U (tp, qn, 1);
}
- else
- {
- hook (ctx, NULL, 0, tp, qn, swapped);
- return an;
- }
+
+ hook (ctx, NULL, 0, tp, qn, swapped);
+ return an;
}
diff -r 71efb0367192 -r f901b9fcfdd3 mpn/generic/gcdext.c
--- a/mpn/generic/gcdext.c Tue May 17 23:05:59 2011 +0200
+++ b/mpn/generic/gcdext.c Thu May 19 20:55:23 2011 +0200
@@ -319,7 +319,7 @@
ctx.un = 1;
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (ap, bp, n, mpn_gcdext_hook, &ctx, tp);
+ n = mpn_gcd_subdiv_step (ap, bp, n, 0, mpn_gcdext_hook, &ctx, tp);
if (n == 0)
{
TMP_FREE;
@@ -374,7 +374,7 @@
ctx.un = un;
/* Temporary storage n */
- n = mpn_gcd_subdiv_step (ap, bp, n, mpn_gcdext_hook, &ctx, tp);
+ n = mpn_gcd_subdiv_step (ap, bp, n, 0, mpn_gcdext_hook, &ctx, tp);
if (n == 0)
{
TMP_FREE;
@@ -386,6 +386,8 @@
}
}
+ ASSERT ( (ap[n-1] | bp[n-1]) > 0);
+
if (UNLIKELY (mpn_cmp (ap, bp, n) == 0))
{
/* Must return the smallest cofactor, +u1 or -u0 */
diff -r 71efb0367192 -r f901b9fcfdd3 mpn/generic/gcdext_lehmer.c
--- a/mpn/generic/gcdext_lehmer.c Tue May 17 23:05:59 2011 +0200
+++ b/mpn/generic/gcdext_lehmer.c Thu May 19 20:55:23 2011 +0200
@@ -35,6 +35,9 @@
{
mp_srcptr up;
+ ASSERT (gn > 0);
+ ASSERT (gp[gn-1] > 0);
+
MPN_COPY (ctx->gp, gp, gn);
ctx->gn = gn;
@@ -216,7 +219,7 @@
/* Temporary storage n for the quotient and ualloc for the
new cofactor. */
- n = mpn_gcd_subdiv_step (ap, bp, n, mpn_gcdext_hook, &ctx, tp);
+ n = mpn_gcd_subdiv_step (ap, bp, n, 0, mpn_gcdext_hook, &ctx, tp);
if (n == 0)
return ctx.gn;
diff -r 71efb0367192 -r f901b9fcfdd3 mpn/generic/jacobi_lehmer.c
--- a/mpn/generic/jacobi_lehmer.c Tue May 17 23:05:59 2011 +0200
+++ b/mpn/generic/jacobi_lehmer.c Thu May 19 20:55:23 2011 +0200
@@ -831,7 +831,7 @@
/* mpn_hgcd2 has failed. Then either one of a or b is very
small, or the difference is very small. Perform one
subtraction followed by one division. */
- n = mpn_gcd_subdiv_step (ap, bp, n, &jacobi_hook, &bits, tp);
+ n = mpn_gcd_subdiv_step (ap, bp, n, 0, &jacobi_hook, &bits, tp);
if (!n)
return bits == BITS_FAIL ? 0 : mpn_jacobi_finish (bits);
}
More information about the gmp-commit
mailing list