21

For very small numbers, say 32 bits unsigned, testing all divisors up to the square root is a very decent approach. Some optimizations can be made to avoid trying all divisors, but these yield marginal improvements. The complexity remains $O(\sqrt n)$.

On the other hand, much faster primality tests are available, but they are pretty sophisticated and deploy their efficiency for much longer numbers.

Is there an intermediate solution, i.e. a relatively simple algorithm, that is of practical use for, say, 64 bits unsigned, with a target running time under 1 ms ?

I am not after micro-optimization of the exhaustive division method. I am after a better working principle, of a reasonable complexity (and of the deterministic type).


Update:

Using a Python version of the Miller-Rabin test from Rosetta code, the time for the prime $2^{64}-59=18446744073709551557$ is $0.7$ ms. (Though this is not a sufficient test because nothing says we are in a worst case.)

http://rosettacode.org/wiki/Miller%E2%80%93Rabin_primality_test#Python:_Proved_correct_up_to_large_N

And I guess that this code can be improved for speed.

  • You can build a sieve filtering out multiples of any primes you have already tried as you go. But it will of course require more memory to keep track of. – mathreadler Oct 20 '17 at 08:37
  • 2
    @Zubzub: for a prime close to 2^64, the function will try close to 2^32 divisions, so no. (By the way, for 18446744073709551557 it takes six minutes). –  Oct 20 '17 at 08:38
  • 2
    @mathreadler: I doubt this is usable for 64 bit integers. And if I am right, the gain will not exceed a factor log(n), not counting overhead. –  Oct 20 '17 at 08:40
  • 100 million divisions instead of 2 billion is better, considering the division is slow. Storing a table requires less than 20 megabytes if we store primality as bits. – mathreadler Oct 20 '17 at 09:06
  • @mathreadler I have a feeling that algorithms such as Miller-Rabin will already be much faster than brute force for 64 bits, and from the question, it seems that even those are too slow. – Qudit Oct 20 '17 at 09:16
  • I think getting down to 1 s should be doable with that approach, but 1ms, you are right... maybe not. – mathreadler Oct 20 '17 at 09:32
  • There are some statements in https://math.stackexchange.com/questions/2450722/whats-the-fastest-and-most-efficent-way-to-find-prime-numbers/2450795#2450795, especially the 64-bit section. Maybe @DanaJ will answer/comment, if not, you can find native C code at his page https://github.com/danaj/Math-Prime-Util. See also https://math.stackexchange.com/a/897442/61216 – gammatester Oct 20 '17 at 09:45
  • Hey, we don't even need to do a division for each new prime, we can also do a "count up, accumulate and find minimum" loop. That would be much nicer for the CPU pipelines and give much lower latencies. – mathreadler Oct 20 '17 at 09:46
  • There exist really efficient SIMD instructions for finding $\displaystyle\min_{i}|a-b_i|$ for like 8 or 16 such expressions a cycle. – mathreadler Oct 20 '17 at 09:52
  • 1
    $1$ second is absolutely no proeblem. Even $100$-digit numbers can be proven to be prime with PARI/GP using the Adleman-Pomerance-Rumely-test within about $200$ milliseconds. But $1ms$ ? Perhaps the BPSW-test is a good idea, which is correct upto at least $2^{64}$. Trial division makes only sense, if we want to test many numbers, to reduce the number of candidates, but verifying $64$-bit numbers via trial division actually is not efficient. But checking the first few primes (lets say upto $100$) could slightly improve the test. – Peter Oct 20 '17 at 10:51
  • I tried a random $20$-digit prime and applied the Adleman-Pomernace-Rumely-test with PARI/GP. The timer displays "$0$", so it seems to take less than a millisecond. – Peter Oct 20 '17 at 10:55
  • I suggest you move your Update to an answer. – lhf Oct 20 '17 at 12:39
  • @lhf: maybe later, if I see that optimizations are possible. In any case I'll keep the accepted answer. –  Oct 20 '17 at 14:47
  • 1
    A deterministic implementation with a maximum of 7 base tests, here. – Brett Hale Sep 05 '18 at 09:33

3 Answers3

20

I think some version of the (deterministic) Miller-Rabin Test should do the trick. It is usually used as a probabilistic test for whopping big numbers, but it can be repurposed as a deterministic test for smaller fixed ranges.

