∑MathIntermediate

Miller-Rabin Primality Test

Key Points

•
Miller–Rabin is a fast primality test that uses modular exponentiation to detect compositeness with very high reliability.
•
It rewrites n − 1 as 2^s × d with d odd and checks if certain powers of a base a are 1 or −1 modulo n.
•
If none of the checks pass for a base a, then a is a witness that n is composite.
•
Running the test with multiple independent bases makes the probability of error drop exponentially, at most $\frac{1}{4}^{k}$ after k rounds.
•
For 64-bit integers, a fixed small set of bases makes the test fully deterministic in practice.
•
Correct implementation in C++ requires overflow-safe modular multiplication (e.g., using __int128).
•
Miller–Rabin is typically combined with small-prime trial division to quickly filter easy composites before the expensive steps.
•
It is a core building block for competitive programming tasks and for algorithms like Pollard’s Rho factorization.

Prerequisites

→Modular arithmetic — Understanding congruences and modulo operations is essential for the test’s conditions and reductions.
→Fast exponentiation (binary exponentiation) — Miller–Rabin relies on computing a^d mod n efficiently in O(log d) steps.
→Bit operations and integer types — Decomposing n − 1 into 2^s × d and preventing overflow require comfort with bit tricks and wide integers.
→Number theory basics (gcd, primes, composites) — Reasoning about primes, gcd, and the multiplicative group modulo n is foundational.
→Probability and error bounds — To choose the number of rounds k, you must interpret and control the error probability (1/4)^k.

Detailed Explanation

Tap terms for definitions

01Overview

The Miller–Rabin primality test is a fast algorithm to decide whether a number n is prime or composite. Unlike simple trial division, which checks many possible factors, Miller–Rabin uses properties of modular arithmetic to find a contradiction that proves compositeness. The algorithm starts by expressing n − 1 as 2^s × d, where d is odd. For a chosen base a (with 2 ≤ a ≤ n − 2), it computes powers of a modulo n. If the results follow a pattern consistent with primes (specifically, hitting 1 or −1 at the right steps), n may be prime; otherwise, a provides evidence (a “witness”) that n is composite. One round is often enough to detect many composites, but we usually repeat with different bases to reduce the chance of falsely classifying a composite as prime. The beautiful part is the error shrinks exponentially with each independent round. For 64-bit integers, mathematicians have identified small fixed sets of bases that make the test deterministic, meaning zero error if implemented correctly. In practical programming, Miller–Rabin is the go-to primality test when n can be as large as 2^64 − 1. It underpins many higher-level algorithms like Pollard’s Rho for factorization, and is widely used in cryptography to quickly identify primes before generating keys (with larger big-integer variants).

02Intuition & Analogies

Imagine checking whether a key fits a lock. A naive way is to try every possible key (trial division), which is slow. Miller–Rabin is like a collection of cleverly designed “test keys” that don’t need to open the lock; they just need to detect that a given blank is definitely wrong. Each test key (base a) interacts with the lock (the number n) in a structured way: you turn it a certain number of times (raising to powers modulo n). If at some stage the motion contradicts the pattern all real locks (primes) must follow, you know immediately the key is fake (n is composite). If the motion always looks plausible, the blank might be a real key (n might be prime), but it could also be a very good fake that mimics the pattern (a pseudoprime). The trick is: most fakes fail quickly for random test keys, and the chance they pass multiple independent tests is tiny. The 2^s × d decomposition is like pulling back a combination lock to the point where important clicks happen: starting at a^d (first click), then doubling the number of turns each time (s more clicks). For a prime, you must either be exactly at the “start” (1) on the first click or land exactly at “−1” at some click; otherwise, the lock is fake. By using several different keys (bases), you layer confidence. For 64-bit locks, there is a known finite set of keys that is guaranteed to tell you the truth every time, making the procedure deterministic. This is why competitive programmers and cryptographers love Miller–Rabin: it’s fast, reliable, and elegant.

03Formal Definition

Let

n > 2

be an odd integer. Write n − 1 as 2^s × d with d odd (i.e., repeatedly divide n − 1 by 2 until it becomes odd). For a base a with 2 ≤

