# The Lucas-Lehmer test for Mersenne primes. #

We define lucasLehmerResidue : Π p : ℕ, ZMod (2^p - 1), and prove lucasLehmerResidue p = 0 → Prime (mersenne p).

We construct a norm_num extension to calculate this residue to certify primality of Mersenne primes using lucas_lehmer_sufficiency.

## TODO #

• Show reverse implication.
• Speed up the calculations using n ≡ (n % 2^p) + (n / 2^p) [MOD 2^p - 1].
• Find some bigger primes!

## History #

This development began as a student project by Ainsley Pahljina, and was then cleaned up for mathlib by Scott Morrison. The tactic for certified computation of Lucas-Lehmer residues was provided by Mario Carneiro. This tactic was ported by Thomas Murrills to Lean 4, and then it was converted to a norm_num extension and made to use kernel reductions by Kyle Miller.

def mersenne (p : ) :

The Mersenne numbers, 2^p - 1.

Equations
Instances For
@[simp]
theorem mersenne_lt_mersenne {p : } {q : } :
p < q
theorem GCongr.mersenne_lt_mersenne {p : } {q : } :
p < q

Alias of the reverse direction of mersenne_lt_mersenne.

@[simp]
theorem mersenne_le_mersenne {p : } {q : } :
p q
theorem GCongr.mersenne_le_mersenne {p : } {q : } :
p q

Alias of the reverse direction of mersenne_le_mersenne.

@[simp]
theorem mersenne_zero :
= 0
@[simp]
theorem mersenne_odd {p : } :
@[simp]
theorem mersenne_pos {p : } :
0 < 0 < p

Alias of the reverse direction of mersenne_pos.

Extension for the positivity tactic: mersenne.

Instances For
@[simp]
theorem one_lt_mersenne {p : } :
1 < 1 < p
@[simp]
theorem succ_mersenne (k : ) :
+ 1 = 2 ^ k

We now define three(!) different versions of the recurrence s (i+1) = (s i)^2 - 2.

These versions take values either in ℤ, in ZMod (2^p - 1), or in ℤ but applying % (2^p - 1) at each step.

They are each useful at different points in the proof, so we take a moment setting up the lemmas relating them.

The recurrence s (i+1) = (s i)^2 - 2 in ℤ.

Equations
Instances For
def LucasLehmer.sZMod (p : ) :
ZMod (2 ^ p - 1)

The recurrence s (i+1) = (s i)^2 - 2 in ZMod (2^p - 1).

Equations
Instances For
def LucasLehmer.sMod (p : ) :

The recurrence s (i+1) = ((s i)^2 - 2) % (2^p - 1) in ℤ.

