# Utility time: find the next prime larger than

**Update 2021-06-26**: extended the utility to work for at least up to 2**81.

Today it's utility time: I'm publishing here a little utility that sits on my disk for many many years now, but hasn't been published for the lack of a more suitable medium.

The sole purpose of it is to find the next primes equal to or larger than a given number. Not primes as in RSA primes (i.e. with a lot of digits), but primes as in hash table sizes or other applications, such as scheduling recurring tasks with a minimum of overlap.

It's not something that I need often, but *occasionally* I just need one or more prime larger than, say, 10000 or so, and searching in the Debian GNU/Linux packages didn't find me an obvious candidate. So I wrote one myself: primenext (well, for that reason, and also because I wanted to learn how `Math::BigInt`

works. I don't like it very much).

## Usage

The usage is very simple:

$ primenext 65536 | head 65537 65539 65543 65551 65557 65563 65579 65581 65587 65599

It works with numbers up to `2**64`

, and although the main goal was for the program to be as simple and short as possible, it is still "fast enough" for finding magic program constants:

$ time primenext 10_000000_000000_000000 | head 10000000000000000051 10000000000000000087 10000000000000000091 10000000000000000097 10000000000000000099 10000000000000000147 10000000000000000169 10000000000000000273 10000000000000000297 10000000000000000307 real 0m0.054s user 0m0.052s sys 0m0.000s

## Theory, errrrr, Practise!

I won't go into the maths behind this program - every child knows the Miller-Rabin primality test which often is used to generate really large random primes for cryptography applications. It's probabilistic, but if you can make the chances of a bad result almost infinitely less likely than getting hit by lightning in the next minute, it's usually fine. Still, I wanted something better, so the program is using a deterministic Miller-Rabin primality test.

For this, you need a list of strong pseudoprimes and the limit to which they can reliably detect primes. Only one set is needed, but for speed reasons, a number of these sets is implemented, with the script choosing the appropriate one:

use Math::BigInt try => "GMP"; my @sprp = ( [ "341531", "9345883071009581737"], [ "520924141", "15, 750068417525532"], [ "109134866497", "2, 45650740, 3722628058"], [ "47636622961201", "2, 2570940, 211991001, 3749873356"], [ "3770579582154547", "2, 2570940, 880937, 610386380, 4130785767"], ["18446744073709551615", "2, 325, 9375, 28178, 450775, 9780504, 1795265022"], # ^^^^^^^^^^^^^^^^^^^^^^ ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ # limit .............. pseudo primes ); $_ = [map +(new Math::BigInt $_), $_->[0], split /,\s*/, $_->[1]] for @sprp;

Wikipedia has a list, but this page has a better one. As you can see, there is still room for improvement, although I would suggest not using `Math::BigInt`

for record searches.

The main program simply starts at $ARGV[0] and iterates, switching to ever larger SPRP sets as required:

my $n = (new Math::BigInt shift) | 1; my @bases; my $max = 0; while () { while ($n > $max) { @sprp or die "range exceeded\n"; ($max, @bases) = @{ shift @sprp }; } print "$n\n" if isprime $n; $n += 2; }

The actual maths is carried out by `isprime`

:

# requires @bases to be set properly before call sub isprime($) { my ($n) = @_; $n % $_ || $n == $_ || return for 3, 5, 7; # 2 not tested return 1 if $n < 121; my ($d, $s) = $n - 1; $s++ until $d & (1<<$s); $d >>= $s; a: for my $a (@bases) { my $x = $a->copy->bmodpow ($d, $n); next if $x == 1 || $x == $n - 1; for (1 .. $s - 1) { $x = $x * $x % $n; return if $x == 1; # warn "$n == 0 mod ", Math::BigInt::bgcd ($a ** ($d * 2 ** $_ / 2) - 1, $n); next a if $x == $n - 1; } return; } 1 }

Anyways, not much else to be said, you can download the `primenext`

source code and be happy with it or not.

Ah, one last thing, the illustrious Steve Worley who is responsible for many of the record entries also has a nifty way to do *really* fast primality tests, for those times where you need it, using a smallish lookup table and a single test.

## Hints welcome

Are there better SPRP lists? Probably not, but if you find one (or an altogether better method), tell me. Likewise, if you find a bug or some material improvement in this script, please tell me as well!

## Responses

As *b_jonas* pointed out on IRC, the `bsdgames`

package contains a `primes` program that is a slight superset of the `primenext` program. Also, it contains `rot13`, which is "why it's installed everywhere".

Also, *Dana Jacobsen* sent a mail with so much additional interesting info that I'll just quote it here in full (any editorial mistakes are mine):

*ntheory (aka Math::Prime::Util) comes with a primes.pl program. It isn't quite the same, but you can do:*

primes.pl 65536 +100

*to get primes between 65536 and 65536+100. It also takes expressions for things like primes.pl 2**9 "nth_prime(105)"*

*You could also use ntheory for the next_prime loop (or forprimes or an iterator), but then you don't have as much fun. You do get support for bigints, results use ES BPSW, and you could check each with is_provable_prime if you really wanted, which shouldn't slow things down until 200 or so digits.*

*For 64-bit, Wojciech's page is the best reference. There are some results for larger than 64-bit now, which is on Wikipedia.*

*Regarding hashing, there are some better results than Worley's. For 32-bit ntheory uses a single-base hashed solution similar to Forisek and Jancina. For larger that 32-bit,*

http://probableprime.org/download/example-primality.c

*shows a 2-base hash solution for 2^49 and 3-base hash for 64-bit. I wrote that almost two years ago, and maybe I should publish it. I've mentioned it on a few forums, but honestly I think BPSW is a better solution for 64-bit.*