a \leq n

− 2 and gcd(a, n) = 1, define

x_{0}

≔

a^{d}

mod n and

x_{r + 1}

≔

x_{r}

^2 mod n for

r = 0

, 1, …, s − 1. The Miller–Rabin test declares n passes for base a (is a “strong probable prime” to base a) if either

x_{0}

≡ 1 (mod n) or there exists r ∈ {0, …, s − 1} such that

x_{r}

≡ −1 (mod n). Otherwise, a is a Miller–Rabin witness for the compositeness of n. If n is prime, it passes for every base a with 1 <

a < n

(by properties following from Fermat’s little theorem and the structure of cyclic groups modulo primes). If n is composite, then at least

\frac{3}{4}

of the bases 2 ≤

a \leq n

− 2 are witnesses; equivalently, at most

\frac{1}{4}

are “strong liars” that mistakenly pass n. Therefore, repeating the test independently with k uniformly random bases yields an error probability bounded above by

\frac{1}{4}^{k}

. For 64-bit integers, specific fixed base sets (e.g., {2, 3, 5, 7, 11, 13, 17} or the stronger set listed in the prompt) are known to detect all composites, making the test deterministic over that domain.

04When to Use

Use Miller–Rabin whenever you need fast primality checks for integers up to 64 bits (and beyond with big-integer support). It excels when trial division is too slow, such as when n has no small prime factors or when you must test many large numbers. In competitive programming, it is ideal for problems requiring checking primality of up to, say, 10^5 numbers near 10^18, finding the next prime ≥ n, or prefiltering numbers before attempting Pollard’s Rho factorization. In cryptographic contexts (with big integers), Miller–Rabin is used with multiple random bases to generate probable primes efficiently; the small chance of error is made negligible by using enough rounds and additional checks (like strong Lucas tests). In mixed strategies, combine a quick small-prime sieve (trial division by primes up to a bound like 10^6) with Miller–Rabin for the remaining hard cases; this speeds up average performance dramatically. Choose a deterministic fixed-base set whenever your inputs are guaranteed to be 64-bit integers, and switch to randomized bases when n can grow beyond that range and you accept a tunable, vanishingly small error probability.

⚠️Common Mistakes

Overflow in modular multiplication: Doing (a * b) % n with 64-bit a, b, n can overflow before the modulo. Use 128-bit intermediates (__int128) or Montgomery/Barrett reduction to avoid undefined behavior.
Forgetting edge cases: Handle n < 2, n = 2, and even n carefully. A common bug is to treat 2 as odd or to skip early returns for small primes.
Using std::pow for modular exponentiation: std::pow works on floating-point and is incorrect for modular arithmetic. Implement fast modular exponentiation by repeated squaring under modulo.
Picking bad bases: For 64-bit determinism, use a proven fixed set of bases. For randomized tests, ensure bases are in [2, n − 2] and coprime to n; if gcd(a, n) > 1, you already found a factor.
Not reusing n − 1 = 2^s × d: Recomputing s and d for every base wastes time. Factor once per n and reuse.
Ignoring small-prime filters: Failing to trial-divide by small primes (e.g., primes ≤ 37 or a small sieve) slows the program on obvious composites.
Misinterpreting the error bound: The (1/4)^k bound assumes independent, appropriately chosen bases; using the same base or a biased generator weakens guarantees.
Returning “prime” on first pass: The test only provides certainty for deterministically proven base sets over bounded domains; otherwise, return “probable prime” unless you used such a set.

Key Formulas

2-adic decomposition

n - 1 = 2^{s} \cdot d, d odd

Explanation: Write n − 1 as a power of two times an odd number. This sets up the sequence of squarings used by the test.

Squaring sequence

x_{0} \equiv a^{d} mod n, x_{r + 1} \equiv x_{r}^{2} mod n

Explanation: Starting from $a^{d}$ mod n, square repeatedly s times. Hitting 1 immediately or −1 at some step is the expected behavior for primes.

Miller–Rabin condition

Pass if a^{d} \equiv 1 mod n or \exists r \in {0, \dots, s - 1} : a^{2^{r} d} \equiv - 1 mod n

