831 Fundamental

Miller-Rabin Primality Test

Functional Programming

Tutorial Video

Text description (accessibility)

This video demonstrates the "Miller-Rabin Primality Test" functional Rust example. Difficulty level: Fundamental. Key concepts covered: Functional Programming. Deterministic primality testing via trial division is O(sqrt(n)), impractical for 64-bit numbers (up to 4 billion operations). Key difference from OCaml: | Aspect | Rust | OCaml |

Tutorial

The Problem

Deterministic primality testing via trial division is O(sqrt(n)), impractical for 64-bit numbers (up to 4 billion operations). Cryptographic applications need to test numbers with hundreds of digits for primality during key generation. Miller-Rabin is a probabilistic primality test that runs in O(k log^2 n) where k is the number of witness rounds — each round reduces error probability to 1/4. With k=25 rounds, the probability of a false positive is negligible. For 64-bit integers, a specific set of deterministic witnesses {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37} makes Miller-Rabin deterministic — no false positives for any n < 3.3 * 10^24.

🎯 Learning Outcomes

• Decompose n-1 as 2^r * d where d is odd, used in the Miller-Rabin witness test

• Implement the witness test: a^d ≠ 1 (mod n) AND a^(2^j * d) ≠ -1 (mod n) for all j means composite

• Understand strong pseudoprimes and why multiple witnesses reduce false positive rate

• Use the deterministic witness set for 64-bit integers to get an exact answer

• Compare with Fermat test: Miller-Rabin has no Carmichael number problem

Code Example

#![allow(clippy::all)]
//! Placeholder

(* Miller-Rabin Primality Test in OCaml *)

(* Modular multiplication using Int64 to avoid overflow *)
let mulmod a b m =
  Int64.(to_int (rem (mul (of_int a) (of_int b)) (of_int m)))

(* Fast modular exponentiation *)
let rec pow_mod base exp m =
  if exp = 0 then 1
  else if exp land 1 = 1 then mulmod base (pow_mod base (exp-1) m) m
  else
    let half = pow_mod base (exp/2) m in
    mulmod half half m