The core of the Miller-Rabin test is the notion of a "probable prime" for some base $a$. Specifically, let $n$ be an odd prime, and write $n - 1 = d \times 2^s$ for some odd $d$. Then, it follows that $a^d \equiv 1 \pmod{n}$ or $a^{d 2^r} \equiv -1 \pmod{n}$ for some $0 \leq r < s$. (The wikipedia page has reasoning for this).

If $n$ is any number, and $a<n$, we could run the same test, and if that test passed we would call $n$ a strong pseudoprime for base $a$. The usual Miller-Rabin test is based on doing this for a lot of different (randomly chosen) $a$, and some number-theoretic argument saying that if $n$ is composite, the probability of not finding an $a$ demonstrating this is vanishingly small (after trying a lot of them).

However, if Wikipedia is to be believed (I haven't followed up the reference), then for $n < 2^{64}$ it is sufficient to test $a \in \{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37\}$. Somehow there are no composite numbers below $2^{64}$ which are a strong pseudoprime for all these $a$.

The above test is very fast, and altogether prime testing would require at most $12 \times 64$ modular exponentiations (this will be well below 1ms). For smaller $n$, you could even use a smaller list of possible $a$'s.

Joppy
  • 13,983
  • I guess that this is exactly what I was looking for: nice tradeoff between complexity and efficiency. –  Oct 20 '17 at 11:01
  • @YvesDaoust For so small numbers, the adleman-pomerance-rumely-test is not much slower (even with the slow program PARI/GP on my also slow computer $64$-bit digit-.numbers require less than $1$ millisecond for a decision. – Peter Oct 20 '17 at 11:12
  • @Peter: The only description of that algorithm I can find goes for pages and pages though. It doesn't look "relatively simple". The algorithm I described I can implement in about 20 lines of Python. – Joppy Oct 20 '17 at 11:24
  • @Joppy Which, knowing what they say about Python, might correspond to 20 pages of code in any other language! – JiK Oct 20 '17 at 11:52
  • 3
    @JiK: Well, I meant reasonably. In another language like C, you might have to write your own modular exponentiation, but that's about it. – Joppy Oct 20 '17 at 11:58
  • 5
    According to http://miller-rabin.appspot.com/, it is sufficient to run the Miller–Rabin test for $a\in{2,325,9375,28178,450775,9780504,1795265022}$; that’s 7 trials instead of 12. – Emil Jeřábek Oct 20 '17 at 15:53
  • In practice then, primes require 7 tests, while most composite numbers just take one test, and few take 3 or more. I'd say the average number of tests is likely less than 1.25. – gnasher729 Oct 20 '17 at 19:51
12

oeis/A014233 says:

The primality of numbers $\lt 2^{64}$ can be determined by asserting strong pseudoprimality to all prime bases $\le 37$.

The reference is the recent paper Strong pseudoprimes to twelve prime bases by Sorenson and Webster.

For code, see Prime64 and also the primes programs in FreeBSD, especially spsp.c.

lhf
  • 221,500
  • This is true, but depends on previous work (i.e. somebody has already checked this works). Does that make it any different to downloading a list of primes? – Henry Oct 20 '17 at 10:36
  • 3
    @Henry, yes, because the list is too long. The OP only wants to test, not to list. – lhf Oct 20 '17 at 10:37
  • @Henry: lhf is right. I am looking for an algorithm that avoids usage of big tables. –  Oct 20 '17 at 10:59
0

Testing for 64-bit integers can be done deterministically in very short time using modifications of the Baille-PSW test or the Miller-Rabin test.

The fastest variant of the Miller-Rabin test involves constructing a list of bases that removes all pseudoprimes within the interval, this can take a large list of bases, 262,144 in the case of a test that uses two bases. A description of such a test is provided by Forisek and Jancina. Alternately one can use a slightly slower modification of the Baille-PSW test that does not require such a large set of bases. The worst case complexity of this test is approximately 2.5 fermat tests.

There is an highly optimised implementation of both of these algorithms in the machine-prime library that is in the U.S public domain so it can be easily copied or translated. Bindings are also provided to multiple languages, including Python, as copying is somewhat tricky and not necessarily even feasible (as many of the optimisations exploit finite ring arithmetic that CPU's use,but may not be available in interpreted languages).

Machine-prime easily satisfies the required speed even in the hardest case, on nearly any machine. On a moderate-end cpu at present day (2024), it is approximately 1 million times faster in the hardest case (2^64-59) than the algorithm in the post. (About 700K times faster for the Python binding)

(Disclaimer: I am the author of machine-prime).

JASory
  • 26