Explanation: If neither condition holds, a is a witness to compositeness. If a condition holds, n is a strong probable prime to base a.

One-round error bound

Pr [MR errs on composite n in one round] \leq \frac{1}{4}

Explanation: For any composite n, at most one quarter of the bases are strong liars. Choosing a base uniformly at random gives at most 25% error for one round.

k-round error bound

Pr [MR errs after k independent rounds] \leq (\frac{1}{4})^{k}

Explanation: Each independent round reduces the error by at least a factor of 4. After k rounds, the error decays exponentially.

Rounds for target error

k \geq ⌈ lo g_{4} (\frac{1}{δ}) ⌉

Explanation: To achieve error at most $δ,$ run at least this many independent bases. For example, δ = 2^{-64} needs $k \geq 32$ .

Fermat’s little theorem

a^{p - 1} \equiv 1 mod p (g cd (a, p) = 1)

Explanation: This underlies the correctness of the passing conditions for primes. It ensures certain powers collapse to 1 modulo a prime.

Modular exponentiation time

T_{modexp} (e) = O (lo g e)

Explanation: Fast exponentiation by squaring uses O(log e) modular multiplications. Each multiplication may be O(1) on a fixed-width machine word model.

Strong liars bound

∣ {a \in [2, n - 2] : n passes MR base a} ∣ \leq \frac{n - 1}{4}

Explanation: At most a quarter of bases can falsely accept a composite. This gives the $\frac{1}{4}$ per-round error bound.

Euler’s totient function

φ (n) = (Z / n Z)^{\times}

Explanation: The size of the multiplicative group modulo n. Its structure explains why primes behave cleanly while composites can misbehave.

Complexity Analysis

Let n be a 64-bit integer and let k be the number of Miller–Rabin bases used. The dominant operation per base is modular exponentiation

a^{d}

mod n, computed via repeated squaring. This takes O(log d) ≈ O(log n) modular multiplications. On a fixed-width word-RAM model (like typical 64-bit CPUs), each modular multiplication and reduction can be treated as O(1) time if we allow 128-bit intermediates (__int128) or use specialized reductions (Barrett or Montgomery). Therefore, the practical running time is O(k log n) word operations. In bit-complexity terms (counting bit operations rather than machine words), a multiplication of two O(log n)-bit integers costs M(log n), where naive multiplication yields O((log n)^2). Hence modular exponentiation costs O(log n · M(log n)) = O((log n)^3) with schoolbook multiplication, and the full test is O(k (log n)^3) bit operations. Space usage is small and essentially constant: we store a handful of 64-bit integers (n, a, d, s, loop temporaries). Thus the auxiliary space is O(1) on fixed-width hardware. If you precompute small primes for trial division (e.g., sieve up to 10^6), the extra space becomes O(B) for the sieve size B, independent of the magnitude of n being tested. In practice, with a deterministic base set for 64-bit inputs (a dozen or fewer bases), the test completes in microseconds per number. Combining quick small-prime filters (constant-time checks for divisibility by the first few primes) often cuts average running time by an order of magnitude, since many composites are eliminated before modular exponentiation begins. Randomized variants allow trading time for an arbitrarily small error probability by selecting k to meet a target δ using k ≥ ⌈log₄(1/

δ)⌉ .

Code Examples

Deterministic Miller–Rabin for 64-bit integers (safe __int128)

