[Gmp-commit] /var/hg/gmp: 2 new changesets
mercurial at gmplib.org
mercurial at gmplib.org
Thu Feb 9 10:47:37 CET 2012
details: /var/hg/gmp/rev/584cba734dfc
changeset: 14609:584cba734dfc
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Thu Feb 09 09:21:03 2012 +0100
description:
Update Copyright year.
details: /var/hg/gmp/rev/159c8b2d8ecd
changeset: 14610:159c8b2d8ecd
user: Marco Bodrato <bodrato at mail.dm.unipi.it>
date: Thu Feb 09 10:33:55 2012 +0100
description:
mpz_oddfac_1: get ready for n!!
diffstat:
ChangeLog | 6 ++++++
doc/gmp.texi | 2 +-
gmp-impl.h | 2 +-
mpz/fac_ui.c | 2 +-
mpz/oddfac_1.c | 38 +++++++++++++++++++++++++++-----------
5 files changed, 36 insertions(+), 14 deletions(-)
diffs (129 lines):
diff -r 37a6d5d3552f -r 159c8b2d8ecd ChangeLog
--- a/ChangeLog Wed Feb 08 09:41:23 2012 +0100
+++ b/ChangeLog Thu Feb 09 10:33:55 2012 +0100
@@ -1,3 +1,9 @@
+2012-02-09 Marco Bodrato <bodrato at mail.dm.unipi.it>
+
+ * mpz/oddfac_1.c (mpz_oddfac_1): Get ready for n!!
+ * gmp-impl.h (mpz_oddfac_1): Update signature.
+ * mpz/fac_ui.c (mpz_fac_ui): Update call to mpz_oddfac_1.
+
2012-02-08 Niels Möller <nisse at lysator.liu.se>
* doc/gmp.texi (mpz_gcdext): Clarified corner cases in cofactor
diff -r 37a6d5d3552f -r 159c8b2d8ecd doc/gmp.texi
--- a/doc/gmp.texi Wed Feb 08 09:41:23 2012 +0100
+++ b/doc/gmp.texi Thu Feb 09 10:33:55 2012 +0100
@@ -15,7 +15,7 @@
arithmetic library, version @value{VERSION}.
Copyright 1991, 1993, 1994, 1995, 1996, 1997, 1998, 1999, 2000, 2001, 2002,
-2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 Free Software
+2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012 Free Software
Foundation, Inc.
Permission is granted to copy, distribute and/or modify this document under
diff -r 37a6d5d3552f -r 159c8b2d8ecd gmp-impl.h
--- a/gmp-impl.h Wed Feb 08 09:41:23 2012 +0100
+++ b/gmp-impl.h Thu Feb 09 10:33:55 2012 +0100
@@ -1573,7 +1573,7 @@
__GMP_DECLSPEC mp_size_t mpz_prodlimbs __GMP_PROTO ((mpz_ptr, mp_ptr, mp_size_t));
#define mpz_oddfac_1 __gmpz_oddfac_1
-__GMP_DECLSPEC void mpz_oddfac_1 __GMP_PROTO ((mpz_ptr, mp_limb_t));
+__GMP_DECLSPEC void mpz_oddfac_1 __GMP_PROTO ((mpz_ptr, mp_limb_t, unsigned));
#define mpz_inp_str_nowhite __gmpz_inp_str_nowhite
#ifdef _GMP_H_HAVE_FILE
diff -r 37a6d5d3552f -r 159c8b2d8ecd mpz/fac_ui.c
--- a/mpz/fac_ui.c Wed Feb 08 09:41:23 2012 +0100
+++ b/mpz/fac_ui.c Thu Feb 09 10:33:55 2012 +0100
@@ -92,7 +92,7 @@
else
{
mp_limb_t count;
- mpz_oddfac_1 (x, n);
+ mpz_oddfac_1 (x, n, 0);
popc_limb (count, n);
mpz_mul_2exp (x, x, n - count);
}
diff -r 37a6d5d3552f -r 159c8b2d8ecd mpz/oddfac_1.c
--- a/mpz/oddfac_1.c Wed Feb 08 09:41:23 2012 +0100
+++ b/mpz/oddfac_1.c Thu Feb 09 10:33:55 2012 +0100
@@ -463,12 +463,18 @@
/* mpz_oddfac_1 computes the odd part of the factorial of the
parameter n. I.e. n! = x 2^a, where x is the returned value: an
odd positive integer.
+
+ If flag != 0 a square is skipped in the DSC part, e.g.
+ if n is odd, n > FAC_DSC_THRESHOLD and flag = 1, x is set to n!!.
+
+ If n is too small, flag is ignored, and an ASSERT can be triggered.
*/
void
-mpz_oddfac_1 (mpz_ptr x, mp_limb_t n)
+mpz_oddfac_1 (mpz_ptr x, mp_limb_t n, unsigned flag)
{
static const mp_limb_t tablef[] = { ONE_LIMB_ODD_FACTORIAL_TABLE };
- static const mp_limb_t tabled[] = { ONE_LIMB_ODD_DOUBLEFACTORIAL_TABLE };
+
+ ASSERT (flag == 0 || (n >= numberof (tablef) && ABOVE_THRESHOLD (n, FAC_DSC_THRESHOLD)));
if (n < numberof (tablef))
{
@@ -484,6 +490,7 @@
s = 0;
{
+ static const mp_limb_t tabled[] = { ONE_LIMB_ODD_DOUBLEFACTORIAL_TABLE };
mp_limb_t tn;
mp_size_t j;
TMP_SDECL;
@@ -546,13 +553,15 @@
TMP_MARK;
- size = n / GMP_NUMB_BITS;
- ASSERT (primesieve_size (n - 1) <= (size + 4) - (size / 2 + 3));
+ flag--;
+ size = n / GMP_NUMB_BITS + 4;
+ ASSERT (primesieve_size (n - 1) <= size - (size / 2 + 1));
/* 2-multiswing(n) < 2^(n-1)*sqrt(n/pi) < 2^(n+GMP_NUMB_BITS);
one more can be overwritten by mul, another for the sieve */
- MPZ_TMP_INIT (mswing, size + 4);
+ MPZ_TMP_INIT (mswing, size);
/* Put the sieve on the second half, it will be overwritten by the last mswing. */
- sieve = PTR(mswing) + size / 2 + 3;
+ sieve = PTR (mswing) + size / 2 + 1;
+ ASSERT ((SIZ (mswing) = 0) || ALLOC (mswing) == size);
size = (bitwise_primesieve (sieve, n - 1) + 1) / log_n_max (n) + 1;
@@ -564,15 +573,22 @@
TMP_DECL;
s--;
+ ASSERT (ABSIZ (mswing) < ALLOC (mswing) / 2); /* Check: sieve has not been overwritten */
mpz_2multiswing_1 (mswing, n >> s, sieve, factors);
TMP_MARK;
nx = SIZ (x);
- size = nx << 1;
- square = TMP_ALLOC_LIMBS (size);
- mpn_sqr (square, PTR (x), nx);
- ns = SIZ(mswing);
- size -= (square[size - 1] == 0);
+ if (s == flag) {
+ size = nx;
+ square = TMP_ALLOC_LIMBS (size);
+ MPN_COPY (square, PTR (x), nx);
+ } else {
+ size = nx << 1;
+ square = TMP_ALLOC_LIMBS (size);
+ mpn_sqr (square, PTR (x), nx);
+ size -= (square[size - 1] == 0);
+ }
+ ns = SIZ (mswing);
nx = size + ns;
MPZ_REALLOC (x, nx);
ASSERT (ns <= size);
More information about the gmp-commit
mailing list