Probabilistic Primality: Miller–Rabin
Every time your laptop opens a secure connection it must, in a few milliseconds, conjure a fresh
prime hundreds of digits long. There is no time to factor anything and no table to look it up in.
Instead the machine plays a game of interrogation: it asks a candidate number a question that
every prime answers correctly, and watches to see whether the candidate stumbles. Composites
almost always stumble — and the Miller–Rabin test is the interrogation sharp enough
that, unlike its predecessor, essentially no composite can bluff its way past.
This is the workhorse of computational number theory: fast, randomised, and wrong only with a
probability you pin as low as you please. To see why it beats the older
Fermat test, we
start with that test's one fatal flaw.
The Fermat test, and the numbers that fool it
\textbf{Fermat's Little Theorem} promises that if
n is prime then, for any base a not divisible
by n,
a^{\,n-1} \equiv 1 \pmod{n}.
So pick a base, compute a^{n-1}\bmod n, and judge: if the answer is not
1, then n is certainly composite
— a proof of compositeness that never exhibits a factor. If the answer is
1, then n is only "probably prime": it passed a
test all primes pass, but so might an impostor. A composite that sneaks through for base
a is a Fermat liar (and a is a
liar's base); a base that exposes it is a witness.
For most composites, witnesses are everywhere and one or two random bases catch them. But a thin,
infinite family of villains ruins the plan: the Carmichael numbers. A Carmichael
number satisfies a^{n-1}\equiv 1 \pmod n for every base
a coprime to n — it is a Fermat pseudoprime to
every base at once. The smallest is
561 = 3 \times 11 \times 17,
followed by 1105, 1729,
2465,\dots — and there are infinitely many of them (Alford,
Granville and Pomerance, 1994). Against a Carmichael number the Fermat test is not merely unreliable;
it is essentially useless, no matter how many bases you try. We need a test that pierces exactly this
disguise.
The Miller–Rabin refinement
The extra idea is small and lethal: in a prime modulus, the only square roots of
1 are \pm 1. (If
x^2\equiv 1 then p \mid (x-1)(x+1), and a prime
dividing a product divides a factor, so x\equiv \pm 1.) Carmichael numbers,
being products of several primes, have other, "non-trivial" square roots of
1 — and Miller–Rabin goes looking for them.
Strip the factors of two out of n-1: for odd
n>2 write
n - 1 = 2^{s}\, d, \qquad d \text{ odd}.
Fermat's condition a^{n-1}\equiv 1 says the sequence
a^{d},\; a^{2d},\; a^{4d},\;\dots,\; a^{2^{s-1}d},\; a^{2^{s}d}=a^{n-1}
ends at 1, each term the square of the one before. If the sequence really
does reach 1, the step just before the first
1 is a square root of 1 — and for a prime that
square root must be -1. So we call n
a strong probable prime to base a if
- either a^{d}\equiv 1 \pmod n,
- or a^{2^{r} d}\equiv -1 \pmod n for some
0 \le r < s.
Every prime passes this for every base. If neither clause holds, we have caught the sequence
arriving at 1 from something other than -1 — a
non-trivial square root of 1, impossible modulo a prime — so
n is certainly composite and a
is a witness. This is precisely the check 561 cannot survive.
Worked example: n = 221, base a = 2
Take n = 221. First peel apart
n - 1 = 220 = 2^{2}\cdot 55, so s = 2 and
d = 55. Now walk the sequence for a = 2, all
arithmetic modulo 221:
| term | exponent | value mod 221 | is it \pm1? |
|---|
| a^{d} | 2^{55} | 128 | no |
| a^{2d} | 2^{110} | 30 | no (and this is a^{n-1}) |
The first clause fails (2^{55}\equiv 128\not\equiv 1), and the only
candidate for the second clause, r = 0, needs
2^{55}\equiv -1\equiv 220, which it isn't. Neither clause holds, so
2 is a witness: 221 is
composite. Indeed 221 = 13\times 17 — but notice the test never went
looking for that factorisation.
Contrast — a Fermat liar. Notice also that
2^{110}\equiv 30\not\equiv 1, so base 2 is a
Fermat witness for 221 too. Now try base
a = 21: one finds 21^{220}\equiv 1 \pmod{221},
so 21 is a Fermat liar — the Fermat test, handed this
base, would shrug and say "probably prime". Miller–Rabin, checking
21^{55}\equiv 200 and 21^{110}\equiv 220\equiv -1,
actually finds -1 and so 21 is a
strong liar here — reminding us that a single base is never enough; the guarantee comes from
the proportion of witnesses, which we turn to next.
Why one random base already wins
The engine of the whole method is a counting theorem due to Monier and Rabin.
- If n is an odd composite, then at least
\tfrac34 of the bases a\in\{1,\dots,n-1\}
are witnesses to its compositeness.
- Equivalently, at most \tfrac14 of the bases are strong liars.
Crucially there is no Carmichael-style exception: the
\tfrac34 bound holds for every odd composite, Carmichael numbers
included. So a single random base catches a composite with probability at least
\tfrac34, and k independent random bases each
miss with probability at most \tfrac14, giving a total error at most
\left(\tfrac14\right)^{k} = 4^{-k}.
With k = 20 rounds that is 4^{-20}\approx 10^{-12};
with k = 40 it is below 10^{-24} — smaller than
the chance of a hardware fault corrupting the answer mid-computation. And this
\tfrac34 is a worst case: for typical composites the real witness density
is far higher, so in practice a handful of rounds is overwhelming.
The test loop, in code
The whole algorithm is a dozen lines. Here is the core, using bigint modular
exponentiation (assume a helper modpow(base, exp, mod)):
function isProbablePrime(n: bigint, bases: bigint[]): boolean {
if (n < 2n) return false;
if (n % 2n === 0n) return n === 2n;
// Write n - 1 = 2^s * d with d odd.
let d = n - 1n, s = 0n;
while (d % 2n === 0n) { d /= 2n; s += 1n; }
witness: for (const a of bases) {
let x = modpow(a, d, n); // a^d mod n
if (x === 1n || x === n - 1n) continue; // clause 1, or -1 at r = 0
for (let r = 1n; r < s; r++) {
x = (x * x) % n; // square: a^(2^r d)
if (x === n - 1n) continue witness; // found -1 -> passes for this base
}
return false; // no clause met: a is a witness, n composite
}
return true; // survived every base: probably prime
}
Feed it random bases for the probabilistic test, or a fixed list for the deterministic variant
below. Each base costs one modular exponentiation —
O(\log n) multiplications — so the whole test is cheap even for
thousand-bit n.
Three snares surround this test.
First: "probable prime" is not "prime". A number Miller–Rabin passes could, in
principle, be composite — this is a randomised test, and its output "probably prime" carries a
residual error of at most 4^{-k}. That error is a dial you control (more
rounds, smaller error), but it is never proven to be exactly zero by the randomised test alone. Only
a witness gives certainty, and a witness only ever certifies compositeness.
Second: Carmichael numbers fool Fermat but NOT Miller–Rabin. The whole point of the
square-root check is that a Carmichael number such as 561, invulnerable to
the Fermat test for every coprime base, still has at least
\tfrac34 of its bases acting as strong witnesses. Do not carry over the
Fermat-test fear of Carmichael numbers — Rabin's bound has no exceptional composites.
Third: the deterministic version needs the Generalised Riemann Hypothesis. Miller's
result — that checking all bases up to 2(\ln n)^2 makes the test exact —
is conditional on the
Generalised
Riemann Hypothesis. If GRH is false, a composite could hide beyond that bound. Only the
fixed small witness sets (proven by direct computation) are unconditional.
Making it deterministic
Randomness is a convenience, not a necessity. Two routes remove it.
Conditional (Miller). Assuming the Generalised Riemann Hypothesis, the smallest
witness of any odd composite n is below
2(\ln n)^2. So testing every base
a \le 2(\ln n)^2 is a deterministic, polynomial-time primality
test — provably correct if GRH holds.
Unconditional (fixed witness sets). For inputs below a fixed bound, a short list of
specific bases has been verified by exhaustive computation to have no common strong liar, so
it decides primality with no hypotheses at all. For example:
| works for all n below | these bases suffice |
|---|
| 2{,}047 | \{2\} |
| 3{,}215{,}031{,}751 | \{2,3,5,7\} |
| 3.3\times 10^{24} | \{2,3,5,7,11,13,17,19,23,29,31,37\} |
The last row — the first twelve primes as bases — settles primality with certainty for every number
below 3.3\times 10^{24}, which comfortably covers 64-bit
and beyond. And for a fully unconditional polynomial-time test with no size limit and no error, the
AKS
algorithm (2002) settled the theoretical question — though Miller–Rabin's speed keeps it
the practical choice.
Yes: the Solovay–Strassen test, historically the first practical randomised
primality test (1977), and the one that first made RSA feasible. It checks Euler's
criterion instead of the square-root ladder: for a prime n and
base a coprime to it,
a^{(n-1)/2} \equiv \left(\tfrac{a}{n}\right) \pmod{n},
where \left(\tfrac{a}{n}\right) is the Jacobi symbol
(a fast-to-compute extension of the Legendre symbol via quadratic reciprocity). A composite passing
this is an Euler pseudoprime. The catch: Solovay–Strassen guarantees only that at
least \tfrac12 of bases are witnesses (error
2^{-k}), and every strong probable prime is also an Euler probable prime —
so Miller–Rabin is strictly stronger and strictly faster in convergence. That is
why Miller–Rabin, not Solovay–Strassen, is what modern libraries actually ship.
Where it lives
Miller–Rabin is the beating heart of key generation: pick a random odd number of the right size, run
the test, and on failure try the next candidate. Because primes of b bits
occur about once every 0.7\,b tries and each test is a few modular
exponentiations, a fresh 2048-bit prime appears in milliseconds. It is a small, self-contained jewel
of computational number theory — a place where a fact about square roots modulo a prime, a counting
bound, and (optionally) the deepest unproven conjecture in analytic number theory all meet inside a
loop your phone runs thousands of times a day.