Equations
Instances For
theorem LucasLehmer.mersenne_int_pos {p : } (hp : p 0) :
0 < 2 ^ p - 1
theorem LucasLehmer.mersenne_int_ne_zero (p : ) (hp : p 0) :
2 ^ p - 1 0
theorem LucasLehmer.sMod_nonneg (p : ) (hp : p 0) (i : ) :
0
theorem LucasLehmer.sMod_mod (p : ) (i : ) :
% (2 ^ p - 1) =
theorem LucasLehmer.sMod_lt (p : ) (hp : p 0) (i : ) :
< 2 ^ p - 1
theorem LucasLehmer.sZMod_eq_s (p' : ) (i : ) :
LucasLehmer.sZMod (p' + 2) i =
theorem LucasLehmer.Int.natCast_pow_pred (b : ) (p : ) (w : 0 < b) :
(b ^ p - 1) = b ^ p - 1
@[deprecated LucasLehmer.Int.natCast_pow_pred]
theorem LucasLehmer.Int.coe_nat_pow_pred (b : ) (p : ) (w : 0 < b) :
(b ^ p - 1) = b ^ p - 1

Alias of LucasLehmer.Int.natCast_pow_pred.

theorem LucasLehmer.Int.coe_nat_two_pow_pred (p : ) :
(2 ^ p - 1) = 2 ^ p - 1

The Lucas-Lehmer residue is s p (p-2) in ZMod (2^p - 1).

Equations
Instances For

Lucas-Lehmer Test: a Mersenne number 2^p-1 is prime if and only if the Lucas-Lehmer residue s p (p-2) % (2^p - 1) is zero.

Equations
Instances For

q is defined as the minimum factor of mersenne p, bundled as an ℕ+.

Equations
Instances For

We construct the ring X q as ℤ/qℤ + √3 ℤ/qℤ.

Equations
Instances For
Equations
Equations
Equations
Equations
theorem LucasLehmer.X.ext_iff {q : ℕ+} {x : } {y : } :
x = y x.1 = y.1 x.2 = y.2
theorem LucasLehmer.X.ext {q : ℕ+} {x : } {y : } (h₁ : x.1 = y.1) (h₂ : x.2 = y.2) :
x = y
@[simp]
theorem LucasLehmer.X.zero_fst {q : ℕ+} :
0.1 = 0
@[simp]
theorem LucasLehmer.X.zero_snd {q : ℕ+} :
0.2 = 0
@[simp]
theorem LucasLehmer.X.add_fst {q : ℕ+} (x : ) (y : ) :
(x + y).1 = x.1 + y.1
@[simp]
theorem LucasLehmer.X.add_snd {q : ℕ+} (x : ) (y : ) :
(x + y).2 = x.2 + y.2
@[simp]
theorem LucasLehmer.X.neg_fst {q : ℕ+} (x : ) :
(-x).1 = -x.1
@[simp]
theorem LucasLehmer.X.neg_snd {q : ℕ+} (x : ) :
(-x).2 = -x.2
instance LucasLehmer.X.instMul {q : ℕ+} :
Equations
• LucasLehmer.X.instMul = { mul := fun (x y : ) => (x.1 * y.1 + 3 * x.2 * y.2, x.1 * y.2 + x.2 * y.1) }
@[simp]
theorem LucasLehmer.X.mul_fst {q : ℕ+} (x : ) (y : ) :
(x * y).1 = x.1 * y.1 + 3 * x.2 * y.2
@[simp]
theorem LucasLehmer.X.mul_snd {q : ℕ+} (x : ) (y : ) :
(x * y).2 = x.1 * y.2 + x.2 * y.1
instance LucasLehmer.X.instOne {q : ℕ+} :
Equations
• LucasLehmer.X.instOne = { one := (1, 0) }
@[simp]
theorem LucasLehmer.X.one_fst {q : ℕ+} :
1.1 = 1
@[simp]
theorem LucasLehmer.X.one_snd {q : ℕ+} :
1.2 = 0
Equations
• LucasLehmer.X.instMonoid = Monoid.mk npowRec
Equations
• LucasLehmer.X.instNatCast = { natCast := fun (n : ) => (n, 0) }
@[simp]
theorem LucasLehmer.X.fst_natCast {q : ℕ+} (n : ) :
(↑n).1 = n
@[simp]
theorem LucasLehmer.X.snd_natCast {q : ℕ+} (n : ) :
(↑n).2 = 0
@[simp]
theorem LucasLehmer.X.ofNat_fst {q : ℕ+} (n : ) [n.AtLeastTwo] :
(OfNat.ofNat n).1 =
@[simp]
theorem LucasLehmer.X.ofNat_snd {q : ℕ+} (n : ) [n.AtLeastTwo] :
(OfNat.ofNat n).2 = 0
Equations
theorem LucasLehmer.X.left_distrib {q : ℕ+} (x : ) (y : ) (z : ) :
x * (y + z) = x * y + x * z
theorem LucasLehmer.X.right_distrib {q : ℕ+} (x : ) (y : ) (z : ) :
(x + y) * z = x * z + y * z
instance LucasLehmer.X.instRing {q : ℕ+} :
Equations
Equations
• LucasLehmer.X.instCommRing =
Equations
• =
@[simp]
theorem LucasLehmer.X.fst_intCast {q : ℕ+} (n : ) :
(↑n).1 = n
@[simp]
theorem LucasLehmer.X.snd_intCast {q : ℕ+} (n : ) :
(↑n).2 = 0
@[deprecated LucasLehmer.X.fst_natCast]
theorem LucasLehmer.X.nat_coe_fst {q : ℕ+} (n : ) :
(↑n).1 = n

Alias of LucasLehmer.X.fst_natCast.

@[deprecated LucasLehmer.X.snd_natCast]
theorem LucasLehmer.X.nat_coe_snd {q : ℕ+} (n : ) :
(↑n).2 = 0

Alias of LucasLehmer.X.snd_natCast.

@[deprecated LucasLehmer.X.fst_intCast]
theorem LucasLehmer.X.int_coe_fst {q : ℕ+} (n : ) :
(↑n).1 = n

Alias of LucasLehmer.X.fst_intCast.

@[deprecated LucasLehmer.X.snd_intCast]
theorem LucasLehmer.X.int_coe_snd {q : ℕ+} (n : ) :
(↑n).2 = 0

Alias of LucasLehmer.X.snd_intCast.

theorem LucasLehmer.X.coe_mul {q : ℕ+} (n : ) (m : ) :
(n * m) = n * m
theorem LucasLehmer.X.coe_natCast {q : ℕ+} (n : ) :
n = n
@[deprecated LucasLehmer.X.coe_natCast]
theorem LucasLehmer.X.coe_nat {q : ℕ+} (n : ) :
n = n

Alias of LucasLehmer.X.coe_natCast.

theorem LucasLehmer.X.card_eq {q : ℕ+} :
= q ^ 2

The cardinality of X is q^2.

theorem LucasLehmer.X.card_units_lt {q : ℕ+} (w : 1 < q) :
< q ^ 2

There are strictly fewer than q^2 units, since 0 is not a unit.

We define ω = 2 + √3.

Equations
• LucasLehmer.X.ω = (2, 1)
Instances For

We define ωb = 2 - √3, which is the inverse of ω.

Equations
• LucasLehmer.X.ωb = (2, -1)
Instances For
theorem LucasLehmer.X.ω_mul_ωb (q : ℕ+) :
LucasLehmer.X.ω * LucasLehmer.X.ωb = 1
theorem LucasLehmer.X.ωb_mul_ω (q : ℕ+) :
LucasLehmer.X.ωb * LucasLehmer.X.ω = 1
theorem LucasLehmer.X.closed_form {q : ℕ+} (i : ) :
= LucasLehmer.X.ω ^ 2 ^ i + LucasLehmer.X.ωb ^ 2 ^ i

A closed form for the recurrence relation.

Here and below, we introduce p' = p - 2, in order to avoid using subtraction in ℕ.

theorem LucasLehmer.two_lt_q (p' : ) :
2 < LucasLehmer.q (p' + 2)

If 1 < p, then q p, the smallest prime factor of mersenne p, is more than 2.

theorem LucasLehmer.ω_pow_formula (p' : ) (h : lucasLehmerResidue (p' + 2) = 0) :
∃ (k : ), LucasLehmer.X.ω ^ 2 ^ (p' + 1) = k * (mersenne (p' + 2)) * LucasLehmer.X.ω ^ 2 ^ p' - 1
theorem LucasLehmer.mersenne_coe_X (p : ) :
(mersenne p) = 0

q is the minimum factor of mersenne p, so M p = 0 in X q.

theorem LucasLehmer.ω_pow_eq_neg_one (p' : ) (h : lucasLehmerResidue (p' + 2) = 0) :
LucasLehmer.X.ω ^ 2 ^ (p' + 1) = -1
theorem LucasLehmer.ω_pow_eq_one (p' : ) (h : lucasLehmerResidue (p' + 2) = 0) :
LucasLehmer.X.ω ^ 2 ^ (p' + 2) = 1

ω as an element of the group of units.

Equations
• = { val := LucasLehmer.X.ω, inv := LucasLehmer.X.ωb, val_inv := , inv_val := }
Instances For
@[simp]
theorem LucasLehmer.ωUnit_coe (p : ) :
= LucasLehmer.X.ω
theorem LucasLehmer.order_ω (p' : ) (h : lucasLehmerResidue (p' + 2) = 0) :
orderOf (LucasLehmer.ωUnit (p' + 2)) = 2 ^ (p' + 2)

The order of ω in the unit group is exactly 2^p.

theorem LucasLehmer.order_ineq (p' : ) (h : lucasLehmerResidue (p' + 2) = 0) :
2 ^ (p' + 2) < (LucasLehmer.q (p' + 2)) ^ 2
theorem lucas_lehmer_sufficiency (p : ) (w : 1 < p) :

