diff --git a/Complexity.lean b/Complexity.lean index e6f116c6e2..3901cfecef 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -88,9 +88,11 @@ public import Complexity.Algorithms.ThreeSum.ChanHe.Sizes public import Complexity.Algorithms.ThreeSum.Convolution public import Complexity.Algorithms.ThreeSum.Defs public import Complexity.Algorithms.ThreeSum.ThreeSets +public import Complexity.Algorithms.WordRAM.AKS.Congruence public import Complexity.Algorithms.WordRAM.AKS.PerfectPower public import Complexity.Algorithms.WordRAM.AKS.PolyArray public import Complexity.Algorithms.WordRAM.AKS.PolyMul +public import Complexity.Algorithms.WordRAM.AKS.PolyPow public import Complexity.Algorithms.WordRAM.Deduplicate.Count public import Complexity.Algorithms.WordRAM.Deduplicate.Model public import Complexity.Algorithms.WordRAM.Deduplicate.Program diff --git a/Complexity/Algorithms/WordRAM/AKS/Congruence.lean b/Complexity/Algorithms/WordRAM/AKS/Congruence.lean new file mode 100644 index 0000000000..9882c7c400 --- /dev/null +++ b/Complexity/Algorithms/WordRAM/AKS/Congruence.lean @@ -0,0 +1,259 @@ +/- +Copyright (c) 2026 Samuel Schlesinger. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Samuel Schlesinger +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.AKS.PolyArray +public import Complexity.Algorithms.WordRAM.AKS.PolyPow + +/-! +# The AKS congruence for one `a` on the word RAM + +Step 5 of the AKS algorithm [Agrawal, Kayal, Saxena, Section 4][AKS2004] tests +`(X + a)^n ≡ X^n + a (mod X^r - 1, n)`. The program writes the coefficients of `X + a` into a base +block `B` and those of `1` into a result block `R`, raises by repeated squaring +(`Cslib.Algorithms.WordRAM.AKS.PolyPow.program`) so that `R` holds `(X + a)^n`, writes the +coefficients of `X^n + a` into the scratch block `T`, and compares `R` with `T`. Both blocks are +reduced, so they agree exactly when the classes agree (`Cslib.AKS.toQuot_eq_toQuot_iff`), which is +the congruence (`Cslib.AKS.congruence_iff_toQuot`). + +## Main definitions + +* `Cslib.Algorithms.WordRAM.AKS.Congruence.program`: the congruence test. +* `Cslib.Algorithms.WordRAM.AKS.Congruence.time`: its instruction bound. + +## Main results + +* `Cslib.Algorithms.WordRAM.AKS.Congruence.exec`: the flag is `1` exactly when + `Cslib.AKS.Congruence n r a` holds, within `time r n` instructions. + +## References + +* [M. Agrawal, N. Kayal, N. Saxena, *PRIMES is in P*][AKS2004] +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.AKS.Congruence + +open Cslib.WordRAM Cslib.AKS Polynomial + +/-- Base address of the result block. -/ +abbrev R : Fin 25 := 0 +/-- Base address of the base block. -/ +abbrev B : Fin 25 := 1 +/-- Base address of the scratch block. -/ +abbrev T : Fin 25 := 2 +/-- The number `r` of coefficients. -/ +abbrev rr : Fin 25 := 3 +/-- The modulus `n`, which is also the exponent. -/ +abbrev nn : Fin 25 := 4 +/-- The exponent of the power. -/ +abbrev e : Fin 25 := 5 +/-- Block argument of the subroutines. -/ +abbrev pa : Fin 25 := 9 +/-- Second block argument of the comparison. -/ +abbrev pb : Fin 25 := 10 +/-- Constant one. -/ +abbrev one : Fin 25 := 19 +/-- The tested value `a`, preserved. -/ +abbrev av : Fin 25 := 21 +/-- The exponent residue for `Fill`. -/ +abbrev er : Fin 25 := 22 +/-- The constant term for `Fill`. -/ +abbrev aa : Fin 25 := 23 +/-- The answer: `1` when the congruence holds. -/ +abbrev flag : Fin 25 := 24 + +/-- The power runs on the first `21` registers. -/ +def powMap : Fin 21 → Fin 25 := Fin.castLE (by norm_num) + +/-- `Fill` writes the block at `pa`, with scratch shared with the power. -/ +def fillMap : Fin 10 → Fin 25 := ![pa, rr, nn, er, aa, 12, 15, 16, one, 20] + +/-- `Compare` compares the blocks at `pa` and `pb`. -/ +def cmpMap : Fin 9 → Fin 25 := ![pa, pb, rr, 12, 15, 16, 17, flag, one] + +theorem powMap_injective : powMap.Injective := Fin.castLE_injective _ +theorem fillMap_injective : fillMap.Injective := by decide +theorem cmpMap_injective : cmpMap.Injective := by decide + +/-- Write `X^e + a` into the block at `pa`. -/ +def fill : Code .standard 25 := Fill.program.rename fillMap + +/-- Raise `R` by `B` to the power in `e`. -/ +def pow : Code .standard 25 := PolyPow.program.rename powMap + +/-- Compare the blocks at `pa` and `pb`. -/ +def cmp : Code .standard 25 := Compare.program.rename cmpMap + +/-- Write `X + a` into `B` and `1` into `R`, raise to the power `n`, write `X^n + a` into `T`, and +compare `R` with `T`. -/ +def program : Code .standard 25 := + .seq (ram[ one ← #1; er ← one % rr; aa ← av; pa ← B ]) (.seq fill + (.seq (ram[ er ← #0; aa ← #0; pa ← R ]) (.seq fill + (.seq (ram[ e ← nn ]) (.seq pow + (.seq (ram[ er ← nn % rr; aa ← av; pa ← T ]) (.seq fill + (.seq (ram[ pa ← R; pb ← T ]) cmp)))))))) + +/-- The instruction bound of `program`. -/ +def time (r n : ℕ) : ℕ := PolyPow.time r n + 26 * r + 57 + +/-- The registers that every subroutine call preserves. -/ +def Keep (j : Fin 25) : Prop := j = R ∨ j = B ∨ j = T ∨ j = rr ∨ j = nn ∨ j = av + +instance : DecidablePred Keep := fun _ => by unfold Keep; infer_instance + +private theorem fill_preserves : ∀ j, Keep j → fill.Preserves j := by decide +private theorem pow_preserves : ∀ j, Keep j → pow.Preserves j := by decide +private theorem cmp_preserves : ∀ j, Keep j → cmp.Preserves j := by decide + +private theorem fill_exec {w : ℕ} (x : Data 25) {bs r n e' A : ℕ} (hb : x.regs pa = bs) + (hr : x.regs rr = r) (hn : x.regs nn = n) (her : x.regs er = e' % r) (ha : x.regs aa = A) + (hr0 : 0 < r) (hn0 : 0 < n) (hfb : Fits w (bs + r)) (hfa : Fits w A) (hfn : Fits w (n + 1)) : + ∃ y, ExecWithin w fill x y (5 * r + 13) ∧ + (∀ k < r, y.mem (bs + k) = powAddCoeff n r e' A k) ∧ + (∀ a', ¬(bs ≤ a' ∧ a' < bs + r) → y.mem a' = x.mem a') ∧ + (∀ j, Keep j → y.regs j = x.regs j) := by + obtain ⟨y, hy, hblk, hframe, -⟩ := Fill.exec (x.comap fillMap) hb hr hn her ha hr0 hn0 hfb hfa hfn + have hy' := hy.rename fillMap_injective + exact ⟨_, hy', hblk, hframe, fun j hj => hy'.regs_eq (fill_preserves j hj)⟩ + +private theorem cmp_exec {w : ℕ} (x : Data 25) {a b r : ℕ} (ha : x.regs pa = a) + (hb : x.regs pb = b) (hr : x.regs rr = r) (hfa : Fits w (a + r)) (hfb : Fits w (b + r)) + (h1 : Fits w 1) : + ∃ y, ExecWithin w cmp x y (11 * r + 5) ∧ y.regs flag ≤ 1 ∧ + (y.regs flag = 1 ↔ ∀ k' < r, x.mem (a + k') = x.mem (b + k')) ∧ + (∀ j, Keep j → y.regs j = x.regs j) := by + obtain ⟨y, hy, hle, hiff, -⟩ := Compare.exec (x.comap cmpMap) ha hb hr hfa hfb h1 + have hy' := hy.rename cmpMap_injective + have hfl : (y.extend cmpMap x).regs flag = y.regs Compare.flag := + cmpMap_injective.extend_apply _ _ Compare.flag + exact ⟨_, hy', by rw [hfl]; exact hle, by rw [hfl]; exact hiff, + fun j hj => hy'.regs_eq (cmp_preserves j hj)⟩ + +private theorem pow_exec {w : ℕ} (x : Data 25) {a b t r n : ℕ} (ha : x.regs R = a) + (hb : x.regs B = b) (ht : x.regs T = t) (hr : x.regs rr = r) (hn : x.regs nn = n) + (hn0 : 0 < n) (hab : a + r ≤ b ∨ b + r ≤ a) (hta : t + r ≤ a ∨ a + r ≤ t) + (htb : t + r ≤ b ∨ b + r ≤ t) (hpa : ∀ i < r, x.mem (a + i) < n) + (hpb : ∀ i < r, x.mem (b + i) < n) (hfa : Fits w (a + r)) (hfb : Fits w (b + r)) + (hft : Fits w (t + r)) (hf2 : Fits w (2 * r)) (hfn : Fits w (n * n + n)) + (hfe : Fits w (x.regs e)) : + ∃ y, ExecWithin w pow x y (PolyPow.time r (x.regs e)) ∧ + toQuot n r (fun i => y.mem (a + i)) = + toQuot n r (fun i => x.mem (a + i)) * toQuot n r (fun i => x.mem (b + i)) ^ x.regs e ∧ + (∀ i < r, y.mem (a + i) < n) ∧ (∀ j, Keep j → y.regs j = x.regs j) := by + obtain ⟨y, hy, hq, hred, -⟩ := PolyPow.exec (x.comap powMap) ha hb ht hr hn hn0 hab hta htb hpa + hpb hfa hfb hft hf2 hfn hfe + have hy' := hy.rename powMap_injective + exact ⟨_, hy', hq, hred, fun j hj => hy'.regs_eq (pow_preserves j hj)⟩ + +/-- The class of `1`, written as `X^0 + 0`. -/ +private theorem toQuot_one {n r : ℕ} (hr : 0 < r) : toQuot n r (powAddCoeff n r 0 0) = 1 := by + rw [toQuot_powAddCoeff hr, pow_zero, Nat.cast_zero, C_0, add_zero, map_one] + +/-- **The AKS congruence for one `a`.** For `r ≥ 1` and `n ≥ 1`, with the result, base and scratch +blocks `[a, a + r)`, `[b, b + r)` and `[t, t + r)` pairwise disjoint, and the word size holding the +block ends, `2 r`, `n² + n` and the tested value `A`, the program leaves `1` in `flag` if +`(X + A)^n ≡ X^n + A (mod X^r - 1, n)` and `0` otherwise, within `time r n` instructions. -/ +theorem exec {w : ℕ} (x : Data 25) {a b t r n A : ℕ} (ha : x.regs R = a) (hb : x.regs B = b) + (ht : x.regs T = t) (hr : x.regs rr = r) (hn : x.regs nn = n) (hA : x.regs av = A) + (hr0 : 0 < r) (hn0 : 0 < n) (hab : a + r ≤ b ∨ b + r ≤ a) (hta : t + r ≤ a ∨ a + r ≤ t) + (htb : t + r ≤ b ∨ b + r ≤ t) (hfa : Fits w (a + r)) (hfb : Fits w (b + r)) + (hft : Fits w (t + r)) (hf2 : Fits w (2 * r)) (hfn : Fits w (n * n + n)) (hfA : Fits w A) : + ∃ y, ExecWithin w program x y (time r n) ∧ y.regs flag ≤ 1 ∧ + (y.regs flag = 1 ↔ Congruence n r A) ∧ + ∀ j, Keep j → y.regs j = x.regs j := by + have h1 : Fits w 1 := hfn.of_le (by nlinarith) + have hfn1 : Fits w (n + 1) := hfn.of_le (by nlinarith) + have hfnn : Fits w n := hfn.of_le (by nlinarith) + have hfr : Fits w r := hf2.of_le (by omega) + -- `B := X + A`. + let x₁ := (((x.set one 1).set er (1 % r)).set aa A).set pa b + have hrun₁ : ExecWithin (B := .standard) w (ram[ one ← #1; er ← one % rr; aa ← av; pa ← B ]) + x x₁ 4 := by + ram_eval [x₁, hr, hA, hb, h1, h1.of_le (Nat.mod_le 1 r)] + have hx₁ (j : Fin 25) (hj : Keep j) : x₁.regs j = x.regs j := by + rcases hj with rfl | rfl | rfl | rfl | rfl | rfl <;> simp +decide [x₁] + obtain ⟨y₁, hrun₂, hB₁, hfr₁, hy₁⟩ := fill_exec (w := w) x₁ (bs := b) (r := r) (n := n) (e' := 1) + (A := A) (by simp +decide [x₁]) (by rw [hx₁ rr (by simp [Keep]), hr]) + (by rw [hx₁ nn (by simp [Keep]), hn]) (by simp +decide [x₁]) (by simp +decide [x₁]) hr0 hn0 + hfb hfA hfn1 + have hy₁' (j : Fin 25) (hj : Keep j) : y₁.regs j = x.regs j := by rw [hy₁ j hj, hx₁ j hj] + -- `R := 1`. + let x₂ := ((y₁.set er 0).set aa 0).set pa a + have hrun₃ : ExecWithin (B := .standard) w (ram[ er ← #0; aa ← #0; pa ← R ]) y₁ x₂ 3 := by + ram_eval [x₂, hy₁' R (by simp [Keep]), ha] + have hx₂ (j : Fin 25) (hj : Keep j) : x₂.regs j = x.regs j := by + rw [← hy₁' j hj] + rcases hj with rfl | rfl | rfl | rfl | rfl | rfl <;> simp +decide [x₂] + obtain ⟨y₂, hrun₄, hR₂, hfr₂, hy₂⟩ := fill_exec (w := w) x₂ (bs := a) (r := r) (n := n) (e' := 0) + (A := 0) (by simp +decide [x₂]) (by rw [hx₂ rr (by simp [Keep]), hr]) + (by rw [hx₂ nn (by simp [Keep]), hn]) (by simp +decide [x₂]) (by simp +decide [x₂]) hr0 hn0 hfa + (fits_zero w) hfn1 + have hy₂' (j : Fin 25) (hj : Keep j) : y₂.regs j = x.regs j := by rw [hy₂ j hj, hx₂ j hj] + have hB₂ (k : ℕ) (hk : k < r) : y₂.mem (b + k) = powAddCoeff n r 1 A k := by + rw [hfr₂ _ (by rcases hab with h | h <;> omega)] + exact hB₁ k hk + -- `R := R · B^n`. + let x₃ := y₂.set e n + have hrun₅ : ExecWithin (B := .standard) w (ram[ e ← nn ]) y₂ x₃ 1 := by + ram_eval [x₃, hy₂' nn (by simp [Keep]), hn] + have hx₃ (j : Fin 25) (hj : Keep j) : x₃.regs j = x.regs j := by + rw [← hy₂' j hj] + rcases hj with rfl | rfl | rfl | rfl | rfl | rfl <;> simp +decide [x₃] + obtain ⟨y₃, hrun₆, hq₃, hred₃, hy₃⟩ := pow_exec (w := w) x₃ (a := a) (b := b) (t := t) (r := r) + (n := n) (by rw [hx₃ R (by simp [Keep]), ha]) (by rw [hx₃ B (by simp [Keep]), hb]) + (by rw [hx₃ T (by simp [Keep]), ht]) (by rw [hx₃ rr (by simp [Keep]), hr]) + (by rw [hx₃ nn (by simp [Keep]), hn]) hn0 hab hta htb + (fun i hi => by rw [show x₃.mem = y₂.mem from rfl, hR₂ i hi]; exact powAddCoeff_lt hn0 _ _ _) + (fun i hi => by rw [show x₃.mem = y₂.mem from rfl, hB₂ i hi]; exact powAddCoeff_lt hn0 _ _ _) + hfa hfb hft hf2 hfn (by simpa [x₃] using hfnn) + have hx₃e : x₃.regs e = n := by simp +decide [x₃] + rw [hx₃e] at hrun₆ hq₃ + have hy₃' (j : Fin 25) (hj : Keep j) : y₃.regs j = x.regs j := by rw [hy₃ j hj, hx₃ j hj] + -- `T := X^n + A`. + let x₄ := ((y₃.set er (n % r)).set aa A).set pa t + have hrun₇ : ExecWithin (B := .standard) w (ram[ er ← nn % rr; aa ← av; pa ← T ]) y₃ x₄ 3 := by + ram_eval [x₄, hy₃' nn (by simp [Keep]), hn, hy₃' rr (by simp [Keep]), hr, + hy₃' av (by simp [Keep]), hA, hy₃' T (by simp [Keep]), ht, + hfnn.of_le (Nat.mod_le _ _)] + have hx₄ (j : Fin 25) (hj : Keep j) : x₄.regs j = x.regs j := by + rw [← hy₃' j hj] + rcases hj with rfl | rfl | rfl | rfl | rfl | rfl <;> simp +decide [x₄] + obtain ⟨y₄, hrun₈, hT₄, hfr₄, hy₄⟩ := fill_exec (w := w) x₄ (bs := t) (r := r) (n := n) (e' := n) + (A := A) (by simp +decide [x₄]) (by rw [hx₄ rr (by simp [Keep]), hr]) + (by rw [hx₄ nn (by simp [Keep]), hn]) (by simp +decide [x₄]) (by simp +decide [x₄]) hr0 hn0 + hft hfA hfn1 + have hy₄' (j : Fin 25) (hj : Keep j) : y₄.regs j = x.regs j := by rw [hy₄ j hj, hx₄ j hj] + have hR₄ (k : ℕ) (hk : k < r) : y₄.mem (a + k) = y₃.mem (a + k) := by + rw [hfr₄ _ (by rcases hta with h | h <;> omega)] + rfl + -- Compare `R` with `T`. + let x₅ := (y₄.set pa a).set pb t + have hrun₉ : ExecWithin (B := .standard) w (ram[ pa ← R; pb ← T ]) y₄ x₅ 2 := by + ram_eval [x₅, hy₄' R (by simp [Keep]), ha, hy₄' T (by simp [Keep]), ht] + have hx₅ (j : Fin 25) (hj : Keep j) : x₅.regs j = x.regs j := by + rw [← hy₄' j hj] + rcases hj with rfl | rfl | rfl | rfl | rfl | rfl <;> simp +decide [x₅] + obtain ⟨y₅, hrun₁₀, hle, hiff, hy₅⟩ := cmp_exec (w := w) x₅ (a := a) (b := t) (r := r) + (by simp +decide [x₅]) (by simp +decide [x₅]) (by rw [hx₅ rr (by simp [Keep]), hr]) hfa hft h1 + refine ⟨y₅, ?_, hle, ?_, fun j hj => by rw [hy₅ j hj, hx₅ j hj]⟩ + · refine (hrun₁.seq (hrun₂.seq (hrun₃.seq (hrun₄.seq (hrun₅.seq (hrun₆.seq (hrun₇.seq + (hrun₈.seq (hrun₉.seq hrun₁₀))))))))).mono ?_ + unfold time + omega + · rw [hiff, congruence_iff_toQuot hr0] + change (∀ k' < r, y₄.mem (a + k') = y₄.mem (t + k')) ↔ _ + have hRred (k : ℕ) (hk : k < r) : y₄.mem (a + k) < n := by rw [hR₄ k hk]; exact hred₃ k hk + have hTred (k : ℕ) (hk : k < r) : y₄.mem (t + k) < n := by + rw [hT₄ k hk]; exact powAddCoeff_lt hn0 _ _ _ + rw [← toQuot_eq_toQuot_iff hn0 hr0 hRred hTred, toQuot_congr hR₄, hq₃, toQuot_congr hT₄, + toQuot_congr (c := fun i => x₃.mem (a + i)) (fun k hk => hR₂ k hk), + toQuot_congr (c := fun i => x₃.mem (b + i)) (fun k hk => hB₂ k hk), toQuot_one hr0, one_mul] + +end Cslib.Algorithms.WordRAM.AKS.Congruence diff --git a/Complexity/Algorithms/WordRAM/AKS/PolyPow.lean b/Complexity/Algorithms/WordRAM/AKS/PolyPow.lean new file mode 100644 index 0000000000..4595ec96e6 --- /dev/null +++ b/Complexity/Algorithms/WordRAM/AKS/PolyPow.lean @@ -0,0 +1,383 @@ +/- +Copyright (c) 2026 Samuel Schlesinger. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Samuel Schlesinger +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.AKS.PolyMul +public import Complexity.Computability.WordRAM.Copy +public import Complexity.Computability.WordRAM.RegisterFrame +public import Mathlib.Data.Nat.Size + +/-! +# Powers modulo `X^r - 1` and `n` on the word RAM + +Raising a polynomial modulo `X^r - 1` and `n` to the power `e` by repeated squaring: with a result +block `R`, a base block `B` and a scratch block `T`, each round multiplies `R` by `B` into `T` and +copies it back when `e` is odd, squares `B` into `T` and copies it back, and halves `e`. The class +of `R · B^e` is invariant, so when `e` reaches `0` the result block holds the class of +`R · B^e` for the original `R`, `B` and `e`. The loop makes `Nat.size e` rounds, each of `O(r²)` +instructions (`Cslib.Algorithms.WordRAM.AKS.PolyMul.program` and +`Cslib.WordRAM.ArrayCopy.program`). This is the polynomial arithmetic of step 5 of the AKS +algorithm [Agrawal, Kayal, Saxena, Section 4][AKS2004]. + +## Main definitions + +* `Cslib.Algorithms.WordRAM.AKS.PolyPow.program`: the repeated-squaring loop. +* `Cslib.Algorithms.WordRAM.AKS.PolyPow.time`: its instruction bound. + +## Main results + +* `Cslib.Algorithms.WordRAM.AKS.PolyPow.exec`: the result block receives the coefficients of + `R · B^e` within `time r e` instructions. + +## References + +* [M. Agrawal, N. Kayal, N. Saxena, *PRIMES is in P*][AKS2004] +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.AKS.PolyPow + +open Cslib.WordRAM Cslib.AKS Polynomial + +/-- Base address of the result block. -/ +abbrev R : Fin 21 := 0 +/-- Base address of the base block. -/ +abbrev B : Fin 21 := 1 +/-- Base address of the scratch block. -/ +abbrev T : Fin 21 := 2 +/-- The number `r` of coefficients. -/ +abbrev rr : Fin 21 := 3 +/-- The modulus `n` of the coefficients. -/ +abbrev nn : Fin 21 := 4 +/-- The exponent, halved every round. -/ +abbrev e : Fin 21 := 5 +/-- The parity of the exponent. -/ +abbrev par : Fin 21 := 6 +/-- Constant two. -/ +abbrev two : Fin 21 := 7 +/-- Constant zero, for the loop test. -/ +abbrev z : Fin 21 := 8 +/-- First argument of the multiplication, and source of the copy. -/ +abbrev pa : Fin 21 := 9 +/-- Second argument of the multiplication, and destination of the copy. -/ +abbrev pb : Fin 21 := 10 +/-- Product block of the multiplication. -/ +abbrev pt : Fin 21 := 11 + +/-- The multiplication runs on `pa`, `pb`, `pt`, `rr`, `nn` and the scratch registers `12, …, 20`. +-/ +def mulMap : Fin 14 → Fin 21 := ![pa, pb, pt, rr, nn, 12, 13, 14, 15, 16, 17, 18, 19, 20] + +/-- The copy runs from `pa` to `pb` over `rr` words, with scratch among `12, …, 20`. -/ +def copyMap : Fin 7 → Fin 21 := ![pa, pb, rr, 12, 15, 16, 19] + +theorem mulMap_injective : mulMap.Injective := by decide +theorem copyMap_injective : copyMap.Injective := by decide + +/-- The multiplication touches no register below `9` other than `rr` and `nn`. -/ +private theorem mulMap_ne {j : Fin 21} (hj : j.val < 9) (hr : j ≠ rr) (hn : j ≠ nn) (a : Fin 14) : + mulMap a ≠ j := by + intro h + subst h + fin_cases a <;> simp_all [mulMap] + +/-- The copy touches no register below `9` other than `rr`. -/ +private theorem copyMap_ne {j : Fin 21} (hj : j.val < 9) (hr : j ≠ rr) (a : Fin 7) : + copyMap a ≠ j := by + intro h + subst h + fin_cases a <;> simp_all [copyMap] + +/-- Multiply the blocks at `x` and `y` into `T`, then copy `T` into the block at `d`. -/ +def mulInto (x y d : Fin 21) : Code .standard 21 := + .seq (ram[ pa ← x; pb ← y; pt ← T ]) (.seq (PolyMul.program.rename mulMap) + (.seq (ram[ pa ← T; pb ← d ]) (ArrayCopy.program.rename copyMap))) + +/-- One round: multiply `R` by `B` if `e` is odd, square `B`, and halve `e`. -/ +def body : Code .standard 21 := ram[ + par ← e % two; if z < par then [ ⟪mulInto R B R⟫ ] else [ par ← par ]; ⟪mulInto B B B⟫; + e ← e / two ] + +/-- Repeated squaring while the exponent is positive. -/ +def program : Code .standard 21 := .seq (ram[ two ← #2; z ← #0 ]) (.loop (.lt z e) body) + +/-- The instruction bound of one multiplication followed by a copy. -/ +def mulTime (r : ℕ) : ℕ := 14 * (r + 1) * (r + 1) + 7 * r + 9 + +/-- The instruction bound of `program`: `Nat.size e` rounds of two multiplications. -/ +def time (r e : ℕ) : ℕ := Nat.size e * (2 * mulTime r + 6) + 4 + +/-- Halving a positive number shortens its binary representation. -/ +theorem size_div_two_lt {m : ℕ} (hm : 0 < m) : Nat.size (m / 2) < Nat.size m := by + have hs : 0 < Nat.size m := Nat.size_pos.mpr hm + have hlt : m < 2 ^ Nat.size m := Nat.lt_size_self m + have : Nat.size (m / 2) ≤ Nat.size m - 1 := by + rw [Nat.size_le] + have h2 : 2 ^ Nat.size m = 2 * 2 ^ (Nat.size m - 1) := by + rw [← pow_succ']; congr 1; omega + omega + omega + +private theorem mulInto_exec {w : ℕ} (X : Data 21) {x y d : Fin 21} (hx : x.val < 9) + (hy : y.val < 9) (hd : d.val < 9) {a b t dd r n : ℕ} (ha : X.regs x = a) (hb : X.regs y = b) + (ht : X.regs T = t) (hdd : X.regs d = dd) (hr : X.regs rr = r) (hn : X.regs nn = n) + (hn0 : 0 < n) (hta : t + r ≤ a ∨ a + r ≤ t) (htb : t + r ≤ b ∨ b + r ≤ t) + (htd : t + r ≤ dd ∨ dd + r ≤ t) + (hpa : ∀ i < r, X.mem (a + i) < n) (hpb : ∀ j < r, X.mem (b + j) < n) + (hfa : Fits w (a + r)) (hfb : Fits w (b + r)) (hft : Fits w (t + r)) (hfd : Fits w (dd + r)) + (hf2 : Fits w (2 * r)) (hfn : Fits w (n * n + n)) : + ∃ Y, ExecWithin w (mulInto x y d) X Y (mulTime r) ∧ + (∀ k < r, Y.mem (dd + k) = + mulCoeff n r (fun i => X.mem (a + i)) (fun j => X.mem (b + j)) k) ∧ + (∀ a', ¬(t ≤ a' ∧ a' < t + r) → ¬(dd ≤ a' ∧ a' < dd + r) → Y.mem a' = X.mem a') ∧ + (∀ j : Fin 21, j.val < 9 → Y.regs j = X.regs j) := by + have hxa : x ≠ pa := fun h => by subst h; simp at hx + have hya : y ≠ pa := fun h => by subst h; simp at hy + have hyb : y ≠ pb := fun h => by subst h; simp at hy + have hda : d ≠ pa := fun h => by subst h; simp at hd + have hdb : d ≠ pb := fun h => by subst h; simp at hd + have h1 : Fits w 1 := hfn.of_le (by nlinarith) + -- Set the arguments of the multiplication. + let X₁ := ((X.set pa a).set pb b).set pt t + have hrun₁ : ExecWithin (B := .standard) w (ram[ pa ← x; pb ← y; pt ← T ]) X X₁ 3 := by + ram_eval [X₁, ha, hb, ht, hxa, hya, hyb] + have hX₁ (j : Fin 21) (hj : j.val < 9) : X₁.regs j = X.regs j := by + have h9 : j ≠ pa ∧ j ≠ pb ∧ j ≠ pt := by + refine ⟨?_, ?_, ?_⟩ <;> (intro h; subst h; simp at hj) + simp [X₁, h9.1, h9.2.1, h9.2.2] + -- Multiply. + obtain ⟨Y₁, hrun₂, hprod, hframe, -, -, -, hY₁r, hY₁n⟩ := + PolyMul.exec (w := w) (X₁.comap mulMap) (a := a) (b := b) (t := t) (r := r) (n := n) + (by change X₁.regs pa = a; simp +decide [X₁]) (by change X₁.regs pb = b; simp +decide [X₁]) + (by change X₁.regs pt = t; simp +decide [X₁]) + (by change X₁.regs rr = r; rw [hX₁ rr (by simp), hr]) + (by change X₁.regs nn = n; rw [hX₁ nn (by simp), hn]) + hn0 hta htb hpa hpb hfa hfb hft hf2 hfn + have hrun₂' := hrun₂.rename mulMap_injective + set Y₁' := Y₁.extend mulMap X₁ with hY₁' + have hY₁'m : Y₁'.mem = Y₁.mem := rfl + have hY₁'r (j : Fin 21) (hj : j.val < 9) : Y₁'.regs j = X.regs j := by + by_cases hjr : j = rr + · subst hjr + rw [hr, hY₁', show rr = mulMap PolyMul.rr from rfl] + exact (mulMap_injective.extend_apply _ _ _).trans hY₁r + by_cases hjn : j = nn + · subst hjn + rw [hn, hY₁', show nn = mulMap PolyMul.nn from rfl] + exact (mulMap_injective.extend_apply _ _ _).trans hY₁n + rw [hrun₂'.regs_eq (Code.preserves_rename_of_forall_ne _ _ _ (mulMap_ne hj hjr hjn)), + hX₁ j hj] + -- Set the arguments of the copy. + let X₂ := (Y₁'.set pa t).set pb dd + have hrun₃ : ExecWithin (B := .standard) w (ram[ pa ← T; pb ← d ]) Y₁' X₂ 2 := by + ram_eval [X₂, hY₁'r T (by simp), ht, hY₁'r d hd, hdd, hda] + have hX₂ (j : Fin 21) (hj : j.val < 9) : X₂.regs j = X.regs j := by + have h9 : j ≠ pa ∧ j ≠ pb := by + refine ⟨?_, ?_⟩ <;> (intro h; subst h; simp at hj) + simp [X₂, h9.1, h9.2, hY₁'r j hj] + -- Copy. + have hXr : X₂.regs rr = r := by rw [hX₂ rr (by simp), hr] + have hsrc : (X₂.comap copyMap).regs ArrayCopy.src = t := by + change X₂.regs pa = t; simp +decide [X₂] + have hdst : (X₂.comap copyMap).regs ArrayCopy.dst = dd := by + change X₂.regs pb = dd; simp +decide [X₂] + have hlen : (X₂.comap copyMap).regs ArrayCopy.len = r := hXr + obtain ⟨Y₂, hrun₄, hY₂m, -, -, hY₂len⟩ := ArrayCopy.exec_of_separated_or_eq (w := w) + (X₂.comap copyMap) h1 (by rw [hsrc, hlen]; exact hft) (by rw [hdst, hlen]; exact hfd) + (by rw [hsrc, hdst, hlen]; omega) + have hrun₄' := hrun₄.rename copyMap_injective + rw [hlen, hsrc, hdst] at hY₂m + rw [hlen] at hrun₄' + refine ⟨Y₂.extend copyMap X₂, ?_, ?_, ?_, ?_⟩ + · exact (hrun₁.seq (hrun₂'.seq (hrun₃.seq hrun₄'))).mono (by unfold mulTime; omega) + · intro k hk + change Y₂.mem (dd + k) = _ + rw [hY₂m, ArrayCopy.copyPrefix_read hk] + exact hprod k hk + · intro a' hat had + change Y₂.mem a' = _ + rw [hY₂m, ArrayCopy.copyPrefix_outside had] + exact hframe a' hat + · intro j hj + by_cases hjr : j = rr + · subst hjr + change Function.extend copyMap Y₂.regs X₂.regs (copyMap ArrayCopy.len) = X.regs rr + rw [copyMap_injective.extend_apply, hY₂len, hlen, hr] + rw [hrun₄'.regs_eq (Code.preserves_rename_of_forall_ne _ _ _ (copyMap_ne hj hjr)), hX₂ j hj] + +/-- The hypotheses of `exec` on the three blocks, the modulus and the word size. -/ +private structure Hyp (w a b t r n e₀ : ℕ) : Prop where + hn0 : 0 < n + hab : a + r ≤ b ∨ b + r ≤ a + hta : t + r ≤ a ∨ a + r ≤ t + htb : t + r ≤ b ∨ b + r ≤ t + hfa : Fits w (a + r) + hfb : Fits w (b + r) + hft : Fits w (t + r) + hf2 : Fits w (2 * r) + hfn : Fits w (n * n + n) + hfe : Fits w e₀ + +/-- The loop invariant: the registers, the reduced blocks and the class of `R · B^e`. -/ +private def Inv (a b t r n e₀ : ℕ) (target : AdjoinRoot (X ^ r - 1 : (ZMod n)[X])) + (y : Data 21) : Prop := + y.regs R = a ∧ y.regs B = b ∧ y.regs T = t ∧ y.regs rr = r ∧ y.regs nn = n ∧ y.regs two = 2 ∧ + y.regs z = 0 ∧ y.regs e ≤ e₀ ∧ (∀ i < r, y.mem (a + i) < n) ∧ (∀ i < r, y.mem (b + i) < n) ∧ + toQuot n r (fun i => y.mem (a + i)) * toQuot n r (fun i => y.mem (b + i)) ^ y.regs e = target + +private theorem body_step {w a b t r n e₀ : ℕ} (H : Hyp w a b t r n e₀) + {target : AdjoinRoot (X ^ r - 1 : (ZMod n)[X])} (y : Data 21) + (hy : Inv a b t r n e₀ target y) (hc : (Cond.lt z e).Holds y) : + ∃ y', ExecWithin w body y y' (2 * mulTime r + 4) ∧ Inv a b t r n e₀ target y' ∧ + Nat.size (y'.regs e) < Nat.size (y.regs e) := by + obtain ⟨hyR, hyB, hyT, hyr, hyn, hy2, hyz, hye, hRred, hBred, hq⟩ := hy + obtain ⟨E, hE⟩ : ∃ E, y.regs e = E := ⟨_, rfl⟩ + have hEpos : 0 < E := by simpa only [Cond.holds_lt, hyz, hE] using hc + rw [hE] at hye hq + have h2 : Fits w 2 := H.hfn.of_le (by have := H.hn0; nlinarith) + have hfE : Fits w E := H.hfe.of_le hye + -- The parity. + let y₁ := y.set par (E % 2) + have hrun₁ : ExecWithin (B := .standard) w (ram[ par ← e % two ]) y y₁ 1 := by + ram_eval [y₁, hE, hy2, h2.of_le (Nat.mod_lt E (by norm_num)).le] + have hy₁ (j : Fin 21) (hj : j ≠ par) : y₁.regs j = y.regs j := by simp [y₁, hj] + have hz₁ : y₁.regs z = 0 := by rw [hy₁ z (by decide), hyz] + have hp₁ : y₁.regs par = E % 2 := by simp [y₁] + -- Multiply `R` by `B` when the exponent is odd. + obtain ⟨y₂, hrun₂, hy₂r, hy₂R, hy₂Rred, hy₂B⟩ : ∃ y₂, ExecWithin w + (.branch (.lt z par) (mulInto R B R) (ram[ par ← par ])) y₁ y₂ (mulTime r + 2) ∧ + (∀ j : Fin 21, j.val < 9 → y₂.regs j = y₁.regs j) ∧ + toQuot n r (fun i => y₂.mem (a + i)) = + toQuot n r (fun i => y.mem (a + i)) * toQuot n r (fun i => y.mem (b + i)) ^ (E % 2) ∧ + (∀ i < r, y₂.mem (a + i) < n) ∧ (∀ i < r, y₂.mem (b + i) = y.mem (b + i)) := by + rcases Nat.mod_two_eq_zero_or_one E with hpar | hpar + · refine ⟨y₁.set par (y₁.regs par), ?_, ?_, ?_, ?_, ?_⟩ + · refine (ExecWithin.branch_neg ?_ ExecWithin.copy).mono (by unfold mulTime; omega) + simp only [Cond.holds_lt, hz₁, hp₁, hpar]; omega + · intro j _; simp + · simp [y₁, hpar] + · simpa [y₁] using hRred + · intro i _; rfl + · have hc' : (Cond.lt z par).Holds y₁ := by + simp only [Cond.holds_lt, hz₁, hp₁, hpar]; omega + obtain ⟨Y, hrun, hprod, hframe, hregs⟩ := mulInto_exec (w := w) y₁ (x := R) (y := B) (d := R) + (by decide) (by decide) (by decide) (a := a) (b := b) (t := t) (dd := a) (r := r) (n := n) + (by rw [hy₁ R (by decide), hyR]) (by rw [hy₁ B (by decide), hyB]) + (by rw [hy₁ T (by decide), hyT]) (by rw [hy₁ R (by decide), hyR]) + (by rw [hy₁ rr (by decide), hyr]) (by rw [hy₁ nn (by decide), hyn]) H.hn0 H.hta H.htb H.hta + (by simpa [y₁] using hRred) (by simpa [y₁] using hBred) H.hfa H.hfb H.hft H.hfa H.hf2 H.hfn + refine ⟨Y, (ExecWithin.branch_pos hc' hrun).mono (by omega), hregs, ?_, ?_, ?_⟩ + · rw [toQuot_congr hprod, toQuot_mulCoeff, hpar, pow_one] + rfl + · intro i hi + rw [hprod i hi] + exact mulCoeff_lt H.hn0 _ _ _ + · intro i hi + refine hframe _ ?_ ?_ + · rcases H.htb with h | h <;> omega + · rcases H.hab with h | h <;> omega + -- Square `B`. + have hy₂' (j : Fin 21) (hj : j.val < 9) (hjp : j ≠ par) : y₂.regs j = y.regs j := by + rw [hy₂r j hj, hy₁ j hjp] + obtain ⟨y₃, hrun₃, hprod₃, hframe₃, hregs₃⟩ := mulInto_exec (w := w) y₂ (x := B) (y := B) + (d := B) (by decide) (by decide) (by decide) (a := b) (b := b) (t := t) (dd := b) (r := r) + (n := n) (by rw [hy₂' B (by decide) (by decide), hyB]) + (by rw [hy₂' B (by decide) (by decide), hyB]) (by rw [hy₂' T (by decide) (by decide), hyT]) + (by rw [hy₂' B (by decide) (by decide), hyB]) (by rw [hy₂' rr (by decide) (by decide), hyr]) + (by rw [hy₂' nn (by decide) (by decide), hyn]) H.hn0 H.htb H.htb H.htb + (fun i hi => by rw [hy₂B i hi]; exact hBred i hi) + (fun i hi => by rw [hy₂B i hi]; exact hBred i hi) + H.hfb H.hfb H.hft H.hfb H.hf2 H.hfn + have hy₃R (i : ℕ) (hi : i < r) : y₃.mem (a + i) = y₂.mem (a + i) := by + refine hframe₃ _ ?_ ?_ + · rcases H.hta with h | h <;> omega + · rcases H.hab with h | h <;> omega + -- Halve the exponent. + have hy₃e : y₃.regs e = E := by rw [hregs₃ e (by decide), hy₂' e (by decide) (by decide), hE] + have hy₃2 : y₃.regs two = 2 := by + rw [hregs₃ two (by decide), hy₂' two (by decide) (by decide), hy2] + let y₄ := y₃.set e (E / 2) + have hrun₄ : ExecWithin (B := .standard) w (ram[ e ← e / two ]) y₃ y₄ 1 := by + ram_eval [y₄, hy₃e, hy₃2, hfE.of_le (Nat.div_le_self E 2)] + have hy₄ (j : Fin 21) (hj : j.val < 9) (hjp : j ≠ par) (hje : j ≠ e) : y₄.regs j = y.regs j := by + simp only [y₄, Data.set_regs, Function.update_of_ne hje] + rw [hregs₃ j hj, hy₂' j hj hjp] + refine ⟨y₄, ?_, ?_, ?_⟩ + · exact (hrun₁.seq (hrun₂.seq (hrun₃.seq hrun₄))).mono (by omega) + · refine ⟨?_, ?_, ?_, ?_, ?_, ?_, ?_, ?_, ?_, ?_, ?_⟩ + · rw [hy₄ R (by decide) (by decide) (by decide), hyR] + · rw [hy₄ B (by decide) (by decide) (by decide), hyB] + · rw [hy₄ T (by decide) (by decide) (by decide), hyT] + · rw [hy₄ rr (by decide) (by decide) (by decide), hyr] + · rw [hy₄ nn (by decide) (by decide) (by decide), hyn] + · rw [hy₄ two (by decide) (by decide) (by decide), hy2] + · rw [hy₄ z (by decide) (by decide) (by decide), hyz] + · simp only [y₄, Data.set_regs, Function.update_self] + exact (Nat.div_le_self E 2).trans hye + · intro i hi + change y₃.mem (a + i) < n + rw [hy₃R i hi] + exact hy₂Rred i hi + · intro i hi + change y₃.mem (b + i) < n + rw [hprod₃ i hi] + exact mulCoeff_lt H.hn0 _ _ _ + · simp only [y₄, Data.set_regs, Function.update_self] + change toQuot n r (fun i => y₃.mem (a + i)) * toQuot n r (fun i => y₃.mem (b + i)) ^ (E / 2) = + target + rw [toQuot_congr hy₃R, toQuot_congr hprod₃, toQuot_mulCoeff, toQuot_congr hy₂B, hy₂R, ← hq, + ← sq, ← pow_mul, mul_assoc, ← pow_add, Nat.mod_add_div] + · simp only [y₄, Data.set_regs, Function.update_self, hE] + exact size_div_two_lt hEpos + +/-- **Powers modulo `X^r - 1` and `n`.** If the result, base and scratch blocks `[a, a + r)`, +`[b, b + r)` and `[t, t + r)` are pairwise disjoint, the result and base coefficients are below +`n ≥ 1`, and the word size holds the block ends, `2 r`, `n² + n` and the exponent `e`, then within +`time r e` instructions the result block receives reduced coefficients whose class is that of +`R · B^e`. -/ +theorem exec {w : ℕ} (x : Data 21) {a b t r n : ℕ} (ha : x.regs R = a) (hb : x.regs B = b) + (ht : x.regs T = t) (hr : x.regs rr = r) (hn : x.regs nn = n) (hn0 : 0 < n) + (hab : a + r ≤ b ∨ b + r ≤ a) (hta : t + r ≤ a ∨ a + r ≤ t) (htb : t + r ≤ b ∨ b + r ≤ t) + (hpa : ∀ i < r, x.mem (a + i) < n) (hpb : ∀ i < r, x.mem (b + i) < n) + (hfa : Fits w (a + r)) (hfb : Fits w (b + r)) (hft : Fits w (t + r)) (hf2 : Fits w (2 * r)) + (hfn : Fits w (n * n + n)) (hfe : Fits w (x.regs e)) : + ∃ y, ExecWithin w program x y (time r (x.regs e)) ∧ + toQuot n r (fun i => y.mem (a + i)) = + toQuot n r (fun i => x.mem (a + i)) * toQuot n r (fun i => x.mem (b + i)) ^ x.regs e ∧ + (∀ i < r, y.mem (a + i) < n) ∧ y.regs R = a ∧ y.regs B = b ∧ y.regs T = t ∧ + y.regs rr = r ∧ y.regs nn = n := by + have H : Hyp w a b t r n (x.regs e) := ⟨hn0, hab, hta, htb, hfa, hfb, hft, hf2, hfn, hfe⟩ + have h2 : Fits w 2 := hfn.of_le (by nlinarith) + let x₀ := (x.set two 2).set z 0 + have hrun₀ : ExecWithin (B := .standard) w (ram[ two ← #2; z ← #0 ]) x x₀ 2 := by + ram_eval [x₀, h2] + have hx₀ (j : Fin 21) (h2 : j ≠ two) (hz : j ≠ z) : x₀.regs j = x.regs j := by + simp [x₀, h2, hz] + have hI : Inv a b t r n (x.regs e) (toQuot n r (fun i => x.mem (a + i)) * + toQuot n r (fun i => x.mem (b + i)) ^ x.regs e) x₀ := by + refine ⟨?_, ?_, ?_, ?_, ?_, by simp +decide [x₀], by simp [x₀], ?_, hpa, hpb, ?_⟩ + · rw [hx₀ R (by decide) (by decide), ha] + · rw [hx₀ B (by decide) (by decide), hb] + · rw [hx₀ T (by decide) (by decide), ht] + · rw [hx₀ rr (by decide) (by decide), hr] + · rw [hx₀ nn (by decide) (by decide), hn] + · rw [hx₀ e (by decide) (by decide)] + · rw [hx₀ e (by decide) (by decide)]; rfl + obtain ⟨y, hrun, hy, hc⟩ := ExecWithin.loop_invariant _ (fun y => Nat.size (y.regs e)) + (body_step H) hI + obtain ⟨hyR, hyB, hyT, hyr, hyn, -, hyz, -, hRred, -, hq⟩ := hy + have hy0 : y.regs e = 0 := by simp only [Cond.holds_lt, hyz] at hc; omega + refine ⟨y, ?_, ?_, hRred, hyR, hyB, hyT, hyr, hyn⟩ + · refine (hrun₀.seq hrun).mono ?_ + rw [hx₀ e (by decide) (by decide), show 2 * mulTime r + 4 + 2 = 2 * mulTime r + 6 by ring] + unfold time + omega + · rw [← hq, hy0, pow_zero, mul_one] + +end Cslib.Algorithms.WordRAM.AKS.PolyPow diff --git a/Complexity/Foundations/NumberTheory/AKS/Coefficients.lean b/Complexity/Foundations/NumberTheory/AKS/Coefficients.lean index 0d4565e023..daed87fd98 100644 --- a/Complexity/Foundations/NumberTheory/AKS/Coefficients.lean +++ b/Complexity/Foundations/NumberTheory/AKS/Coefficients.lean @@ -92,6 +92,11 @@ theorem toQuot_eq_sum (n r : ℕ) (c : ℕ → ℕ) : simp only [toQuot, coeffPoly, map_sum, map_mul, map_pow, AdjoinRoot.mk_X, AdjoinRoot.mk_C] rfl +/-- The class of `∑_{i < r} c i X^i` depends only on the coefficients below `r`. -/ +theorem toQuot_congr {c d : ℕ → ℕ} (h : ∀ i < r, c i = d i) : toQuot n r c = toQuot n r d := by + rw [toQuot_eq_sum, toQuot_eq_sum] + exact sum_congr rfl fun i hi => by rw [h i (mem_range.mp hi)] + /-- **Coefficients of products.** -/ theorem toQuot_mulCoeff (n r : ℕ) (p q : ℕ → ℕ) : toQuot n r (mulCoeff n r p q) = toQuot n r p * toQuot n r q := by diff --git a/ComplexityTests/WordRAMAKSCongruence.lean b/ComplexityTests/WordRAMAKSCongruence.lean new file mode 100644 index 0000000000..6d7e61b348 --- /dev/null +++ b/ComplexityTests/WordRAMAKSCongruence.lean @@ -0,0 +1,71 @@ +/- +Copyright (c) 2026 Samuel Schlesinger. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Samuel Schlesinger +-/ + +import Complexity.Algorithms.WordRAM.AKS.Congruence + +/-! +# Powers and the AKS congruence on the word RAM + +Concrete runs of repeated squaring modulo `X^r - 1` and `n` and of the congruence test, and +instances of their correctness theorems, including a congruence refuted by running the program. +-/ + +namespace ComplexityTests.WordRAMAKSCongruence + +open Cslib.WordRAM Cslib.Algorithms.WordRAM + +/-! ## Repeated squaring -/ + +/-- `R = 1` at `0`, `B = X + 1` at `2`, scratch at `4`, with `r = 2`, `n = 3` and exponent `e`. -/ +private def powInput (e : ℕ) : Data 21 := + ⟨fun j => if j = AKS.PolyPow.B then 2 else if j = AKS.PolyPow.T then 4 else + if j = AKS.PolyPow.rr then 2 else if j = AKS.PolyPow.nn then 3 else + if j = AKS.PolyPow.e then e else 0, + fun a => if a = 0 ∨ a = 2 ∨ a = 3 then 1 else if a < 4 then 0 else 5⟩ + +/-- `(X + 1)^2 = X² + 2X + 1 ≡ 2X + 2` modulo `X² - 1` and `3`. -/ +example : (Code.eval 16 4 AKS.PolyPow.program (powInput 2)).map + (fun yn => [yn.1.mem 0, yn.1.mem 1]) = some [2, 2] := by + decide +kernel + +/-- The correctness theorem at the same input: the result block is reduced modulo `7`. -/ +example : ∃ y, ExecWithin 16 AKS.PolyPow.program (powInput 5) y (AKS.PolyPow.time 2 5) ∧ + ∀ i < 2, y.mem i < 3 := by + obtain ⟨y, hy, -, hred, -⟩ := AKS.PolyPow.exec (w := 16) (powInput 5) (a := 0) (b := 2) (t := 4) + (r := 2) (n := 3) rfl rfl rfl rfl rfl (by decide) (by decide) (by decide) (by decide) + (by decide) (by decide) (by decide) (by decide) (by decide) (by decide) (by decide) + (by decide) + exact ⟨y, hy, fun i hi => by simpa using hred i hi⟩ + +/-! ## The congruence test -/ + +/-- Blocks at `0`, `r` and `2r`, with `r`, `n` and `a` in their registers, over dirty memory. -/ +private def congInput (r n a : ℕ) : Data 25 := + ⟨fun j => if j = AKS.Congruence.B then r else if j = AKS.Congruence.T then 2 * r else + if j = AKS.Congruence.rr then r else if j = AKS.Congruence.nn then n else + if j = AKS.Congruence.av then a else 0, + fun _ => 3⟩ + +/-- `(X + 1)^3 ≡ X^3 + 1` modulo `X - 1` and `3`: both sides evaluate to `2` at `X = 1`. -/ +example : (Code.eval 16 4 AKS.Congruence.program (congInput 1 3 1)).map + (·.1.regs AKS.Congruence.flag) = some 1 := by + decide +kernel + +/-- The correctness theorem and a run refute `(X + 1)^4 ≡ X^4 + 1 (mod X - 1, 4)`: at `X = 1` +the sides are `16 ≡ 0` and `2`. -/ +example : ¬ Cslib.AKS.Congruence 4 1 1 := by + obtain ⟨y, ⟨m, -, hy⟩, -, hiff, -⟩ := AKS.Congruence.exec (w := 16) (congInput 1 4 1) + (a := 0) (b := 1) (t := 2) (r := 1) (n := 4) (A := 1) rfl rfl rfl rfl rfl rfl (by decide) + (by decide) (by decide) (by decide) (by decide) (by decide) (by decide) (by decide) + (by decide) (by decide) (by decide) + have hrun : (Code.eval 16 4 AKS.Congruence.program (congInput 1 4 1)).map + (·.1.regs AKS.Congruence.flag) = some 0 := by decide +kernel + obtain ⟨⟨y', k⟩, heval, hflag⟩ := Option.map_eq_some_iff.mp hrun + obtain ⟨rfl, -⟩ := Exec.unique hy (Exec.of_eval heval) + rw [← hiff, hflag] + decide + +end ComplexityTests.WordRAMAKSCongruence