From 11da511c784eca003deb90c23570f0873954e0de Mon Sep 17 00:00:00 2001 From: Duncan Wilkie Date: Sat, 18 Nov 2023 06:11:09 -0600 Subject: Initial commit. --- gmp-6.3.0/tests/devel/primes.c | 424 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 424 insertions(+) create mode 100644 gmp-6.3.0/tests/devel/primes.c (limited to 'gmp-6.3.0/tests/devel/primes.c') diff --git a/gmp-6.3.0/tests/devel/primes.c b/gmp-6.3.0/tests/devel/primes.c new file mode 100644 index 0000000..8e58962 --- /dev/null +++ b/gmp-6.3.0/tests/devel/primes.c @@ -0,0 +1,424 @@ +/* +Copyright 2018-2020 Free Software Foundation, Inc. + +This file is part of the GNU MP Library test suite. + +The GNU MP Library test suite 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. + +The GNU MP Library test suite 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 +the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */ + +/* Usage: + + ./primes [p|c] [n0] + + Checks mpz_probab_prime_p(n, r) exhaustively, starting from n=n0 + up to nMax. + If n0 * n0 > nMax, the intervall is sieved piecewise, else the + full intervall [0..nMax] is sieved at once. + With the parameter "p" (or nothing), tests all numbers. With "c" + only composites are tested. + + ./primes n|N [n0] + + Checks mpz_nextprime() exhaustively, starting from n=n0 up to + nMax. With "n", only the sequence of primes is checked, with "N" + the function is tested on every number in the interval. + + WARNING: The full intervall [0..nMax] is sieved at once, even if + only a piece is needed. This may require a lot of memory! + + */ + +#include +#include +#include "gmp-impl.h" +#include "longlong.h" +#include "tests.h" +#define STOP(x) return (x) +/* #define STOP(x) x */ +#define REPS 10 +/* #define TRACE(x,n) if ((n)>1) {x;} */ +#define TRACE(x,n) + +/* The full primesieve.c is included, just for block_resieve, that + is not exported ... */ +#undef gmp_primesieve +#include "../../primesieve.c" + +#ifndef BLOCK_SIZE +#define BLOCK_SIZE 2048 +#endif + +/*********************************************************/ +/* Section sieve: sieving functions and tools for primes */ +/*********************************************************/ + +static mp_size_t +primesieve_size (mp_limb_t n) { return n_fto_bit(n) / GMP_LIMB_BITS + 1; } + +/*************************************************************/ +/* Section macros: common macros, for swing/fac/bin (&sieve) */ +/*************************************************************/ + +#define LOOP_ON_SIEVE_CONTINUE(prime,end) \ + __max_i = (end); \ + \ + do { \ + ++__i; \ + if ((*__sieve & __mask) == 0) \ + { \ + mp_limb_t prime; \ + prime = id_to_n(__i) + +#define LOOP_ON_SIEVE_BEGIN(prime,start,end,off,sieve) \ + do { \ + mp_limb_t __mask, *__sieve, __max_i, __i; \ + \ + __i = (start)-(off); \ + __sieve = (sieve) + __i / GMP_LIMB_BITS; \ + __mask = CNST_LIMB(1) << (__i % GMP_LIMB_BITS); \ + __i += (off); \ + \ + LOOP_ON_SIEVE_CONTINUE(prime,end) + +#define LOOP_ON_SIEVE_STOP \ + } \ + __mask = __mask << 1 | __mask >> (GMP_LIMB_BITS-1); \ + __sieve += __mask & 1; \ + } while (__i <= __max_i) + +#define LOOP_ON_SIEVE_END \ + LOOP_ON_SIEVE_STOP; \ + } while (0) + +mpz_t g; + +int something_wrong (mpz_t er, int exp) +{ + fprintf (stderr, "value = %lu , expected = %i\n", mpz_get_ui (er), exp); + return -1; +} + +int +check_pprime (unsigned long begin, unsigned long end, int composites) +{ + begin = (begin / 6U) * 6U; + for (;(begin < 2) & (begin <= end); ++begin) + { + *(g->_mp_d) = begin; + TRACE(printf ("-%li ", begin),1); + if (mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 0)); + } + for (;(begin < 4) & (begin <= end); ++begin) + { + *(g->_mp_d) = begin; + TRACE(printf ("+%li ", begin),2); + if (!composites && !mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 1)); + } + if (end > 4) { + if ((end > 10000) && (begin > end / begin)) + { + mp_limb_t *sieve, *primes; + mp_size_t size_s, size_p, off; + unsigned long start; + + mpz_set_ui (g, end); + mpz_sqrt (g, g); + start = mpz_get_ui (g) + GMP_LIMB_BITS; + size_p = primesieve_size (start); + + primes = __GMP_ALLOCATE_FUNC_LIMBS (size_p); + gmp_primesieve (primes, start); + + size_s = BLOCK_SIZE * 2; + sieve = __GMP_ALLOCATE_FUNC_LIMBS (size_s); + off = n_cto_bit(begin); + + do { + TRACE (printf ("off =%li\n", off),3); + block_resieve (sieve, BLOCK_SIZE, off, primes); + TRACE (printf ("LOOP =%li - %li\n", id_to_n (off+1), id_to_n (off + BLOCK_SIZE * GMP_LIMB_BITS)),3); + LOOP_ON_SIEVE_BEGIN (prime, off, off + BLOCK_SIZE * GMP_LIMB_BITS - 1, + off, sieve); + + do { + *(g->_mp_d) = begin; + TRACE(printf ("-%li ", begin),1); + if (mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 0)); + if ((begin & 0xff) == 0) + { + spinner(); + if ((begin & 0xfffffff) == 0) + printf ("%li (0x%lx)\n", begin, begin); + } + } while (++begin < prime); + + *(g->_mp_d) = begin; + TRACE(printf ("+%li ", begin),2); + if (!composites && ! mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 1)); + ++begin; + + LOOP_ON_SIEVE_END; + off += BLOCK_SIZE * GMP_LIMB_BITS; + } while (begin < end); + + __GMP_FREE_FUNC_LIMBS (sieve, size_s); + __GMP_FREE_FUNC_LIMBS (primes, size_p); + } + else + { + mp_limb_t *sieve; + mp_size_t size; + unsigned long start; + + size = primesieve_size (end); + + sieve = __GMP_ALLOCATE_FUNC_LIMBS (size); + gmp_primesieve (sieve, end); + start = MAX (begin, 5) | 1; + LOOP_ON_SIEVE_BEGIN (prime, n_cto_bit(start), + n_fto_bit (end), 0, sieve); + + do { + *(g->_mp_d) = begin; + TRACE(printf ("-%li ", begin),1); + if (mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 0)); + if ((begin & 0xff) == 0) + { + spinner(); + if ((begin & 0xfffffff) == 0) + printf ("%li (0x%lx)\n", begin, begin); + } + } while (++begin < prime); + + *(g->_mp_d) = begin; + TRACE(printf ("+%li ", begin),2); + if (!composites && ! mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 1)); + ++begin; + + LOOP_ON_SIEVE_END; + + __GMP_FREE_FUNC_LIMBS (sieve, size); + } + } + + for (;begin < end; ++begin) + { + *(g->_mp_d) = begin; + TRACE(printf ("-%li ", begin),1); + if (mpz_probab_prime_p (g, REPS)) + STOP (something_wrong (g, 0)); + } + + gmp_printf ("%Zd\n", g); + return 0; +} + +int +check_nprime (unsigned long begin, unsigned long end) +{ + if (begin < 2) + { + *(g->_mp_d) = begin; + g->_mp_size = begin; + TRACE(printf ("%li ", begin),1); + mpz_nextprime (g, g); + if (mpz_cmp_ui (g, 2) != 0) + STOP (something_wrong (g, 2)); + begin = mpz_get_ui (g); + } + if (begin < 3) + { + *(g->_mp_d) = begin; + TRACE(printf ("%li ", begin),1); + mpz_nextprime (g, g); + if (mpz_cmp_ui (g, 3) != 0) + STOP (something_wrong (g, 3)); + begin = mpz_get_ui (g); + } + if (end > 4) + { + mp_limb_t *sieve; + mp_size_t size; + unsigned long start; + + size = primesieve_size (end); + + sieve = __GMP_ALLOCATE_FUNC_LIMBS (size); + gmp_primesieve (sieve, end); + start = MAX (begin, 5) | 1; + *(g->_mp_d) = begin; + LOOP_ON_SIEVE_BEGIN (prime, n_cto_bit(start), + n_fto_bit (end), 0, sieve); + + mpz_nextprime (g, g); + if (mpz_cmp_ui (g, prime) != 0) + STOP (something_wrong (g, prime)); + + if (prime - start > 200) + { + start = prime; + spinner(); + if (prime - begin > 0xfffffff) + { + begin = prime; + printf ("%li (0x%lx)\n", begin, begin); + } + } + + LOOP_ON_SIEVE_END; + + __GMP_FREE_FUNC_LIMBS (sieve, size); + } + + if (mpz_cmp_ui (g, end) < 0) + { + mpz_nextprime (g, g); + if (mpz_cmp_ui (g, end) <= 0) + STOP (something_wrong (g, -1)); + } + + gmp_printf ("%Zd\n", g); + return 0; +} + +int +check_Nprime (unsigned long begin, unsigned long end) +{ + mpz_t op; + mpz_init_set_ui (op, end); + + for (;begin < 2; ++begin) + { + *(op->_mp_d) = begin; + op->_mp_size = begin; + TRACE(printf ("%li ", begin),1); + mpz_nextprime (g, op); + if (mpz_cmp_ui (g, 2) != 0) + STOP (something_wrong (g, 2)); + } + if (begin < 3) + { + *(op->_mp_d) = begin; + TRACE(printf ("%li ", begin),1); + mpz_nextprime (g, op); + if (mpz_cmp_ui (g, 3) != 0) + STOP (something_wrong (g, 3)); + begin = 3; + } + if (end > 4) + { + mp_limb_t *sieve; + mp_size_t size; + unsigned long start; + unsigned long opl; + + size = primesieve_size (end); + + sieve = __GMP_ALLOCATE_FUNC_LIMBS (size); + gmp_primesieve (sieve, end); + start = MAX (begin, 5) | 1; + opl = begin; + LOOP_ON_SIEVE_BEGIN (prime, n_cto_bit(start), + n_fto_bit (end), 0, sieve); + + do { + *(op->_mp_d) = opl; + mpz_nextprime (g, op); + if (mpz_cmp_ui (g, prime) != 0) + STOP (something_wrong (g, prime)); + ++opl; + } while (opl < prime); + + if (prime - start > 200) + { + start = prime; + spinner(); + if (prime - begin > 0xfffffff) + { + begin = prime; + printf ("%li (0x%lx)\n", begin, begin); + } + } + + LOOP_ON_SIEVE_END; + + __GMP_FREE_FUNC_LIMBS (sieve, size); + } + + if (mpz_cmp_ui (g, end) < 0) + { + mpz_nextprime (g, g); + if (mpz_cmp_ui (g, end) <= 0) + STOP (something_wrong (g, -1)); + } + + gmp_printf ("%Zd\n", g); + return 0; +} + +int +main (int argc, char **argv) +{ + int ret, mode = 0; + unsigned long begin = 0, end = 0; + + for (;argc > 1;--argc,++argv) + switch (*argv[1]) { + case 'p': + mode = 0; + break; + case 'c': + mode = 2; + break; + case 'n': + mode = 1; + break; + case 'N': + mode = 3; + break; + default: + begin = end; + end = atol (argv[1]); + } + + if (begin >= end) + { + fprintf (stderr, "usage: primes [N|n|p|c] [n0] \n"); + exit (1); + } + + mpz_init_set_ui (g, ULONG_MAX); + + switch (mode) { + case 1: + ret = check_nprime (begin, end); + break; + case 3: + ret = check_Nprime (begin, end); + break; + default: + ret = check_pprime (begin, end, mode); + } + + mpz_clear (g); + + if (ret == 0) + printf ("Prime tests checked in [%lu - %lu] [0x%lx - 0x%lx].\n", begin, end, begin, end); + return ret; +} -- cgit v1.2.3