El Pop

Unlike the Sieve of Eratosthenes, which finds all primes in a range, Miller-Rabin is used to check if a single, specific number is prime. It is a probabilistic test, meaning it can tell you for certain if a number is composite, but it can only tell you that a number is “probably” prime.

It is the industry standard for generating large primes (e.g., for RSA encryption) because it is computationally efficient even for numbers with hundreds of digits.

The Mathematical Basis

It is based on an extension of Fermat’s Little Theorem, which states that if $n$ is prime, then for any a where $1 < a < n-1: a^{n-1} \equiv 1 \pmod n$

Miller-Rabin goes further by checking for “non-trivial square roots of unity.” In a prime field, the only square roots of $1 \pmod n$ are 1 and -1 (which is $n-1$). If we find any other square root of 1, the number is definitely composite.

How the Algorithm Works

To test an odd number $n > 2$:

1. Decomposition: Write $n-1 as 2^s \cdot d$, where d is odd. (Essentially, you keep dividing $n-1$ by 2 until you can’t anymore).
2. Witness Picking: Pick a random integer a (the “witness”) in the range [2, n-2].
3. Initial Compute: Compute $x = a^d \pmod n$.
   • If $x = 1 or x = n-1, n$ passes this round. It is “probably prime.”
4. Squaring Loop: If the previous step didn’t pass, square x repeatedly s-1 times:
   • Compute x = $x^2 \pmod n$.
   • If $x = n-1$, $n$ passes. It is “probably prime.”
   • If $x = 1$ (without having hit $n-1$ first), it’s a non-trivial square root of unity $\rightarrow$ Definitely Composite.
5. Failure: If the loop finishes and you never hit n-1, the number is Definitely Composite.

Operational Reality & Trade-offs

1. The Error Rate

Each round with a different random a reduces the probability that a composite number is misidentified as prime. The probability of a composite number passing k rounds is less than $(1/4)^k$. With 40 rounds, the chance of error is lower than the chance of a cosmic ray flipping a bit in your CPU and causing a crash.

2. Deterministic Miller-Rabin

For smaller numbers, we don’t need randomness. It has been proven that if you test a specific set of bases, the test becomes 100% deterministic:

• If $n < 4,759,123,141$, you only need to test a = 2, 7, 61.
• If $n < 3,317,044,064,679,887,361$, you only need to test $a = 2, 3, 5, 7, 11, 13, 17, 19, 23.$

3. Complexity

The time complexity is $O(k \log^3 n)$, where k is the number of bases tested. This is incredibly fast compared to trial division or sieving for a single large number.

Summary Comparison

Feature	Sieve (Basic/Segmented)	Miller-Rabin
Goal	Find all primes in a range	Check one specific number
Certainty	100% Deterministic	Probabilistic (by default)
Memory	$O(\sqrt{n})$(Segmented)	O(1)
Speed for 10^{12}	Very Slow (must sieve all)	Instant (modular exponentiation)
Best use case	Pre-computing prime tables	Generating huge crypto keys

C program implementation

To implement Miller-Rabin for 64-bit integers in C, the biggest challenge is overflow during modular multiplication. When you calculate $(a \times b) \pmod n$, if both a and b are near $2^{64}$, the result exceeds the capacity of a uint64_t.

I have used __int128_t (a GCC/Clang extension) to handle these intermediate products safely. I have also used a deterministic set of bases, meaning this program is 100% accurate for any number up to $2^{64}-1$; it is not probabilistic for this range.

miller_rabin.c

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>

typedef __uint128_t uint128;

// Modular exponentiation: (base^exp) % mod
uint64_t power(uint64_t base, uint64_t exp, uint64_t mod) {
    uint64_t res = 1;
    base %= mod;
    while (exp > 0) {
        if (exp % 2 == 1) res = (uint64_t)((uint128)res * base % mod);
        base = (uint64_t)((uint128)base * base % mod);
        exp /= 2;
    }
    return res;
}

// The Miller-Rabin test for a specific witness 'a'
int miller_test(uint64_t n, uint64_t a) {
    if (a % n == 0) return 1; // a is a multiple of n, not a useful witness

    uint64_t d = n - 1;
    int s = 0;
    while (d % 2 == 0) {
        d /= 2;
        s++;
    }

    uint64_t x = power(a, d, n);
    if (x == 1 || x == n - 1) return 1;

    for (int r = 1; r < s; r++) {
        x = (uint64_t)((uint128)x * x % n);
        if (x == n - 1) return 1;
    }
    return 0; // Definitely composite
}

// Deterministic primality check for uint64_t
int is_prime(uint64_t n) {
    if (n < 2) return 0;
    if (n == 2 || n == 3) return 1;
    if (n % 2 == 0) return 0;

    // Deterministic bases for n < 2^64
    // These bases are proven to be sufficient for all 64-bit integers
    uint64_t bases[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
    int num_bases = sizeof(bases) / sizeof(bases[0]);

    for (int i = 0; i < num_bases; i++) {
        if (n == bases[i]) return 1;
        if (!miller_test(n, bases[i])) return 0;
    }
    return 1;
}

int main(int argc, char *argv[]) {
    if (argc < 2) {
        printf("Usage: %s <number>\n", argv[0]);
        return 1;
    }

    uint64_t n = strtoull(argv[1], NULL, 10);
    if (is_prime(n)) {
        printf("%llu is PRIME\n", (unsigned long long)n);
    } else {
        printf("%llu is COMPOSITE\n", (unsigned long long)n);
    }

    return 0;
}Code language: C++ (cpp)

Compilation and Execution

gcc -O3 miller_rabin.c -o miller_rabin
./miller_rabin 18446744073709551557  # Largest 64-bit prime

Engineering Breakdown:

1. __uint128_t: This is the most pragmatic way to implement Miller-Rabin in C. Without it, you would have to implement a custom mul_mod function that performs multiplication via addition (similar to binary exponentiation), which is significantly slower.
2. Deterministic Bases: I used the bases $\{2, 3, 5, \dots, 37\}$. This specific set is proven to be deterministic for all $n < 2^{64}$. You don’t need to worry about “probably prime” vs “definitely prime” for standard 64-bit unsigned integers.
3. Complexity: The time complexity is $O(k \log^3 n)$, where k is the number of bases (12 in this case). The check is virtually instantaneous even for the largest possible uint64_t.
4. Edge Case Handling:
   • $n < 2$: Not prime.
   • $n=2, 3$: Prime.
   • Even numbers > 2: Composite.
   • Base a being a multiple of n: Handled in miller_test().