[Gmp-commit] /home/hgfiles/gmp: Added testcases for bdiv, and fixed a few bugs.
mercurial at gmplib.org
mercurial at gmplib.org
Tue Dec 8 15:55:40 CET 2009
details: /home/hgfiles/gmp/rev/0fe3db5757a6
changeset: 13006:0fe3db5757a6
user: Niels M?ller <nisse at lysator.liu.se>
date: Tue Dec 08 15:48:45 2009 +0100
description:
Added testcases for bdiv, and fixed a few bugs.
diffstat:
ChangeLog | 12 ++++
mpn/generic/bdiv_q.c | 12 +--
mpn/generic/bdiv_qr.c | 12 ++--
tests/mpn/Makefile.am | 2 +-
tests/mpn/t-bdiv.c | 167 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
5 files changed, 191 insertions(+), 14 deletions(-)
diffs (273 lines):
diff -r 415fa264bb46 -r 0fe3db5757a6 ChangeLog
--- a/ChangeLog Tue Dec 08 10:54:47 2009 +0100
+++ b/ChangeLog Tue Dec 08 15:48:45 2009 +0100
@@ -1,3 +1,15 @@
+2009-12-08 Niels Möller <nisse at lysator.liu.se>
+
+ * tests/mpn/Makefile.am (check_PROGRAMS): Added t-bdiv.
+
+ * tests/mpn/t-bdiv.c: New file.
+
+ * mpn/generic/bdiv_q.c (mpn_bdiv_q): Fixed bad quotient length,
+ should have qn == nn.
+
+ * mpn/generic/bdiv_qr.c (mpn_bdiv_qr): Pass correct nn length to
+ the lower-level functions.
+
2009-12-08 Torbjorn Granlund <tege at gmplib.org>
* tune/speed.c (routine): Add mpn_invert.
diff -r 415fa264bb46 -r 0fe3db5757a6 mpn/generic/bdiv_q.c
--- a/mpn/generic/bdiv_q.c Tue Dec 08 10:54:47 2009 +0100
+++ b/mpn/generic/bdiv_q.c Tue Dec 08 15:48:45 2009 +0100
@@ -37,24 +37,22 @@
mp_ptr tp)
{
mp_limb_t di;
- mp_size_t qn = nn - dn + 1;
- dn = MIN (qn, dn);
if (BELOW_THRESHOLD (dn, DC_BDIV_Q_THRESHOLD))
{
- MPN_COPY (tp, np, qn);
+ MPN_COPY (tp, np, nn);
binvert_limb (di, dp[0]); di = -di;
- mpn_sbpi1_bdiv_q (qp, tp, qn, dp, dn, di);
+ mpn_sbpi1_bdiv_q (qp, tp, nn, dp, dn, di);
}
else if (BELOW_THRESHOLD (dn, MU_BDIV_Q_THRESHOLD))
{
- MPN_COPY (tp, np, qn);
+ MPN_COPY (tp, np, nn);
binvert_limb (di, dp[0]); di = -di;
- mpn_dcpi1_bdiv_q (qp, tp, qn, dp, dn, di);
+ mpn_dcpi1_bdiv_q (qp, tp, nn, dp, dn, di);
}
else
{
- mpn_mu_bdiv_q (qp, np, qn, dp, dn, tp);
+ mpn_mu_bdiv_q (qp, np, nn, dp, dn, tp);
}
return;
}
diff -r 415fa264bb46 -r 0fe3db5757a6 mpn/generic/bdiv_qr.c
--- a/mpn/generic/bdiv_qr.c Tue Dec 08 10:54:47 2009 +0100
+++ b/mpn/generic/bdiv_qr.c Tue Dec 08 15:48:45 2009 +0100
@@ -46,16 +46,16 @@
MPN_COPY (tp, np, nn);
tp[nn] = 0;
binvert_limb (di, dp[0]); di = -di;
- rh = mpn_sbpi1_bdiv_qr (qp, tp, nn + 1, dp, dn, di);
- MPN_COPY (rp, tp + nn + 1 - dn, dn);
+ rh = mpn_sbpi1_bdiv_qr (qp, tp, nn, dp, dn, di);
+ MPN_COPY (rp, tp + nn - dn, dn);
}
else if (BELOW_THRESHOLD (dn, MU_BDIV_QR_THRESHOLD))
{
MPN_COPY (tp, np, nn);
tp[nn] = 0;
binvert_limb (di, dp[0]); di = -di;
- rh = mpn_dcpi1_bdiv_qr (qp, tp, nn + 1, dp, dn, di);
- MPN_COPY (rp, tp + nn + 1 - dn, dn);
+ rh = mpn_dcpi1_bdiv_qr (qp, tp, nn, dp, dn, di);
+ MPN_COPY (rp, tp + nn - dn, dn);
}
else
{
@@ -63,11 +63,11 @@
mp_ptr scratch_out;
TMP_DECL;
TMP_MARK;
- itch = mpn_mu_bdiv_qr_itch (nn + 1, dn);
+ itch = mpn_mu_bdiv_qr_itch (nn, dn);
scratch_out = TMP_BALLOC_LIMBS (itch);
MPN_COPY (tp, np, nn);
tp[nn] = 0;
- rh = mpn_mu_bdiv_qr (qp, rp, tp, nn + 1, dp, dn, scratch_out);
+ rh = mpn_mu_bdiv_qr (qp, rp, tp, nn, dp, dn, scratch_out);
TMP_FREE;
}
diff -r 415fa264bb46 -r 0fe3db5757a6 tests/mpn/Makefile.am
--- a/tests/mpn/Makefile.am Tue Dec 08 10:54:47 2009 +0100
+++ b/tests/mpn/Makefile.am Tue Dec 08 15:48:45 2009 +0100
@@ -25,7 +25,7 @@
t-instrument t-iord_u t-mp_bases t-perfsqr t-scan \
t-toom22 t-toom32 t-toom33 t-toom42 t-toom43 t-toom44 \
t-toom52 t-toom53 t-toom62 \
- t-hgcd t-matrix22 t-mullo t-mulmod_bnm1
+ t-bdiv t-hgcd t-matrix22 t-mullo t-mulmod_bnm1
EXTRA_DIST = toom-shared.h
diff -r 415fa264bb46 -r 0fe3db5757a6 tests/mpn/t-bdiv.c
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/tests/mpn/t-bdiv.c Tue Dec 08 15:48:45 2009 +0100
@@ -0,0 +1,167 @@
+/* Test bdiv functions
+
+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 <stdio.h>
+#include <stdlib.h>
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "tests.h"
+
+static int
+bdiv_q_valid_p (mp_srcptr np, mp_size_t nn, mp_srcptr dp, mp_size_t dn,
+ mp_srcptr qp)
+{
+ mp_ptr tp;
+ int c;
+ TMP_DECL;
+
+ TMP_MARK;
+ tp = TMP_ALLOC_LIMBS (dn + nn);
+ refmpn_mul_basecase (tp, qp, nn, dp, dn);
+ MPN_CMP (c, tp, np, nn);
+ TMP_FREE;
+
+ return c == 0;
+}
+
+static int
+bdiv_qr_valid_p (mp_srcptr np, mp_size_t nn, mp_srcptr dp, mp_size_t dn,
+ mp_srcptr qp, mp_srcptr rp, mp_limb_t rh)
+{
+ mp_size_t qn;
+ mp_ptr tp;
+ mp_size_t cy;
+ int c;
+ TMP_DECL;
+
+ TMP_MARK;
+ qn = nn - dn;
+ tp = TMP_ALLOC_LIMBS (nn);
+
+ if (qn >= dn)
+ refmpn_mul_basecase (tp, qp, qn, dp, dn);
+ else
+ refmpn_mul_basecase (tp, dp, dn, qp, qn);
+
+ MPN_CMP (c, tp, np, qn);
+ if (c)
+ {
+ TMP_FREE;
+ return 0;
+ }
+
+ cy = refmpn_sub_nc (tp + qn, np + qn, tp + qn, dn, 0);
+ MPN_CMP (c, tp + qn, rp, dn);
+ TMP_FREE;
+
+ return cy == rh && c == 0;
+}
+
+#define SIZE_LOG 10
+#define MAX_DN (1L << SIZE_LOG)
+#define MAX_NN (1L << (SIZE_LOG + 1))
+
+#define COUNT 1000
+int
+main (int argc, char **argv)
+{
+ mp_ptr np, dp, qp, rp, scratch;
+ mp_size_t itch;
+ int count = COUNT;
+ int test;
+ gmp_randstate_ptr rands;
+ TMP_DECL;
+ TMP_MARK;
+
+ 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;
+ }
+ }
+
+ tests_start ();
+ rands = RANDS;
+
+ dp = TMP_ALLOC_LIMBS (MAX_DN + 1);
+ np = TMP_ALLOC_LIMBS (MAX_NN + 1);
+ qp = TMP_ALLOC_LIMBS (MAX_NN + 1);
+ rp = TMP_ALLOC_LIMBS (MAX_DN + 1);
+
+ /* Claim larger nn that we'll use, to allow for a large quotient.
+ FIXME: This probably isn't the right way. */
+ itch = mpn_bdiv_qr_itch (MAX_NN + MAX_DN, MAX_DN);
+ scratch = TMP_ALLOC_LIMBS (itch);
+
+ for (test = 0; test < count; test++)
+ {
+ unsigned size_range;
+ mp_size_t nn, dn;
+ mp_limb_t rh;
+
+ /* We generate dn in the range 1 <= dn <= (1 << size_range). */
+ size_range = 1
+ + gmp_urandomm_ui (rands, SIZE_LOG);
+
+ dn = 1
+ + gmp_urandomm_ui (rands, (1L << size_range));
+ nn = dn + 1 + gmp_urandomm_ui (rands, MAX_NN - dn);
+
+ mpn_random2 (dp, dn);
+ mpn_random2 (np, nn);
+ dp[0] |= 1;
+
+ ASSERT_ALWAYS (mpn_bdiv_q_itch (nn, dn) <= itch);
+
+ mpn_bdiv_q (qp, np, nn, dp, dn, scratch);
+ if (!bdiv_q_valid_p (np, nn, dp, dn, qp))
+ {
+ gmp_printf ("bdiv_q test %d failed, nn = %d, dn = %d:\n"
+ "n = %Nx\n"
+ "d = %Nx\n"
+ "q = %Nx\n",
+ test, nn, dn,
+ np, nn, dp, dn, qp, nn);
+ abort();
+ }
+
+ ASSERT_ALWAYS (mpn_bdiv_qr_itch (nn, dn) <= itch);
+
+ rh = mpn_bdiv_qr (qp, rp, np, nn, dp, dn, scratch);
+ if (!bdiv_qr_valid_p (np,nn, dp, dn, qp, rp, rh))
+ {
+ gmp_printf ("bdiv_qr test %d failed, nn = %d, dn = %d:\n"
+ "n = %Nx\n"
+ "d = %Nx\n"
+ "q = %Nx\n"
+ "r = %Mx %Nx\n",
+ test, nn, dn,
+ np, nn, dp, dn, qp, nn - dn,
+ rh, rp, dn);
+ abort();
+ }
+ }
+ TMP_FREE;
+ return 0;
+}
More information about the gmp-commit
mailing list