1 #include <bits/stdc++.h>
2 using namespace std;
3 
4 using u128 = unsigned __int128;
5 using u64 = uint64_t;
6 
7 // Safe (a * b) % mod using 128-bit intermediate
8 static inline u64 mul_mod(u64 a, u64 b, u64 mod) {
9     return (u64)((u128)a * b % mod);
10 }
11 
12 // Fast modular exponentiation: (a^e) % mod
13 static inline u64 pow_mod(u64 a, u64 e, u64 mod) {
14     u64 res = 1;
15     a %= mod;
16     while (e) {
17         if (e & 1) res = mul_mod(res, a, mod);
18         a = mul_mod(a, a, mod);
19         e >>= 1;
20     }
21     return res;
22 }
23 
24 // One Miller–Rabin round for a given base 'a'; assumes n > 2 and odd
25 static bool miller_rabin_round(u64 n, u64 a, u64 d, int s) {
26     if (a % n == 0) return true; // base equals n: trivial pass
27     u64 x = pow_mod(a, d, n);
28     if (x == 1 || x == n - 1) return true;
29     for (int r = 1; r < s; ++r) {
30         x = mul_mod(x, x, n);
31         if (x == n - 1) return true;
32     }
33     return false; // 'a' is a witness: composite
34 }
35 
36 // Deterministic for 64-bit using a strong base set.
37 bool isPrime64(u64 n) {
38     if (n < 2) return false;
39     // Small primes quick checks
40     static const u64 small_primes[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
41     for (u64 p : small_primes) {
42         if (n == p) return true;
43         if (n % p == 0) return (n == p);
44     }
45     if ((n & 1ULL) == 0ULL) return false; // even numbers > 2 are composite
46 
47     // Write n-1 = 2^s * d with d odd
48     u64 d = n - 1;
49     int s = 0;
50     while ((d & 1ULL) == 0ULL) { d >>= 1; ++s; }
51 
52     // Deterministic base set for 64-bit (superset; safe and simple)
53     static const u64 bases[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
54     for (u64 a : bases) {
55         if (!miller_rabin_round(n, a, d, s)) return false; // witness found
56     }
57     return true; // passed all rounds
58 }
59 
60 int main() {
61     ios::sync_with_stdio(false);
62     cin.tie(nullptr);
63 
64     vector<unsigned long long> tests = {
65         2ULL, 3ULL, 4ULL, 5ULL, 17ULL, 19ULL,
66         561ULL,      // Carmichael (composite)
67         1000000007ULL, // prime
68         18446744073709551557ULL // a known 64-bit prime near 2^64
69     };
70 
71     for (auto n : tests) {
72         cout << n << ": " << (isPrime64(n) ? "prime" : "composite") << "\n";
73     }
74     return 0;
75 }
76

This program implements a deterministic Miller–Rabin test for unsigned 64-bit integers. It first filters out trivial cases with small primes and even numbers, then writes n − 1 = 2^s × d. For each base in a proven set for 64-bit inputs, it runs one Miller–Rabin round using overflow-safe modular arithmetic with 128-bit intermediates. Passing all rounds implies primality for all 64-bit n with this base selection.

Time: O(B log n), where B is the number of bases (here, 12)Space: O(1)

Randomized Miller–Rabin with tunable error (probable prime)

1 #include <bits/stdc++.h>
2 using namespace std;
3 using u128 = unsigned __int128;
4 using u64 = uint64_t;
5 
6 static inline u64 mul_mod(u64 a, u64 b, u64 mod) {
7     return (u64)((u128)a * b % mod);
8 }
9 
10 static inline u64 pow_mod(u64 a, u64 e, u64 mod) {
11     u64 res = 1;
12     a %= mod;
13     while (e) {
14         if (e & 1ULL) res = mul_mod(res, a, mod);
15         a = mul_mod(a, a, mod);
16         e >>= 1ULL;
17     }
18     return res;
19 }
20 
21 static bool miller_rabin_round(u64 n, u64 a, u64 d, int s) {
22     if (a % n == 0) return true;
23     u64 x = pow_mod(a, d, n);
24     if (x == 1 || x == n - 1) return true;
25     for (int r = 1; r < s; ++r) {
26         x = mul_mod(x, x, n);
27         if (x == n - 1) return true;
28     }
29     return false;
30 }
31 
32 bool isProbablePrime(u64 n, int k = 8) {
33     if (n < 2) return false;
34     // Quick filters
35     static const u64 small_primes[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
36     for (u64 p : small_primes) {
37         if (n == p) return true;
38         if (n % p == 0) return (n == p);
39     }
40     if ((n & 1ULL) == 0ULL) return false;
41 
42     // n-1 = 2^s * d
43     u64 d = n - 1;
44     int s = 0;
45     while ((d & 1ULL) == 0ULL) { d >>= 1; ++s; }
46 
47     // Random bases; each round reduces error by at least 4x
48     static thread_local mt19937_64 rng((uint64_t)chrono::high_resolution_clock::now().time_since_epoch().count());
49     uniform_int_distribution<u64> dist(2ULL, n - 2ULL);
50 
51     for (int i = 0; i < k; ++i) {
52         u64 a = dist(rng);
53         if (std::gcd(a, n) != 1ULL) return false; // found a nontrivial gcd
54         if (!miller_rabin_round(n, a, d, s)) return false;
55     }
56     return true; // probable prime with error <= (1/4)^k
57 }
58 
59 int main() {
60     ios::sync_with_stdio(false);
61     cin.tie(nullptr);
62 
63     vector<unsigned long long> tests = {3ULL, 4ULL, 5ULL, 17ULL, 341ULL, 561ULL, 1105ULL, 1'000'000'007ULL};
64     for (auto n : tests) {
65         bool pr = isProbablePrime(n, 10); // set rounds k
66         cout << n << ": " << (pr ? "probable prime" : "composite") << "\n";
67     }
68     return 0;
69 }
70

This variant selects k random bases in [2, n − 2] and runs Miller–Rabin for each, short-circuiting on any witness. Each independent round reduces the error probability by at least a factor of 4, so k = 10 gives error ≤ 1/4^{10} ≈ 1e−6. It includes small-prime filters and a quick gcd(a, n) check that may immediately find a factor.

Time: O(k log n) word operations (or O(k (log n)^3) bit complexity)Space: O(1)

Next prime >= n using small-prime filters + Miller–Rabin + wheel stepping

1 #include <bits/stdc++.h>
2 using namespace std;
3 using u128 = unsigned __int128;
4 using u64 = uint64_t;
5 
6 static inline u64 mul_mod(u64 a, u64 b, u64 mod) { return (u64)((u128)a * b % mod); }
7 static inline u64 pow_mod(u64 a, u64 e, u64 mod) {
8     u64 r = 1; a %= mod;
9     while (e) { if (e & 1ULL) r = mul_mod(r, a, mod); a = mul_mod(a, a, mod); e >>= 1ULL; }
10     return r;
11 }
12 static bool mr_round(u64 n, u64 a, u64 d, int s) {
13     if (a % n == 0) return true;
14     u64 x = pow_mod(a, d, n);
15     if (x == 1 || x == n - 1) return true;
16     for (int r = 1; r < s; ++r) { x = mul_mod(x, x, n); if (x == n - 1) return true; }
17     return false;
18 }
19 static bool isPrime64(u64 n) {
20     if (n < 2) return false;
21     static const u64 sp[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
22     for (u64 p : sp) { if (n == p) return true; if (n % p == 0) return (n == p); }
23     if ((n & 1ULL) == 0ULL) return false;
24     u64 d = n - 1; int s = 0; while ((d & 1ULL) == 0ULL) { d >>= 1; ++s; }
25     for (u64 a : sp) if (!mr_round(n, a, d, s)) return false;
26     return true;
27 }
28 
29 // Wheel skipping multiples of 2,3,5: step pattern of length 8 with jumps {2,4,2,4,6,2,6,4}
30 static const int WHEEL[8] = {2,4,2,4,6,2,6,4};
31 
32 unsigned long long nextPrime(unsigned long long n) {
33     if (n <= 2ULL) return 2ULL;
34     // Make n odd and not divisible by 3 or 5 by finding the next candidate on the wheel
35     // Start from the next odd
36     if ((n & 1ULL) == 0ULL) ++n;
37     // Align to the wheel position based on n % 30
38     static const bool isBadMod30[30] = {
39         // mark residues divisible by 2,3,5
40         true,true,false,true,false,true,false,true,false,true,true,true,false,true,false,true,true,true,false,true,false,true,true,true,false,true,true,true,false,true
41     };
42     while (n % 30ULL < 30ULL && isBadMod30[n % 30ULL]) ++n; // move to a 30k+{1,7,11,13,17,19,23,29}
43 
44     // Iterate candidates using the wheel steps
45     int idx = 0;
46     while (true) {
47         if (isPrime64(n)) return n;
48         n += WHEEL[idx];
49         idx = (idx + 1) & 7; // modulo 8
50     }
51 }
52 
53 int main(){
54     ios::sync_with_stdio(false);
55     cin.tie(nullptr);
56 
57     unsigned long long start = 18446744073709551500ULL; // near 2^64
58     auto p = nextPrime(start);
59     cout << "nextPrime(" << start << ") = " << p << "\n";
60 
61     // Demonstrate for several values
62     for (unsigned long long x : {1ULL, 2ULL, 14ULL, 1000ULL, 1'000'000ULL}) {
63         cout << "nextPrime(" << x << ") = " << nextPrime(x) << "\n";
64     }
65     return 0;
66 }
67

This utility finds the next prime ≥ n efficiently. It uses a 2·3·5 wheel to skip obvious composites and a deterministic 64-bit Miller–Rabin for correctness. The wheel reduces the density of candidates by filtering numbers divisible by 2, 3, or 5 before invoking the heavier primality test.

Time: Average O((log n) · 1/π-density) per candidate; overall near O(log n) per tested candidate, with candidate density ≈ 8/30.Space: O(1)

1	#include <bits/stdc++.h>
2	using namespace std;
3
4	using u128 = unsigned __int128;
5	using u64 = uint64_t;
6
7	// Safe (a * b) % mod using 128-bit intermediate
8	static inline u64 mul_mod(u64 a, u64 b, u64 mod) {
9	return (u64)((u128)a * b % mod);
10	}
11
12	// Fast modular exponentiation: (a^e) % mod
13	static inline u64 pow_mod(u64 a, u64 e, u64 mod) {
14	u64 res = 1;
15	a %= mod;
16	while (e) {
17	if (e & 1) res = mul_mod(res, a, mod);
18	a = mul_mod(a, a, mod);
19	e >>= 1;
20	}
21	return res;
22	}
23
24	// One Miller–Rabin round for a given base 'a'; assumes n > 2 and odd
25	static bool miller_rabin_round(u64 n, u64 a, u64 d, int s) {
26	if (a % n == 0) return true; // base equals n: trivial pass
27	u64 x = pow_mod(a, d, n);
28	if (x == 1 \|\| x == n - 1) return true;
29	for (int r = 1; r < s; ++r) {
30	x = mul_mod(x, x, n);
31	if (x == n - 1) return true;
32	}
33	return false; // 'a' is a witness: composite
34	}
35
36	// Deterministic for 64-bit using a strong base set.
37	bool isPrime64(u64 n) {
38	if (n < 2) return false;
39	// Small primes quick checks
40	static const u64 small_primes[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
41	for (u64 p : small_primes) {
42	if (n == p) return true;
43	if (n % p == 0) return (n == p);
44	}
45	if ((n & 1ULL) == 0ULL) return false; // even numbers > 2 are composite
46
47	// Write n-1 = 2^s * d with d odd
48	u64 d = n - 1;
49	int s = 0;
50	while ((d & 1ULL) == 0ULL) { d >>= 1; ++s; }
51
52	// Deterministic base set for 64-bit (superset; safe and simple)
53	static const u64 bases[] = {2ULL,3ULL,5ULL,7ULL,11ULL,13ULL,17ULL,19ULL,23ULL,29ULL,31ULL,37ULL};
54	for (u64 a : bases) {
55	if (!miller_rabin_round(n, a, d, s)) return false; // witness found
56	}
57	return true; // passed all rounds
58	}
59
60	int main() {
61	ios::sync_with_stdio(false);
62	cin.tie(nullptr);
63
64	vector<unsigned long long> tests = {
65	2ULL, 3ULL, 4ULL, 5ULL, 17ULL, 19ULL,
66	561ULL, // Carmichael (composite)
67	1000000007ULL, // prime
68	18446744073709551557ULL // a known 64-bit prime near 2^64
69	};
70
71	for (auto n : tests) {
72	cout << n << ": " << (isPrime64(n) ? "prime" : "composite") << "\n";
73	}
74	return 0;
75	}
76