Documentation

Mathlib.NumberTheory.LucasLehmer

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 #

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
    theorem mersenne_pos {p : } (h : 0 < p) :
    theorem one_lt_mersenne {p : } (hp : 1 < p) :
    @[simp]
    theorem succ_mersenne (k : ) :
    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 : ) :
          theorem LucasLehmer.sMod_lt (p : ) (hp : p 0) (i : ) :
          LucasLehmer.sMod p i < 2 ^ p - 1
          theorem LucasLehmer.sZMod_eq_s (p' : ) (i : ) :
          theorem LucasLehmer.Int.coe_nat_pow_pred (b : ) (p : ) (w : 0 < b) :
          (b ^ p - 1) = b ^ p - 1
          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 {q : ℕ+} {x : LucasLehmer.X q} {y : LucasLehmer.X q} (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 : LucasLehmer.X q) (y : LucasLehmer.X q) :
                  (x + y).1 = x.1 + y.1
                  @[simp]
                  theorem LucasLehmer.X.add_snd {q : ℕ+} (x : LucasLehmer.X q) (y : LucasLehmer.X q) :
                  (x + y).2 = x.2 + y.2
                  @[simp]
                  theorem LucasLehmer.X.neg_fst {q : ℕ+} (x : LucasLehmer.X q) :
                  (-x).1 = -x.1
                  @[simp]
                  theorem LucasLehmer.X.neg_snd {q : ℕ+} (x : LucasLehmer.X q) :
                  (-x).2 = -x.2
                  Equations
                  • LucasLehmer.X.instMulX = { mul := fun (x y : LucasLehmer.X q) => (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 : LucasLehmer.X q) (y : LucasLehmer.X q) :
                  (x * y).1 = x.1 * y.1 + 3 * x.2 * y.2
                  @[simp]
                  theorem LucasLehmer.X.mul_snd {q : ℕ+} (x : LucasLehmer.X q) (y : LucasLehmer.X q) :
                  (x * y).2 = x.1 * y.2 + x.2 * y.1
                  Equations
                  • LucasLehmer.X.instOneX = { one := (1, 0) }
                  @[simp]
                  theorem LucasLehmer.X.one_fst {q : ℕ+} :
                  1.1 = 1
                  @[simp]
                  theorem LucasLehmer.X.one_snd {q : ℕ+} :
                  1.2 = 0
                  Equations
                  Equations
                  • LucasLehmer.X.instNatCastX = { natCast := fun (n : ) => (n, 0) }
                  @[simp]
                  theorem LucasLehmer.X.nat_coe_fst {q : ℕ+} (n : ) :
                  (n).1 = n
                  @[simp]
                  theorem LucasLehmer.X.nat_coe_snd {q : ℕ+} (n : ) :
                  (n).2 = 0
                  @[simp]
                  theorem LucasLehmer.X.ofNat_snd {q : ℕ+} (n : ) [Nat.AtLeastTwo n] :
                  (OfNat.ofNat n).2 = 0
                  Equations
                  • One or more equations did not get rendered due to their size.
                  theorem LucasLehmer.X.left_distrib {q : ℕ+} (x : LucasLehmer.X q) (y : LucasLehmer.X q) (z : LucasLehmer.X q) :
                  x * (y + z) = x * y + x * z
                  theorem LucasLehmer.X.right_distrib {q : ℕ+} (x : LucasLehmer.X q) (y : LucasLehmer.X q) (z : LucasLehmer.X q) :
                  (x + y) * z = x * z + y * z
                  Equations
                  • One or more equations did not get rendered due to their size.
                  Equations
                  Equations
                  • =
                  @[simp]
                  theorem LucasLehmer.X.int_coe_fst {q : ℕ+} (n : ) :
                  (n).1 = n
                  @[simp]
                  theorem LucasLehmer.X.int_coe_snd {q : ℕ+} (n : ) :
                  (n).2 = 0
                  theorem LucasLehmer.X.coe_mul {q : ℕ+} (n : ) (m : ) :
                  (n * m) = n * m
                  theorem LucasLehmer.X.coe_nat {q : ℕ+} (n : ) :
                  n = n

                  The cardinality of X is q^2.

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

                  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.s 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
                      • LucasLehmer.ωUnit p = { val := LucasLehmer.X.ω, inv := LucasLehmer.X.ωb, val_inv := , inv_val := }
                      Instances For
                        @[simp]
                        theorem LucasLehmer.ωUnit_coe (p : ) :
                        (LucasLehmer.ωUnit 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

                        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]