(* Single Miller-Rabin witness test: is 'a' a witness to n's primality? *)
(* Returns true if n is probably prime w.r.t. witness a *)
let miller_witness (n : int) (d : int) (s : int) (a : int) : bool =
  let x = ref (pow_mod a d n) in
  if !x = 1 || !x = n - 1 then true
  else begin
    let composite = ref true in
    for _ = 1 to s - 1 do
      x := mulmod !x !x n;
      if !x = n - 1 then composite := false
    done;
    not !composite
  end

(* Deterministic Miller-Rabin for n < 3.2*10^18 *)
let is_prime (n : int) : bool =
  if n < 2 then false
  else if n = 2 || n = 3 || n = 5 || n = 7 then true
  else if n mod 2 = 0 || n mod 3 = 0 then false
  else begin
    (* Factor n-1 = 2^s * d, d odd *)
    let d = ref (n - 1) and s = ref 0 in
    while !d mod 2 = 0 do d := !d / 2; incr s done;
    (* Witnesses sufficient for n < 3.2*10^18 *)
    let witnesses = [2; 3; 5; 7; 11; 13; 17; 19; 23; 29; 31; 37] in
    List.for_all (fun a ->
      if a >= n then true  (* witness >= n means n is prime by definition *)
      else miller_witness n !d !s a
    ) witnesses
  end

let () =
  let primes   = [2; 3; 5; 17; 97; 997; 1_000_003; 999_999_937] in
  let composites = [4; 9; 15; 100; 1001; 1_000_000] in
  List.iter (fun n ->
    Printf.printf "is_prime(%d) = %b  (expected true)\n" n (is_prime n)
  ) primes;
  List.iter (fun n ->
    Printf.printf "is_prime(%d) = %b  (expected false)\n" n (is_prime n)
  ) composites;
  (* Carmichael number: 561 = 3*11*17, fools Fermat test *)
  Printf.printf "is_prime(561) = %b  (Carmichael, expected false)\n" (is_prime 561)

Key Differences

Aspect	Rust	OCaml
64-bit modular mul	`as u128` intermediate	`Int64` or `Zarith`
Witness iteration	`for &a in &[...]`	`List.for_all`
Overflow safety	`u128` prevents overflow	`Zarith` arbitrary precision
Deterministic range	n < 3.3 * 10^24	Same witnesses, same range
False positives	None with chosen witnesses	None
Probabilistic mode	Add random witnesses	`Random.int64` witnesses

OCaml Approach

OCaml's Miller-Rabin decomposes n-1 using a recursive/iterative bit-strip. The witness list is a let witnesses = [2; 3; 5; 7; 11; 13; 17; 19; 23; 29; 31; 37]. The test uses List.for_all (fun a -> not (is_composite n a d r)) witnesses. OCaml's Int64 or Zarith handles modular multiplication for large n. The 63-bit native int suffices for n < 2^62; for full 64-bit range, Int64 or Zarith.mpz is needed.

Full Source

#![allow(clippy::all)]
//! Placeholder

(* Miller-Rabin Primality Test in OCaml *)

(* Modular multiplication using Int64 to avoid overflow *)
let mulmod a b m =
  Int64.(to_int (rem (mul (of_int a) (of_int b)) (of_int m)))

(* Fast modular exponentiation *)
let rec pow_mod base exp m =
  if exp = 0 then 1
  else if exp land 1 = 1 then mulmod base (pow_mod base (exp-1) m) m
  else
    let half = pow_mod base (exp/2) m in
    mulmod half half m

(* Single Miller-Rabin witness test: is 'a' a witness to n's primality? *)
(* Returns true if n is probably prime w.r.t. witness a *)
let miller_witness (n : int) (d : int) (s : int) (a : int) : bool =
  let x = ref (pow_mod a d n) in
  if !x = 1 || !x = n - 1 then true
  else begin
    let composite = ref true in
    for _ = 1 to s - 1 do
      x := mulmod !x !x n;
      if !x = n - 1 then composite := false
    done;
    not !composite
  end

(* Deterministic Miller-Rabin for n < 3.2*10^18 *)
let is_prime (n : int) : bool =
  if n < 2 then false
  else if n = 2 || n = 3 || n = 5 || n = 7 then true
  else if n mod 2 = 0 || n mod 3 = 0 then false
  else begin
    (* Factor n-1 = 2^s * d, d odd *)
    let d = ref (n - 1) and s = ref 0 in
    while !d mod 2 = 0 do d := !d / 2; incr s done;
    (* Witnesses sufficient for n < 3.2*10^18 *)
    let witnesses = [2; 3; 5; 7; 11; 13; 17; 19; 23; 29; 31; 37] in
    List.for_all (fun a ->
      if a >= n then true  (* witness >= n means n is prime by definition *)
      else miller_witness n !d !s a
    ) witnesses
  end

let () =
  let primes   = [2; 3; 5; 17; 97; 997; 1_000_003; 999_999_937] in
  let composites = [4; 9; 15; 100; 1001; 1_000_000] in
  List.iter (fun n ->
    Printf.printf "is_prime(%d) = %b  (expected true)\n" n (is_prime n)
  ) primes;
  List.iter (fun n ->
    Printf.printf "is_prime(%d) = %b  (expected false)\n" n (is_prime n)
  ) composites;
  (* Carmichael number: 561 = 3*11*17, fools Fermat test *)
  Printf.printf "is_prime(561) = %b  (Carmichael, expected false)\n" (is_prime 561)

Deep Comparison

OCaml vs Rust: Miller Rabin Primality

Overview

See the example.rs and example.ml files for detailed implementations.

Key Differences

Aspect	OCaml	Rust
Type system	Hindley-Milner	Ownership + traits
Memory	GC	Zero-cost abstractions
Mutability	Explicit ref	mut keyword
Error handling	Option/Result	Result<T, E>

See README.md for detailed comparison.

Exercises

Implement the probabilistic version with k random witnesses and estimate false positive rate empirically.

Use Miller-Rabin to generate large primes: generate random odd n until Miller-Rabin returns true.

Compare Miller-Rabin with the Fermat primality test on Carmichael numbers (which fool Fermat but not Miller-Rabin).

Implement Baillie-PSW primality test (Miller-Rabin base 2 + strong Lucas) which has no known pseudoprimes.

Benchmark Miller-Rabin vs. trial division for n up to 10^12 and measure crossover point.

Open Source Repos

functional-rust

View the source for this example on GitHub — OCaml and Rust side by side in the repo.

Rust