### norm_num extension #

Next we define a norm_num extension that calculates LucasLehmerTest p for 1 < p. It makes use of a version of sMod that is specifically written to be reducible by the Lean 4 kernel, which has the capability of efficiently reducing natural number expressions. With this reduction in hand, it's a simple matter of applying the lemma LucasLehmer.residue_eq_zero_iff_sMod_eq_zero.

See [Archive/Examples/MersennePrimes.lean] for certifications of all Mersenne primes up through mersenne 4423.

Version of sMod that is ℕ-valued. One should have q = 2 ^ p - 1. This can be reduced by the kernel.

Equations
Instances For

Calculate LucasLehmer.LucasLehmerTest p for 2 ≤ p by using kernel reduction for the sMod' function.

Equations
• One or more equations did not get rendered due to their size.
Instances For

This implementation works successfully to prove (2^4423 - 1).Prime, and all the Mersenne primes up to this point appear in [Archive/Examples/MersennePrimes.lean]. These can be calculated nearly instantly, and (2^9689 - 1).Prime only fails due to deep recursion.

(Note by kmill: the following notes were for the Lean 3 version. They seem like they could still be useful, so I'm leaving them here.)

There's still low hanging fruit available to do faster computations based on the formula

n ≡ (n % 2^p) + (n / 2^p) [MOD 2^p - 1]


and the fact that % 2^p and / 2^p can be very efficient on the binary representation. Someone should do this, too!

theorem modEq_mersenne (n : ) (k : ) :
k k / 2 ^ n + k % 2 ^ n [MOD 2 ^ n - 1]