[Gmp-commit] /home/hgfiles/gmp: 4 new changesets
mercurial at gmplib.org
mercurial at gmplib.org
Tue Dec 22 22:32:16 CET 2009
details: /home/hgfiles/gmp/rev/807e9d2eff12
changeset: 13186:807e9d2eff12
user: Torbjorn Granlund <tege at gmplib.org>
date: Tue Dec 22 22:20:58 2009 +0100
description:
(mpn_mu_div_q_itch): New function.
details: /home/hgfiles/gmp/rev/30b41ab5301b
changeset: 13187:30b41ab5301b
user: Torbjorn Granlund <tege at gmplib.org>
date: Tue Dec 22 22:28:43 2009 +0100
description:
Handle quotient overflow.
details: /home/hgfiles/gmp/rev/eb5530e19ae2
changeset: 13188:eb5530e19ae2
user: Torbjorn Granlund <tege at gmplib.org>
date: Tue Dec 22 22:30:53 2009 +0100
description:
New test.
details: /home/hgfiles/gmp/rev/c795ab883a97
changeset: 13189:c795ab883a97
user: Torbjorn Granlund <tege at gmplib.org>
date: Tue Dec 22 22:31:42 2009 +0100
description:
Call tests_end().
diffstat:
ChangeLog | 9 +
gmp-impl.h | 2 +
mpn/generic/mu_div_q.c | 10 +
mpn/generic/mu_divappr_q.c | 14 +-
tests/mpn/Makefile.am | 2 +-
tests/mpn/t-bdiv.c | 2 +
tests/mpn/t-div.c | 397 +++++++++++++++++++++++++++++++++++++++++++++
tests/mpn/t-hgcd.c | 7 +
tests/mpn/t-invert.c | 1 +
tests/mpn/t-matrix22.c | 1 +
tests/mpn/t-mullo.c | 1 +
tests/mpn/t-mulmod_bnm1.c | 1 +
tests/mpn/t-sqrmod_bnm1.c | 1 +
tests/mpn/toom-shared.h | 2 +
tests/mpz/t-remove.c | 1 +
15 files changed, 444 insertions(+), 7 deletions(-)
diffs (truncated from 624 to 300 lines):
diff -r 6eb2023eb8f4 -r c795ab883a97 ChangeLog
--- a/ChangeLog Tue Dec 22 14:15:35 2009 +0100
+++ b/ChangeLog Tue Dec 22 22:31:42 2009 +0100
@@ -1,3 +1,12 @@
+2009-12-22 Torbjorn Granlund <tege at gmplib.org>
+
+ * tests/mpn/t-div.c: New file.
+ * tests/mpn/Makefile.am: Compile it.
+
+ * mpn/generic/mu_divappr_q.c: Handle quotient overflow.
+
+ * mpn/generic/mu_div_q.c (mpn_mu_div_q_itch): New function.
+
2009-12-22 Niels Möller <<nisse at lysator.liu.se>>
* mpn/generic/sbpi1_div_q.c (mpn_sbpi1_div_q): Use udiv_qr_3by2.
diff -r 6eb2023eb8f4 -r c795ab883a97 gmp-impl.h
--- a/gmp-impl.h Tue Dec 22 14:15:35 2009 +0100
+++ b/gmp-impl.h Tue Dec 22 22:31:42 2009 +0100
@@ -1202,6 +1202,8 @@
#define mpn_mu_div_q __MPN(mu_div_q)
__GMP_DECLSPEC mp_limb_t mpn_mu_div_q __GMP_PROTO ((mp_ptr, mp_ptr, mp_size_t, mp_srcptr, mp_size_t, mp_ptr));
+#define mpn_mu_div_q_itch __MPN(mu_div_q_itch)
+__GMP_DECLSPEC mp_size_t mpn_mu_div_q_itch __GMP_PROTO ((mp_size_t, mp_size_t, int));
#define mpn_invert __MPN(invert)
__GMP_DECLSPEC void mpn_invert __GMP_PROTO ((mp_ptr, mp_srcptr, mp_size_t, mp_ptr));
diff -r 6eb2023eb8f4 -r c795ab883a97 mpn/generic/mu_div_q.c
--- a/mpn/generic/mu_div_q.c Tue Dec 22 14:15:35 2009 +0100
+++ b/mpn/generic/mu_div_q.c Tue Dec 22 22:31:42 2009 +0100
@@ -154,3 +154,13 @@
TMP_FREE;
return 0;
}
+
+mp_size_t
+mpn_mu_div_q_itch (mp_size_t nn, mp_size_t dn, int mua_k)
+{
+ mp_size_t itch1, itch2;
+
+ itch1 = mpn_mu_divappr_q_itch (nn, dn, mua_k);
+ itch2 = mpn_mu_div_qr_itch (nn, dn, mua_k);
+ return MAX (itch1, itch2);
+}
diff -r 6eb2023eb8f4 -r c795ab883a97 mpn/generic/mu_divappr_q.c
--- a/mpn/generic/mu_divappr_q.c Tue Dec 22 14:15:35 2009 +0100
+++ b/mpn/generic/mu_divappr_q.c Tue Dec 22 22:31:42 2009 +0100
@@ -11,7 +11,7 @@
SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES. IN FACT, IT IS ALMOST
GUARANTEED THAT THEY WILL CHANGE OR DISAPPEAR IN A FUTURE GMP RELEASE.
-Copyright 2005, 2006, 2007 Free Software Foundation, Inc.
+Copyright 2005, 2006, 2007, 2009 Free Software Foundation, Inc.
This file is part of the GNU MP Library.
@@ -256,7 +256,7 @@
mp_ptr scratch)
{
mp_ptr rp;
- mp_size_t qn;
+ mp_size_t qn, qnown;
mp_limb_t cy;
mp_ptr tp;
mp_limb_t r;
@@ -274,6 +274,7 @@
MPN_COPY (rp, np, dn);
+ qnown = 0; /* how many limbs we've developed thus far */
while (qn > 0)
{
if (qn < in)
@@ -291,6 +292,7 @@
ASSERT_ALWAYS (cy == 0); /* FIXME */
qn -= in;
+ qnown += in;
if (qn == 0)
break;
@@ -349,16 +351,16 @@
/* We loop 0 times with about 69% probability, 1 time with about 31%
probability, 2 times with about 0.6% probability, if inverse is
computed as recommended. */
- mpn_incr_u (qp, 1);
cy = mpn_sub_n (rp, rp, dp, dn);
r -= cy;
+ cy = mpn_add_1 (qp, qp, qnown, 1);
STAT (err++);
}
if (mpn_cmp (rp, dp, dn) >= 0)
{
/* This is executed with about 76% probability. */
- mpn_incr_u (qp, 1);
- cy = mpn_sub_n (rp, rp, dp, dn);
+ mpn_sub_n (rp, rp, dp, dn);
+ cy = mpn_add_1 (qp, qp, qnown, 1);
STAT (err++);
}
@@ -380,7 +382,7 @@
quotient. For now, just make sure the returned quotient is >= the real
quotient. */
qn = nn - dn;
- cy = mpn_add_1 (qp, qp, qn, 3);
+ cy += mpn_add_1 (qp, qp, qn, 3);
if (cy != 0)
{
MPN_ZERO (qp, qn);
diff -r 6eb2023eb8f4 -r c795ab883a97 tests/mpn/Makefile.am
--- a/tests/mpn/Makefile.am Tue Dec 22 14:15:35 2009 +0100
+++ b/tests/mpn/Makefile.am Tue Dec 22 22:31:42 2009 +0100
@@ -26,7 +26,7 @@
t-toom22 t-toom32 t-toom33 t-toom42 t-toom43 t-toom44 \
t-toom52 t-toom53 t-toom62 t-toom63 t-toom6h \
t-bdiv t-hgcd t-matrix22 t-mullo t-mulmod_bnm1 t-sqrmod_bnm1 \
- t-invert
+ t-invert t-div
EXTRA_DIST = toom-shared.h
diff -r 6eb2023eb8f4 -r c795ab883a97 tests/mpn/t-bdiv.c
--- a/tests/mpn/t-bdiv.c Tue Dec 22 14:15:35 2009 +0100
+++ b/tests/mpn/t-bdiv.c Tue Dec 22 22:31:42 2009 +0100
@@ -204,5 +204,7 @@
}
}
TMP_FREE;
+
+ tests_end ();
return 0;
}
diff -r 6eb2023eb8f4 -r c795ab883a97 tests/mpn/t-div.c
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/tests/mpn/t-div.c Tue Dec 22 22:31:42 2009 +0100
@@ -0,0 +1,397 @@
+/* Copyright 2006, 2007, 2009 Free Software Foundation, Inc.
+
+This program is free software; you can redistribute it and/or modify it under
+the terms of the GNU General Public License as published by the Free Software
+Foundation; either version 3 of the License, or (at your option) any later
+version.
+
+This program 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 General Public License for more details.
+
+You should have received a copy of the GNU General Public License along with
+this program. If not, see http://www.gnu.org/licenses/. */
+
+
+#include <stdlib.h> /* for strtol */
+#include <stdio.h> /* for printf */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+#include "tests/tests.h"
+
+
+static void
+dumpy (mp_srcptr p, mp_size_t n)
+{
+ mp_size_t i;
+ if (n > 20)
+ {
+ for (i = n - 1; i >= n - 4; i--)
+ {
+ printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
+ printf (" ");
+ }
+ printf ("... ");
+ for (i = 3; i >= 0; i--)
+ {
+ printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
+ printf (" " + (i == 0));
+ }
+ }
+ else
+ {
+ for (i = n - 1; i >= 0; i--)
+ {
+ printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
+ printf (" " + (i == 0));
+ }
+ }
+ puts ("");
+}
+
+static unsigned long test;
+
+static void
+check_one (mp_srcptr qrefp, mp_srcptr rrefp, mp_ptr qp, mp_srcptr rp,
+ mp_srcptr np, mp_size_t nn, mp_srcptr dp, mp_size_t dn,
+ char *fname, mp_limb_t q_allowed_err)
+{
+ mp_size_t qn = nn - dn + 1;
+ int qcmp, rcmp;
+ mp_ptr tp;
+ mp_size_t i;
+ TMP_DECL;
+ TMP_MARK;
+
+ qcmp = mpn_cmp (qp, qrefp, qn);
+ if (qcmp != 0 && rp == rrefp) /* FIXME: rp=rrefp indicates we check a approx-q fn */
+ {
+ mp_limb_t qerr;
+
+ tp = TMP_ALLOC_LIMBS (qn);
+ qerr = mpn_sub_n (tp, qp, qrefp, qn);
+
+ for (i = 1; i < qn; i++)
+ qerr |= tp[i];
+
+ qcmp = qerr | (tp[0] > q_allowed_err);
+ }
+
+ rcmp = mpn_cmp (rp, rrefp, dn);
+
+ if (qcmp != 0 || rcmp != 0)
+ {
+ mp_size_t lo, hi;
+ printf ("\r*******************************************************************************\n");
+ printf ("%s and refmpn_tdiv_qr disagree in test %lu\n", fname, test);
+ printf ("N= "); dumpy (np, nn);
+ printf ("D= "); dumpy (dp, dn);
+ printf ("Q= "); dumpy (qp, qn);
+ printf ("refQ="); dumpy (qrefp, qn);
+ printf ("R= "); dumpy (rp, dn);
+ printf ("refR="); dumpy (rrefp, dn);
+ printf ("nn = %ld, dn = %ld, qn = %ld", nn, dn, nn-dn+1);
+ if (qcmp != 0)
+ {
+ for (lo = 0; qp[lo] == qrefp[lo]; lo++);
+ for (hi = nn - dn; qp[hi] == qrefp[hi]; hi--);
+ printf (", quotient bad %ld--%ld", hi, lo);
+ }
+ if (rcmp != 0)
+ {
+ for (lo = 0; rp[lo] == rrefp[lo]; lo++);
+ for (hi = dn - 1; rp[hi] == rrefp[hi]; hi--);
+ printf (", remainder bad %ld--%ld", hi, lo);
+ }
+ printf ("\n*******************************************************************************\n");
+ abort ();
+ }
+
+ TMP_FREE;
+}
+
+
+/* These are *bit* sizes. */
+#define SIZE_LOG 16
+#define MAX_DN (1L << SIZE_LOG)
+#define MAX_NN (1L << (SIZE_LOG + 1))
+
+#define COUNT 100
+
+mp_limb_t
+random_word (gmp_randstate_ptr rs)
+{
+ mpz_t x;
+ mp_limb_t r;
+ TMP_DECL;
+ TMP_MARK;
+
+ MPZ_TMP_INIT (x, 2);
+ mpz_urandomb (x, rs, 32);
+ r = mpz_get_ui (x);
+ TMP_FREE;
+ return r;
+}
+
+int
+main (int argc, char **argv)
+{
+ gmp_randstate_ptr rands;
+ unsigned long maxnbits, maxdbits, nbits, dbits;
+ mpz_t n, d, tz;
+ mp_size_t maxnn, maxdn, nn, dn, clearn, i;
+ mp_ptr np, dp, qp, rp, qrefp, rrefp;
+ mp_limb_t t;
+ gmp_pi1_t dinv;
+ int count = COUNT;
+ mp_ptr scratch;
+ mp_limb_t ran;
+ mp_size_t alloc, itch;
+ mp_limb_t rran0, rran1, qran0, qran1;
+ TMP_DECL;
+
+ if (argc > 1)
+ {
+ char *end;
+ count = strtol (argv[1], &end, 0);
+ if (*end || count <= 0)
+ {
+ fprintf (stderr, "Invalid test count: %s.\n", argv[1]);
+ return 1;
More information about the gmp-commit
mailing list