Fast Deterministic Prime Test for n Less Than Quintillion
Learn about the Miller-Rabin test and explore a Python implementation.
Join the DZone community and get the full member experience.Join For Free
The most common way to test whether a large number is prime is the Miller-Rabin test. If the test says a number is composite, it's definitely composite. Otherwise, the number is very likely, but not certain, to be prime. A pseudoprime is a composite number that slips past the Miller-Rabin test. (Actually, a strong pseudoprime. More on that below.)
The Miller-Rabin test is actually a sequence of tests, one for each prime number. First, you run the test associated with 2, then the test associated with 3, then the one associated with 5, etc. If we knew the smallest numbers for which these tests fail, then for smaller numbers, we know for certain that they're prime if they pass. In other words, we can turn the Miller-Rabin test for probable primes into a test for provable primes.
Lower Bound on Failure
A recent result by Yupeng Jiang and Yingpu Deng finds the smallest number for which the Miller-Rabin test fails for the first nine primes. This number is:
N = 3,825,123,056,546,413,051
or more than 3.8 quintillion. So, if a number passes the first nine Miller-Rabin tests, and it's less than N, then it's prime. Not just a probable prime, but definitely prime. For a number n < N, this will be more efficient than running previously-known deterministic primality tests on n.
Let's play with this in Python. The SymPy library implements the Miller-Rabin test in a function
The following shows that N is composite and that it is a false positive for the first nine Miller-Rabin tests.
from sympy.ntheory.primetest import mr N = 3825123056546413051 assert(N == 149491*747451*34233211) ps = [2, 3, 5, 7, 11, 13, 17, 19, 23] print( mr(N, ps) )
This doesn't prove that N is the smallest number with these properties; we need the proof of Jiang and Deng for that. But assuming their result is right, here's an efficient deterministic primality test that works for all n less than N.
def is_prime(n): N = 3825123056546413051 assert(n < N) ps = [2, 3, 5, 7, 11, 13, 17, 19, 23] return mr(n, ps)
Jiang and Deng assert that N is also the smallest composite number to slip by the first 10 and 11 Miller-Rabin tests. We can show that N is indeed a strong pseudoprime for the 10th and 11th primes, but not for the 12th prime.
print( mr(N, [29, 31]) ) print( mr(N, ) )
This code prints
True for the first test and
False for the second. That is, N is a strong pseudoprime for bases 29 and 31, but not for 37.
Pseudoprimes and Strong Pseudoprimes
Fermat's little theorem says that if n is prime, then:
an-1 = 1 mod n
for all 0 < a < n. This gives a necessary but not sufficient test for primality. A (Fermat) pseudoprime for base a is a composite number n such that the above holds, an example of where the test is not sufficient.
The Miller-Rabin test refines Fermat's test by looking at additional necessary conditions for a number being prime. Often a composite number will fail one of these conditions, but not always. The composite numbers that slip by are called strong pseudoprimes or sometimes Miller-Rabin pseudoprimes.
Miller and Rabin's extra testing starts by factoring n-1 into 2 d where d is odd. If n is prime, then for all 0 < a < n either:
ad = 1 mod n
a2kd = -1 mod n
for all k satisfying 0 ≤ k < s. If one of these two conditions holds for a particular a, then n passes the Miller-Rabin test for the base a.
It wouldn't be hard to write your own implementation of the Miller-Rabin test. You'd need a way to work with large integers and to compute modular exponents, both of which are included in Python without having to use SymPy.
Let us know your thoughts and questions in the comments.
Published at DZone with permission of John Cook, DZone MVB. See the original article here.
Opinions expressed by DZone contributors are their own.