diff options
author | christos <christos@NetBSD.org> | 2018-02-03 15:40:29 +0000 |
---|---|---|
committer | christos <christos@NetBSD.org> | 2018-02-03 15:40:29 +0000 |
commit | e5c17e55df9ac3ba0cc4e6836b1296ca4e09c863 (patch) | |
tree | 583cfb88d5033a8358b72877f951c45d1307673f /primes | |
parent | 18958ee2a20b22fe016fc80f46385a0c7cbe04f1 (diff) | |
download | bsdgames-darwin-e5c17e55df9ac3ba0cc4e6836b1296ca4e09c863.tar.gz bsdgames-darwin-e5c17e55df9ac3ba0cc4e6836b1296ca4e09c863.tar.zst bsdgames-darwin-e5c17e55df9ac3ba0cc4e6836b1296ca4e09c863.zip |
PR/52976: Eitan Adler: handle larger primes
Using results from
J. Sorenson and J. Webster, Strong pseudoprimes to twelve prime
bases, Math. Comp. 86(304):985-1003, 2017.
teach primes(6) to enumerate primes up to 2^64 - 1. Until Sorenson
and Webster's paper, we did not know how many strong speudoprime tests
were required when testing alleged primes between 3825123056546413051
and 2^64 - 1.
Adapted from: FreeBSD
Diffstat (limited to 'primes')
-rw-r--r-- | primes/primes.6 | 13 | ||||
-rw-r--r-- | primes/primes.c | 9 | ||||
-rw-r--r-- | primes/primes.h | 5 | ||||
-rw-r--r-- | primes/spsp.c | 42 |
4 files changed, 38 insertions, 31 deletions
diff --git a/primes/primes.6 b/primes/primes.6 index 36c76cd1..aa000ef4 100644 --- a/primes/primes.6 +++ b/primes/primes.6 @@ -1,4 +1,4 @@ -.\" $NetBSD: primes.6,v 1.5 2014/10/04 13:15:50 wiz Exp $ +.\" $NetBSD: primes.6,v 1.6 2018/02/03 15:40:29 christos Exp $ .\" .\" Copyright (c) 1989, 1993 .\" The Regents of the University of California. All rights reserved. @@ -35,7 +35,7 @@ .\" .\" By Landon Curt Noll, http://www.isthe.com/chongo/index.html /\oo/\ .\" -.Dd February 3, 2008 +.Dd February 2, 2018 .Dt PRIMES 6 .Os .Sh NAME @@ -100,14 +100,7 @@ Originally by .An Landon Curt Noll , extended to some 64-bit primes by .An Colin Percival . -.Sh CAVEATS +.Sh BUGS This .Nm program won't get you a world record. -.Pp -The program is not able to list primes between -3825123056546413050 and 18446744073709551615 (2^64 -- 1) as it relies on strong pseudoprime tests after -sieving, and it is yet unknown how many of those -tests are needed to prove primality for integers -larger than 3825123056546413050. diff --git a/primes/primes.c b/primes/primes.c index 7c874711..a022aa96 100644 --- a/primes/primes.c +++ b/primes/primes.c @@ -1,4 +1,4 @@ -/* $NetBSD: primes.c,v 1.21 2014/10/04 13:15:50 wiz Exp $ */ +/* $NetBSD: primes.c,v 1.22 2018/02/03 15:40:29 christos Exp $ */ /* * Copyright (c) 1989, 1993 @@ -42,7 +42,7 @@ __COPYRIGHT("@(#) Copyright (c) 1989, 1993\ #if 0 static char sccsid[] = "@(#)primes.c 8.5 (Berkeley) 5/10/95"; #else -__RCSID("$NetBSD: primes.c,v 1.21 2014/10/04 13:15:50 wiz Exp $"); +__RCSID("$NetBSD: primes.c,v 1.22 2018/02/03 15:40:29 christos Exp $"); #endif #endif /* not lint */ @@ -118,7 +118,7 @@ main(int argc, char *argv[]) argv += optind; start = 0; - stop = SPSPMAX; + stop = (uint64_t)(-1); /* * Convert low and high args. Strtoumax(3) sets errno to @@ -145,9 +145,6 @@ main(int argc, char *argv[]) err(1, "%s", argv[1]); if (*p != '\0') errx(1, "%s: illegal numeric format.", argv[1]); - if (stop > SPSPMAX) - errx(1, "%s: stop value too large (>%" PRIu64 ").", - argv[1], (uint64_t) SPSPMAX); break; case 1: /* Start on the command line. */ diff --git a/primes/primes.h b/primes/primes.h index 20da4ef6..2c004465 100644 --- a/primes/primes.h +++ b/primes/primes.h @@ -1,4 +1,4 @@ -/* $NetBSD: primes.h,v 1.6 2014/10/02 21:36:37 ast Exp $ */ +/* $NetBSD: primes.h,v 1.7 2018/02/03 15:40:29 christos Exp $ */ /* * Copyright (c) 1989, 1993 @@ -69,6 +69,3 @@ extern const size_t pattern_size; /* length of pattern array */ /* Test for primality using strong pseudoprime tests. */ int isprime(uint64_t); - -/* Maximum value which the SPSP code can handle. */ -#define SPSPMAX 3825123056546413050ULL diff --git a/primes/spsp.c b/primes/spsp.c index 7d1b0841..998751c3 100644 --- a/primes/spsp.c +++ b/primes/spsp.c @@ -1,4 +1,4 @@ -/* $NetBSD: spsp.c,v 1.1 2014/10/02 21:36:37 ast Exp $ */ +/* $NetBSD: spsp.c,v 1.2 2018/02/03 15:40:29 christos Exp $ */ /*- * Copyright (c) 2014 Colin Percival @@ -36,7 +36,7 @@ __COPYRIGHT("@(#) Copyright (c) 1989, 1993\ #if 0 static char sccsid[] = "@(#)primes.c 8.5 (Berkeley) 5/10/95"; #else -__RCSID("$NetBSD: spsp.c,v 1.1 2014/10/02 21:36:37 ast Exp $"); +__RCSID("$NetBSD: spsp.c,v 1.2 2018/02/03 15:40:29 christos Exp $"); #endif #endif /* not lint */ @@ -46,23 +46,33 @@ __RCSID("$NetBSD: spsp.c,v 1.1 2014/10/02 21:36:37 ast Exp $"); #include "primes.h" -/* Return a * b % n, where 0 <= a, b < 2^63, 0 < n < 2^63. */ +/* Return a * b % n, where 0 <= n. */ static uint64_t mulmod(uint64_t a, uint64_t b, uint64_t n) { uint64_t x = 0; + uint64_t an = a % n; while (b != 0) { - if (b & 1) - x = (x + a) % n; - a = (a + a) % n; + if (b & 1) { + x += an; + if ((x < an) || (x >= n)) + x -= n; + } + if (an + an < an) + an = an + an - n; + else if (an + an >= n) + an = an + an - n; + else + an = an + an; + b >>= 1; } return (x); } -/* Return a^r % n, where 0 <= a < 2^63, 0 < n < 2^63. */ +/* Return a^r % n, where 0 < n. */ static uint64_t powmod(uint64_t a, uint64_t r, uint64_t n) { @@ -186,10 +196,20 @@ isprime(uint64_t _n) return (0); if (n < 3825123056546413051) return (1); + /* + * Value from: + * J. Sorenson and J. Webster, Strong pseudoprimes to twelve prime + * bases, Math. Comp. 86(304):985-1003, 2017. + */ - /* We can't handle values larger than this. */ - assert(n <= SPSPMAX); + /* No SPSPs to bases 2..37 less than 318665857834031151167461. */ + if (!spsp(n, 29)) + return (0); + if (!spsp(n, 31)) + return (0); + if (!spsp(n, 37)) + return (0); - /* UNREACHABLE */ - return (0); + /* All 64-bit values are less than 318665857834031151167461. */ + return (1); } |