From 6094ec9cc975f8e361f8ea49146ac4bea36aa115 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 05:13:43 -0400 Subject: [PATCH 01/11] feat(exactTriangle): relate lopsided queries and matrix products Adapt the finite counting interface and signed digit expansion from anthropics/formal-math@16a5ce41304fcdaea4e0be7c132e1c361dc26da5 for milestone #1438. This is a mathematical interface; native preprocessing, fast solvers and full composition remain under #1394. --- Complexity.lean | 2 + .../Algorithms/ExactTriangle/Lopsided.lean | 265 ++++++++++++++++++ .../Foundations/Data/Matrix/IntegerBits.lean | 112 ++++++++ ComplexityTests/Lopsided.lean | 163 +++++++++++ NOTICE | 5 +- ROADMAP.md | 8 +- UPSTREAM.md | 11 + 7 files changed, 561 insertions(+), 5 deletions(-) create mode 100644 Complexity/Algorithms/ExactTriangle/Lopsided.lean create mode 100644 Complexity/Foundations/Data/Matrix/IntegerBits.lean create mode 100644 ComplexityTests/Lopsided.lean diff --git a/Complexity.lean b/Complexity.lean index bcc3d2303d..d661ead1c1 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -1,6 +1,7 @@ module -- shake: keep-all --deprecated_module: ignore public import Complexity.Algorithms.ExactTriangle.Defs +public import Complexity.Algorithms.ExactTriangle.Lopsided public import Complexity.Algorithms.Lean.FusionTree.Capacity public import Complexity.Algorithms.Lean.FusionTree.Compression public import Complexity.Algorithms.Lean.FusionTree.Correction @@ -2084,6 +2085,7 @@ public import Complexity.Foundations.Data.List.OfFn public import Complexity.Foundations.Data.List.Range public import Complexity.Foundations.Data.List.SelectDistinct public import Complexity.Foundations.Data.Matrix.Binary +public import Complexity.Foundations.Data.Matrix.IntegerBits public import Complexity.Foundations.Data.Matrix.Permanent public import Complexity.Foundations.Data.Matrix.Reindex public import Complexity.Foundations.Data.Matrix.TotallyRegular diff --git a/Complexity/Algorithms/ExactTriangle/Lopsided.lean b/Complexity/Algorithms/ExactTriangle/Lopsided.lean new file mode 100644 index 0000000000..d52e779437 --- /dev/null +++ b/Complexity/Algorithms/ExactTriangle/Lopsided.lean @@ -0,0 +1,265 @@ +/- +Copyright (c) 2026 Anthropic, PBC, Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Anthropic, PBC, Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Foundations.Data.Matrix.IntegerBits +public import Mathlib.Data.Fintype.EquivFin +import Mathlib.Tactic.Ring + +/-! +# Finite lopsided triangle queries as matrix products + +Two outer parts have `n` vertices each; the finite middle part embeds into `Fin D`. +The matrices count common neighbors at every outer pair. Only pairs in `W` are graph edges, +so only those counts are triangle counts. Answers outside `W` are unconstrained. + +This mathematical interface gives no native program, preprocessing cost, sparsity bound, +fast solver or subquadratic running time. + +## Main definitions + +* `LopsidedGraph`: two adjacency relations, a finite middle part and queried outer edges. +* `LopsidedGraph.matX`, `LopsidedGraph.matY`: padded biadjacency matrices. +* `LopsidedGraph.signedBitInstance`: four sign choices for each pair of binary digits. + +## Main results + +* `LopsidedGraph.matX_mul_matY_apply`: entries count common middle neighbors. +* `LopsidedGraph.matX_mul_matY_apply_ne_zero_iff`: nonzero entries detect common neighbors. +* `LopsidedGraph.zero_one`: the converse for zero/one matrices. +* `LopsidedGraph.entries_from_counts`: bounded signed products from `4 * β²` counts. + +## References + +Adapted from [Anthropic's MatrixLanguage.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/MatrixLanguage.lean) +and the lopsided definitions in [PaperStatements.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/PaperStatements.lean#L798-L852). +The adaptation uses `Set.EqOn`, pointwise signed parts, and direct finite cardinalities. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.ExactTriangle + +/-- An unweighted tripartite graph whose outer parts each have `n` vertices. -/ +structure LopsidedGraph (n : ℕ) where + /-- The middle vertices. -/ + M : Type + /-- The middle part is finite. -/ + [fintypeM : Fintype M] + /-- Edges from the first outer part to the middle. -/ + adjA : Fin n → M → Prop + /-- Edges from the middle to the second outer part. -/ + adjB : M → Fin n → Prop + /-- Queried outer edges. -/ + W : Finset (Fin n × Fin n) + +attribute [instance] LopsidedGraph.fintypeM + +namespace LopsidedGraph + +/-- On a pair in `W`, a common middle neighbor closes a triangle. +For a pair outside `W`, this predicate does not assert an outer edge. -/ +def HasCommonNeighbor {n : ℕ} (I : LopsidedGraph n) (a b : Fin n) : Prop := + ∃ v : I.M, I.adjA a v ∧ I.adjB v b + +open Classical in +/-- Number of common middle neighbors. At pairs in `W` this is the triangle count. -/ +noncomputable def commonNeighborCount {n : ℕ} (I : LopsidedGraph n) (a b : Fin n) : ℕ := + (Finset.univ.filter fun v : I.M => I.adjA a v ∧ I.adjB v b).card + +/-- Two thin matrices describe adjacency where their entries equal one. +The counting theorem `zero_one` requires all entries to be zero or one. -/ +def ofMatrices {n D : ℕ} (X : Matrix (Fin n) (Fin D) ℤ) + (Y : Matrix (Fin D) (Fin n) ℤ) (W : Finset (Fin n × Fin n)) : LopsidedGraph n where + M := Fin D + adjA a v := X a v = 1 + adjB v b := Y v b = 1 + W := W + +open Classical in +/-- The count for a matrix-defined graph uses the ordinary finite middle indices. -/ +theorem commonNeighborCount_ofMatrices {n D : ℕ} (X : Matrix (Fin n) (Fin D) ℤ) + (Y : Matrix (Fin D) (Fin n) ℤ) (W : Finset (Fin n × Fin n)) (a b : Fin n) : + (ofMatrices X Y W).commonNeighborCount a b = + (Finset.univ.filter fun v : Fin D => X a v = 1 ∧ Y v b = 1).card := by + unfold ofMatrices commonNeighborCount + congr 1 + ext v + simp + +open Classical in +/-- The first biadjacency matrix, with zero columns outside the injection range. -/ +noncomputable def matX {n D : ℕ} (I : LopsidedGraph n) (e : I.M ↪ Fin D) : + Matrix (Fin n) (Fin D) ℤ := + Matrix.of fun a j => if ∃ v : I.M, e v = j ∧ I.adjA a v then 1 else 0 + +open Classical in +/-- The second biadjacency matrix, with zero rows outside the injection range. -/ +noncomputable def matY {n D : ℕ} (I : LopsidedGraph n) (e : I.M ↪ Fin D) : + Matrix (Fin D) (Fin n) ℤ := + Matrix.of fun j b => if ∃ v : I.M, e v = j ∧ I.adjB v b then 1 else 0 + +variable {n D : ℕ} (I : LopsidedGraph n) (e : I.M ↪ Fin D) + +/-- Every first biadjacency entry is zero or one. -/ +theorem matX_eq_zero_or_one (a : Fin n) (j : Fin D) : + I.matX e a j = 0 ∨ I.matX e a j = 1 := by + classical + simp only [matX, Matrix.of_apply] + split_ifs <;> simp + +/-- Every second biadjacency entry is zero or one. -/ +theorem matY_eq_zero_or_one (j : Fin D) (b : Fin n) : + I.matY e j b = 0 ∨ I.matY e j b = 1 := by + classical + simp only [matY, Matrix.of_apply] + split_ifs <;> simp + +/-- A middle part of cardinality at most `D` admits the required numbering. -/ +theorem nonempty_embedding (h : Fintype.card I.M ≤ D) : Nonempty (I.M ↪ Fin D) := + Function.Embedding.nonempty_of_card_le (by simpa using h) + +open Classical in +/-- At a middle vertex, the first matrix records its adjacency. -/ +theorem matX_apply_emb (a : Fin n) (v : I.M) : + I.matX e a (e v) = if I.adjA a v then 1 else 0 := by + simp [matX] + +open Classical in +/-- At a middle vertex, the second matrix records its adjacency. -/ +theorem matY_apply_emb (v : I.M) (b : Fin n) : + I.matY e (e v) b = if I.adjB v b then 1 else 0 := by + simp [matY] + +/-- An unused column is zero. -/ +theorem matX_apply_of_notMem_range (a : Fin n) {j : Fin D} (hj : j ∉ Set.range e) : + I.matX e a j = 0 := + ite_eq_right fun ⟨v, hv, _⟩ => hj ⟨v, hv⟩ + +/-- An unused row is zero. -/ +theorem matY_apply_of_notMem_range {j : Fin D} (hj : j ∉ Set.range e) (b : Fin n) : + I.matY e j b = 0 := + ite_eq_right fun ⟨v, hv, _⟩ => hj ⟨v, hv⟩ + +/-- A pair has a common neighbor if and only if its count is nonzero. -/ +theorem commonNeighborCount_ne_zero_iff (a b : Fin n) : + I.commonNeighborCount a b ≠ 0 ↔ I.HasCommonNeighbor a b := by + classical + rw [commonNeighborCount, Finset.card_ne_zero, Finset.filter_nonempty_iff] + simp [HasCommonNeighbor] + +/-- Every product entry counts common middle neighbors, independently of the numbering. -/ +theorem matX_mul_matY_apply (a b : Fin n) : + (I.matX e * I.matY e) a b = (I.commonNeighborCount a b : ℤ) := by + classical + rw [Matrix.mul_apply, commonNeighborCount, Finset.card_filter, Nat.cast_sum] + -- Only the columns of the middle vertices matter. + rw [← Finset.sum_subset (Finset.subset_univ (Finset.univ.map e)) fun j _ hj => by + rw [I.matX_apply_of_notMem_range e a (by simpa using hj), zero_mul], Finset.sum_map] + refine Finset.sum_congr rfl fun v _ => ?_ + rw [I.matX_apply_emb, I.matY_apply_emb] + by_cases hA : I.adjA a v <;> by_cases hB : I.adjB v b <;> simp [hA, hB] + +/-- A nonzero product entry detects a common middle neighbor. -/ +theorem matX_mul_matY_apply_ne_zero_iff (a b : Fin n) : + (I.matX e * I.matY e) a b ≠ 0 ↔ I.HasCommonNeighbor a b := by + rw [matX_mul_matY_apply, Ne, Nat.cast_eq_zero] + exact I.commonNeighborCount_ne_zero_iff a b + +/-- On queried outer edges, correct natural counts are exactly the integer product entries. +Ordinary `Set.EqOn` leaves all unqueried entries unconstrained. -/ +theorem counting_eqOn_iff (ans : Fin n × Fin n → ℕ) : + Set.EqOn ans (fun q => I.commonNeighborCount q.1 q.2) (I.W : Set _) ↔ + Set.EqOn (fun q => (ans q : ℤ)) (fun q => (I.matX e * I.matY e) q.1 q.2) + (I.W : Set _) := by + simp only [Set.EqOn, matX_mul_matY_apply, Nat.cast_inj] + +/-- Counts on queried edges provide their detection answers. -/ +theorem detection_eqOn_of_counting (ans : Fin n × Fin n → ℕ) + (h : Set.EqOn ans (fun q => I.commonNeighborCount q.1 q.2) (I.W : Set _)) : + Set.EqOn (fun q => ans q ≠ 0) (fun q => I.HasCommonNeighbor q.1 q.2) (I.W : Set _) := by + intro q hq + change (ans q ≠ 0) = I.HasCommonNeighbor q.1 q.2 + apply propext + rw [h hq, commonNeighborCount_ne_zero_iff] + +/-- For zero/one matrices, each product entry is the count in the graph they describe. -/ +theorem zero_one {n D : ℕ} (X : Matrix (Fin n) (Fin D) ℤ) (Y : Matrix (Fin D) (Fin n) ℤ) + (W : Finset (Fin n × Fin n)) (hX : ∀ a j, X a j = 0 ∨ X a j = 1) + (hY : ∀ j b, Y j b = 0 ∨ Y j b = 1) (a b : Fin n) : + (X * Y) a b = ((ofMatrices X Y W).commonNeighborCount a b : ℤ) := by + classical + rw [Matrix.mul_apply, commonNeighborCount, Finset.card_filter, Nat.cast_sum] + refine Finset.sum_congr rfl fun j _ => ?_ + rcases hX a j with hx | hx <;> rcases hY j b with hy | hy <;> + simp [ofMatrices, hx, hy] + +/-- One of four sign choices for a pair of binary digits; `true` selects the negative part. +Every instance has exactly `D` middle vertices and the original queried edges `W`. -/ +def signedBitInstance (X : Matrix (Fin n) (Fin D) ℤ) (Y : Matrix (Fin D) (Fin n) ℤ) + (W : Finset (Fin n × Fin n)) (negativeX negativeY : Bool) (i j : ℕ) : LopsidedGraph n := + ofMatrices (IntMatrix.bit (if negativeX then (X.map negPart) else (X.map posPart)) i) + (IntMatrix.bit (if negativeY then (Y.map negPart) else (Y.map posPart)) j) W + +/-- The finite family of sign/digit instances for a magnitude width `β`. +Different slots may describe the same graph. -/ +def signedBitInstances (X : Matrix (Fin n) (Fin D) ℤ) (Y : Matrix (Fin D) (Fin n) ℤ) + (W : Finset (Fin n × Fin n)) (β : ℕ) : + (Bool × Bool × Fin β × Fin β) → LopsidedGraph n := + fun q => signedBitInstance X Y W q.1 q.2.1 q.2.2.1.val q.2.2.2.val + +/-- The indexed family has `4 * β²` slots, allowing repeated or zero instances. -/ +theorem card_signedBit_indices (β : ℕ) : + Fintype.card (Bool × Bool × Fin β × Fin β) = 4 * β ^ 2 := by + simp [pow_two, Nat.mul_assoc, Nat.mul_left_comm, Nat.mul_comm] + +/-- Each nonnegative product entry is a weighted sum of digit-graph counts. -/ +private theorem mul_apply_eq_sum_counts {n D : ℕ} (X : Matrix (Fin n) (Fin D) ℤ) + (Y : Matrix (Fin D) (Fin n) ℤ) (W : Finset (Fin n × Fin n)) (β : ℕ) + (hX : ∀ a j, 0 ≤ X a j ∧ X a j < 2 ^ β) (hY : ∀ j b, 0 ≤ Y j b ∧ Y j b < 2 ^ β) (a b : Fin n) : + (X * Y) a b = ∑ i ∈ Finset.range β, ∑ j ∈ Finset.range β, (2 : ℤ) ^ (i + j) * + ((ofMatrices (IntMatrix.bit X i) (IntMatrix.bit Y j) W).commonNeighborCount a b : ℤ) := by + rw [IntMatrix.product X Y β hX hY, Matrix.sum_apply] + refine Finset.sum_congr rfl fun i _ => ?_ + rw [Matrix.sum_apply] + refine Finset.sum_congr rfl fun j _ => ?_ + rw [Matrix.smul_apply, smul_eq_mul, zero_one _ _ W (IntMatrix.bit_eq_zero_or_one X i) + (IntMatrix.bit_eq_zero_or_one Y j)] + +/-- Signed bounded product entries are weighted sums of four count instances per digit +pair, with signs `+, −, −, +`. Only counts at the same queried edge are needed. The absolute +bounds are strict; `β = 0` therefore permits only zero entries. -/ +theorem entries_from_counts {n D : ℕ} (X : Matrix (Fin n) (Fin D) ℤ) + (Y : Matrix (Fin D) (Fin n) ℤ) (W : Finset (Fin n × Fin n)) (β : ℕ) + (hX : ∀ a j, |X a j| < 2 ^ β) (hY : ∀ j b, |Y j b| < 2 ^ β) (a b : Fin n) : + (X * Y) a b = ∑ i ∈ Finset.range β, ∑ j ∈ Finset.range β, (2 : ℤ) ^ (i + j) * + (((signedBitInstance X Y W false false i j).commonNeighborCount + a b : ℤ) + - ((signedBitInstance X Y W false true i j).commonNeighborCount + a b : ℤ) + - ((signedBitInstance X Y W true false i j).commonNeighborCount + a b : ℤ) + + ((signedBitInstance X Y W true true i j).commonNeighborCount + a b : ℤ)) := by + simp only [signedBitInstance, Bool.false_eq_true, ↓reduceIte] + obtain ⟨hXpos, hXneg⟩ := IntMatrix.parts_lt X hX + obtain ⟨hYpos, hYneg⟩ := IntMatrix.parts_lt Y hY + -- `XY` is a signed sum of the four products of a part of `X` and a part of `Y`. + conv_lhs => rw [IntMatrix.parts X, IntMatrix.parts Y] + rw [Matrix.sub_mul, Matrix.mul_sub, Matrix.mul_sub] + simp only [Matrix.sub_apply] + rw [mul_apply_eq_sum_counts _ _ W β hXpos hYpos, + mul_apply_eq_sum_counts _ _ W β hXpos hYneg, + mul_apply_eq_sum_counts _ _ W β hXneg hYpos, + mul_apply_eq_sum_counts _ _ W β hXneg hYneg] + simp only [mul_sub, mul_add, Finset.sum_sub_distrib, Finset.sum_add_distrib] + ring + +end LopsidedGraph + +end Cslib.Algorithms.ExactTriangle diff --git a/Complexity/Foundations/Data/Matrix/IntegerBits.lean b/Complexity/Foundations/Data/Matrix/IntegerBits.lean new file mode 100644 index 0000000000..5083b77db3 --- /dev/null +++ b/Complexity/Foundations/Data/Matrix/IntegerBits.lean @@ -0,0 +1,112 @@ +/- +Copyright (c) 2026 Anthropic, PBC, Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Anthropic, PBC, Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Init +public import Mathlib.Algebra.Order.Group.PosPart +public import Mathlib.Data.Matrix.Basic +public import Mathlib.Data.Nat.Log +import Mathlib.Tactic.Positivity +import Mathlib.Tactic.Ring + +/-! +# Binary digits of integer matrices + +Finite binary expansions and their products, using Mathlib's pointwise positive and negative +parts. These are algebraic identities, with no algorithm or running-time assertion. + +## Main definitions + +* `IntMatrix.bit`: the binary digit matrix of a given weight. + +## Main results + +* `IntMatrix.binary`: expansion of a nonnegative bounded matrix. +* `IntMatrix.product`: expansion of a product into products of digit matrices. +* `IntMatrix.abs_lt_two_pow_log`: an explicit integer bit bound, including zero. + +## References + +Adapted from [Anthropic's MatrixLanguage.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/MatrixLanguage.lean), +`matBit` and `Footnote8.binary/product`. General index types and the integer bit bound are original +adaptation code. +-/ + +@[expose] public section + +namespace Cslib.IntMatrix + +/-- The binary digit of weight `2^i`, entry by entry. For negative entries this is the +Euclidean remainder of the quotient; reconstruction below assumes nonnegative entries. -/ +def bit {m n : Type*} (X : Matrix m n ℤ) (i : ℕ) : Matrix m n ℤ := + Matrix.of fun a j => X a j / 2 ^ i % 2 + +/-- `x mod 2^β` is the sum of the last `β` binary digits of `x`, with their weights. -/ +private theorem emod_two_pow_eq_sum_digits (x : ℤ) (β : ℕ) : + x % 2 ^ β = ∑ i ∈ Finset.range β, 2 ^ i * (x / 2 ^ i % 2) := by + induction β with + | zero => simp + | succ β ih => + have h : x / (2 ^ β * 2) = x / 2 ^ β / 2 := (Int.ediv_ediv_of_nonneg (by positivity)).symm + rw [Finset.sum_range_succ, ← ih, pow_succ, Int.emod_def, Int.emod_def x (2 ^ β), + Int.emod_def (x / 2 ^ β), h] + ring + +/-- A matrix with entries in `[0, 2^β)` is the weighted sum of its first `β` digit matrices. -/ +theorem binary {m n : Type*} (X : Matrix m n ℤ) (β : ℕ) + (hX : ∀ a j, 0 ≤ X a j ∧ X a j < 2 ^ β) : + X = ∑ i ∈ Finset.range β, (2 : ℤ) ^ i • bit X i := by + ext a j + rw [Matrix.sum_apply] + simp only [Matrix.smul_apply, bit, Matrix.of_apply, smul_eq_mul] + rw [← emod_two_pow_eq_sum_digits, Int.emod_eq_of_lt (hX a j).1 (hX a j).2] + +/-- Every digit matrix has entries zero or one, even for signed input. -/ +theorem bit_eq_zero_or_one {m n : Type*} (X : Matrix m n ℤ) (i : ℕ) (a : m) + (j : n) : bit X i a j = 0 ∨ bit X i a j = 1 := by + simp only [bit, Matrix.of_apply] + omega + +/-- Expand a product of nonnegative bounded matrices into products of binary digits. -/ +theorem product {m d n : Type*} [Fintype d] (X : Matrix m d ℤ) (Y : Matrix d n ℤ) + (β : ℕ) (hX : ∀ a j, 0 ≤ X a j ∧ X a j < 2 ^ β) (hY : ∀ j b, 0 ≤ Y j b ∧ Y j b < 2 ^ β) : + X * Y = ∑ i ∈ Finset.range β, ∑ j ∈ Finset.range β, + (2 : ℤ) ^ (i + j) • (bit X i * bit Y j) := by + conv_lhs => rw [binary X β hX, binary Y β hY] + rw [Matrix.sum_mul] + refine Finset.sum_congr rfl fun i _ => ?_ + rw [Matrix.mul_sum] + refine Finset.sum_congr rfl fun j _ => ?_ + rw [Matrix.smul_mul, Matrix.mul_smul, smul_smul, pow_add] + +/-- Entrywise positive and negative parts reconstruct an integer matrix. +`Matrix.map` reuses scalar signed parts without adding an ordering instance to matrices. -/ +theorem parts {m n : Type*} (X : Matrix m n ℤ) : + X = X.map posPart - X.map negPart := by + ext a j + exact (posPart_sub_negPart (X a j)).symm + +/-- Both pointwise signed parts of an absolutely bounded matrix satisfy the same nonnegative +bound. This uses the existing ordered-group operations, without introducing matrix variants. -/ +theorem parts_lt {m n : Type*} (X : Matrix m n ℤ) {β : ℕ} + (hX : ∀ a j, |X a j| < 2 ^ β) : + (∀ a j, 0 ≤ (X.map posPart) a j ∧ (X.map posPart) a j < 2 ^ β) ∧ + ∀ a j, 0 ≤ (X.map negPart) a j ∧ (X.map negPart) a j < 2 ^ β := by + have hpow : (0 : ℤ) < 2 ^ β := by positivity + exact ⟨fun a j => ⟨le_max_right _ _, max_lt ((le_abs_self _).trans_lt (hX a j)) hpow⟩, + fun a j => ⟨le_max_right _ _, max_lt ((neg_le_abs _).trans_lt (hX a j)) hpow⟩⟩ + +/-- An integer bounded in magnitude by `U` needs at most `log₂ U + 1` magnitude bits. +This also holds for `U = 0`; no real logarithms or asymptotic claim are involved. -/ +theorem abs_lt_two_pow_log {U : ℕ} (x : ℤ) (hx : |x| ≤ U) : + |x| < 2 ^ (Nat.log 2 U + 1) := by + have h := Nat.lt_pow_succ_log_self (by decide : 1 < 2) U + have hcast : (U : ℤ) < 2 ^ (Nat.log 2 U + 1) := by exact_mod_cast h + exact hx.trans_lt hcast + +end Cslib.IntMatrix diff --git a/ComplexityTests/Lopsided.lean b/ComplexityTests/Lopsided.lean new file mode 100644 index 0000000000..737cbf386a --- /dev/null +++ b/ComplexityTests/Lopsided.lean @@ -0,0 +1,163 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +import Complexity.Algorithms.ExactTriangle.Lopsided +import Mathlib.LinearAlgebra.Matrix.Notation + +/-! Concrete tests for padding, counts, queried answers and signed reconstruction. -/ + +open Cslib.Algorithms.ExactTriangle +open LopsidedGraph + +private abbrev twoNeighbors : LopsidedGraph 2 where + M := Fin 2 + adjA a _ := a = 0 + adjB _ b := b = 1 + W := {(0, 1)} + +private def padded : Fin 2 ↪ Fin 3 := ⟨fun v => ⟨v.val, by omega⟩, + fun _ _ h => Fin.ext (show _ from congrArg (fun v : Fin 3 => v.val) h)⟩ + +example : twoNeighbors.commonNeighborCount 0 1 = 2 := by + norm_num [commonNeighborCount, twoNeighbors, Finset.univ_fin2] + rfl + +example : twoNeighbors.HasCommonNeighbor 0 1 := ⟨0, rfl, rfl⟩ + +example : ¬twoNeighbors.HasCommonNeighbor 1 1 := by + simp [HasCommonNeighbor, twoNeighbors] + +example : twoNeighbors.matX padded 0 2 = 0 ∧ twoNeighbors.matY padded 2 1 = 0 := by + have h : (2 : Fin 3) ∉ Set.range padded := by + rintro ⟨v, hv⟩ + have hv' : v.val = 2 := congrArg (fun w : Fin 3 => w.val) hv + have := v.isLt + omega + exact ⟨matX_apply_of_notMem_range _ _ _ h, matY_apply_of_notMem_range _ _ h _⟩ + +example : (twoNeighbors.matX padded * twoNeighbors.matY padded) 0 1 = 2 := by + rw [matX_mul_matY_apply] + norm_num [commonNeighborCount, twoNeighbors, Finset.univ_fin2] + rfl + +-- Values away from the sole query are deliberately incorrect and still allowed. +example : Set.EqOn (fun q : Fin 2 × Fin 2 => if q = (0, 1) then 2 else 99) + (fun q => twoNeighbors.commonNeighborCount q.1 q.2) (twoNeighbors.W : Set _) := by + intro q hq + have h : q = (0, 1) := by simpa [twoNeighbors] using hq + subst q + norm_num [commonNeighborCount, twoNeighbors, Finset.univ_fin2] + rfl + +private abbrev emptyMiddle : LopsidedGraph 1 where + M := Fin 0 + adjA _ _ := True + adjB _ _ := True + W := {(0, 0)} + +example : emptyMiddle.commonNeighborCount 0 0 = 0 := by + simp [commonNeighborCount, emptyMiddle] + +example : ¬emptyMiddle.HasCommonNeighbor 0 0 := by + simp [HasCommonNeighbor, emptyMiddle] + +example : (emptyMiddle.matX (Function.Embedding.refl _) * + emptyMiddle.matY (Function.Embedding.refl _)) 0 0 = 0 := by + rw [matX_mul_matY_apply] + simp [commonNeighborCount, emptyMiddle] + +private def signedX : Matrix (Fin 1) (Fin 2) ℤ := !![-3, 2] +private def signedY : Matrix (Fin 2) (Fin 1) ℤ := !![2; -1] + +example : (signedX * signedY) 0 0 = -8 := by + norm_num [signedX, signedY, Matrix.mul_apply, Fin.sum_univ_two] + +example : (signedBitInstance signedX signedY {(0, 0)} true false 0 1).commonNeighborCount + 0 0 = 1 := by + classical + simp only [signedBitInstance, commonNeighborCount_ofMatrices, ↓reduceIte] + norm_num [signedX, signedY, + Cslib.IntMatrix.bit, posPart_def, negPart_def, Finset.univ_fin2, + Finset.filter_insert, Finset.filter_singleton] + +example : (signedBitInstance signedX signedY {(0, 0)} false true 1 0).commonNeighborCount + 0 0 = 1 := by + classical + simp only [signedBitInstance, commonNeighborCount_ofMatrices, ↓reduceIte] + norm_num [signedX, signedY, + Cslib.IntMatrix.bit, posPart_def, negPart_def, Finset.univ_fin2, + Finset.filter_insert, Finset.filter_singleton] + +-- Apply the signed reconstruction theorem to an actual product, including all four signs. +example : (∑ i ∈ Finset.range 2, ∑ j ∈ Finset.range 2, (2 : ℤ) ^ (i + j) * + (((signedBitInstance signedX signedY {(0, 0)} false false i j).commonNeighborCount 0 0 : ℤ) + - ((signedBitInstance signedX signedY {(0, 0)} false true i j).commonNeighborCount 0 0 : ℤ) + - ((signedBitInstance signedX signedY {(0, 0)} true false i j).commonNeighborCount 0 0 : ℤ) + + ((signedBitInstance signedX signedY {(0, 0)} true true i j).commonNeighborCount 0 0 : ℤ))) + = -8 := by + rw [← entries_from_counts signedX signedY {(0, 0)} 2 + (by intro a j; fin_cases a; fin_cases j <;> norm_num [signedX]) + (by intro j b; fin_cases j <;> fin_cases b <;> norm_num [signedY])] + norm_num [signedX, signedY, Matrix.mul_apply, Fin.sum_univ_two] + +example : Cslib.IntMatrix.bit (!![-1, 0, 3] : Matrix (Fin 1) (Fin 3) ℤ) 1 0 0 = 1 := by + norm_num [Cslib.IntMatrix.bit] + +example : |(-7 : ℤ)| < 2 ^ (Nat.log 2 7 + 1) := by + exact Cslib.IntMatrix.abs_lt_two_pow_log (-7) (by norm_num) + +example : |(0 : ℤ)| < 2 ^ (Nat.log 2 0 + 1) := by + exact Cslib.IntMatrix.abs_lt_two_pow_log 0 (by norm_num) + +example : (signedX * (-signedY)) 0 0 = 8 := by + norm_num [signedX, signedY, Matrix.mul_apply, Fin.sum_univ_two] + +example : (ofMatrices (!![1, 1] : Matrix (Fin 1) (Fin 2) ℤ) !![1; 1] + {(0, 0)}).commonNeighborCount 0 0 = 2 := by + have h := zero_one (!![1, 1] : Matrix (Fin 1) (Fin 2) ℤ) !![1; 1] {(0, 0)} + (by intro a j; fin_cases a; fin_cases j <;> simp) + (by intro j b; fin_cases j <;> fin_cases b <;> simp) 0 0 + norm_num [Matrix.mul_apply, Fin.sum_univ_two] at h + exact_mod_cast h.symm + +example : (signedBitInstance signedX (-signedY) {(0, 0)} true true 0 1).commonNeighborCount + 0 0 = 1 := by + classical + simp only [signedBitInstance, commonNeighborCount_ofMatrices, ↓reduceIte] + norm_num [signedX, signedY, + Cslib.IntMatrix.bit, posPart_def, negPart_def, Finset.univ_fin2, + Finset.filter_insert, Finset.filter_singleton] + +example : (signedBitInstances signedX signedY {(0, 0)} 2 (true, false, 0, 1)).W = + {(0, 0)} := rfl + +example : Fintype.card (Bool × Bool × Fin 2 × Fin 2) = 16 := by decide + +-- Empty outer parts with a genuinely nonempty middle part still give an empty product. +example : (0 : Matrix (Fin 0) (Fin 3) ℤ) * (0 : Matrix (Fin 3) (Fin 0) ℤ) = 0 := by + ext a b + exact Fin.elim0 a + +example : Set.EqOn (fun _ : Fin 0 × Fin 0 => (99 : ℕ)) + (fun q => + (ofMatrices (0 : Matrix (Fin 0) (Fin 3) ℤ) 0 ∅).commonNeighborCount q.1 q.2) + (∅ : Set (Fin 0 × Fin 0)) := by + simp + +-- With beta zero, strict magnitude bounds are satisfied by zero matrices of nonempty shape. +example : ((0 : Matrix (Fin 1) (Fin 2) ℤ) * (0 : Matrix (Fin 2) (Fin 1) ℤ)) 0 0 = 0 := by + have h := entries_from_counts (0 : Matrix (Fin 1) (Fin 2) ℤ) 0 {(0, 0)} 0 + (by simp) (by simp) 0 0 + simpa only [Finset.range_zero, Finset.sum_empty] using h + +-- A queried edge with disjoint middle neighborhoods is a negative triangle query. +example : let I := ofMatrices (!![1, 0] : Matrix (Fin 1) (Fin 2) ℤ) !![0; 1] {(0, 0)} + (0, 0) ∈ I.W ∧ I.commonNeighborCount 0 0 = 0 := by + dsimp only + constructor + · simp [ofMatrices] + · rw [commonNeighborCount_ofMatrices] + norm_num [Finset.univ_fin2, Finset.filter_insert, Finset.filter_singleton] diff --git a/NOTICE b/NOTICE index 81584c2411..3ab01a2df3 100644 --- a/NOTICE +++ b/NOTICE @@ -1,8 +1,9 @@ The 3SUM/Exact Triangle reductions, their finite counting lemmas and native matrix builder in -Complexity/Algorithms/ExactTriangle/Defs.lean, +Complexity/Algorithms/ExactTriangle/{Defs,Lopsided}.lean, Complexity/Algorithms/ThreeSum/, Complexity/Algorithms/WordRAM/ThreeSum/, -Complexity/Foundations/Data/Finset/Average.lean and +Complexity/Foundations/Data/Finset/Average.lean, +Complexity/Foundations/Data/Matrix/IntegerBits.lean and Complexity/Foundations/NumberTheory/PrimeCounting.lean are adapted from anthropics/formal-math, 3sum-apsp, revision 16a5ce41304fcdaea4e0be7c132e1c361dc26da5. diff --git a/ROADMAP.md b/ROADMAP.md index 8d53f5250f..5b3b230013 100644 --- a/ROADMAP.md +++ b/ROADMAP.md @@ -421,8 +421,10 @@ each statement against the numbered result. Organize modules by subject, not cha transfer of Anthropic's modular machine and 3SUM/APSP algorithms; direct 3SUM adaptation ([#1394]) has the Convolution-3SUM/Exact Triangle correctness reduction, a charged native single-instance matrix builder and the complete finite mathematical Chan–He reduction ([#1435]), including - ordinary distinct-position inputs and explicit size bounds. It still needs native Chan–He - preprocessing, the reduction host, the fast Exact Triangle solver and full input/time/width/ + ordinary distinct-position inputs and explicit size bounds. The finite lopsided matrix/counting + interface ([#1438]) adds padding, detection and signed digit reconstruction; it still needs + modular Exact Triangle hashing/false-positive splitting, tensor/sparse-query programs, native + Chan–He preprocessing, the reduction host, the fast Exact Triangle solver and full input/time/width/ asymptotic composition, then the improved Section 4 route to exponent 1.9992; size-based width selection with an input-fit guarantee and explicit overhead; expected-time rules for unbounded almost-surely terminating structured loops; random-tape equivalence and BPP/RP/ZPP simulation @@ -546,6 +548,6 @@ not planned, not as proved. [#1155]: https://github.com/cslib-community/complexitylib/issues/1155 [#1351]: https://github.com/cslib-community/complexitylib/issues/1351 [#1394]: https://github.com/cslib-community/complexitylib/issues/1394 - [#1435]: https://github.com/cslib-community/complexitylib/issues/1435 +[#1438]: https://github.com/cslib-community/complexitylib/issues/1438 [#1449]: https://github.com/cslib-community/complexitylib/issues/1449 diff --git a/UPSTREAM.md b/UPSTREAM.md index f465fed043..87d9cb5288 100644 --- a/UPSTREAM.md +++ b/UPSTREAM.md @@ -165,6 +165,17 @@ Do not add a ledger row for an ordinary candidate or repeat its source marker he length, count and magnitude bounds; `scanPairs` counts candidate-prime/element pairs and is not a charged running time. Generic finite-average and prime-counting lemmas live in Foundations. See milestone [#1435](https://github.com/cslib-community/complexitylib/issues/1435). + `Algorithms/ExactTriangle/Lopsided.lean` adapts the finite unweighted lopsided query/counting + interface from the same pin's `Sec3/MatrixLanguage.lean`: arbitrary-injection zero padding, + common-neighbor counts as product entries, detection, the zero/one converse, and signed products + from `4 * β²` count instances with the original queried edges. Its digit identities live in + `Foundations/Data/Matrix/IntegerBits.lean`, reuse Mathlib's pointwise positive/negative parts, + and give an explicit `Nat.log 2 U + 1` magnitude-bit bound including zero. Query answers use + `Set.EqOn`; no additional weighted-triangle predicate or source machine is introduced. + This is the mathematical interface milestone + [#1438](https://github.com/cslib-community/complexitylib/issues/1438), with no sparsity assumption, + preprocessing/time claim or solver. Modular hashing/false-positive splitting and tensor/sparse + query programs remain in #1394, alongside the native and asymptotic composition below. The existing signed and unsigned APIs were reused; no source machine or separate run semantics were copied. The full reduction host, native Chan–He preprocessing, fast Exact Triangle and full input/time/width/asymptotic composition remain in From dec963a5dd581ae8acecb846e7b87b1bf7e5c78e Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 06:02:25 -0400 Subject: [PATCH 02/11] feat(exactTriangle): prove modular counts and the selected-prime bound --- Complexity.lean | 4 + .../Algorithms/ExactTriangle/Hashing.lean | 482 ++++++++++++++++++ .../Algorithms/ExactTriangle/Weights.lean | 61 +++ .../Foundations/NumberTheory/PrimeWindow.lean | 346 +++++++++++++ .../Foundations/RingTheory/CyclicIntRing.lean | 84 +++ ComplexityTests/ExactTriangleHashing.lean | 115 +++++ NOTICE | 9 +- ROADMAP.md | 6 +- UPSTREAM.md | 12 +- 9 files changed, 1113 insertions(+), 6 deletions(-) create mode 100644 Complexity/Algorithms/ExactTriangle/Hashing.lean create mode 100644 Complexity/Algorithms/ExactTriangle/Weights.lean create mode 100644 Complexity/Foundations/NumberTheory/PrimeWindow.lean create mode 100644 Complexity/Foundations/RingTheory/CyclicIntRing.lean create mode 100644 ComplexityTests/ExactTriangleHashing.lean diff --git a/Complexity.lean b/Complexity.lean index d661ead1c1..e6ac945545 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -1,7 +1,9 @@ module -- shake: keep-all --deprecated_module: ignore public import Complexity.Algorithms.ExactTriangle.Defs +public import Complexity.Algorithms.ExactTriangle.Hashing public import Complexity.Algorithms.ExactTriangle.Lopsided +public import Complexity.Algorithms.ExactTriangle.Weights public import Complexity.Algorithms.Lean.FusionTree.Capacity public import Complexity.Algorithms.Lean.FusionTree.Compression public import Complexity.Algorithms.Lean.FusionTree.Correction @@ -2190,6 +2192,7 @@ public import Complexity.Foundations.NumberTheory.AKS.Introspective public import Complexity.Foundations.NumberTheory.AKS.SmallR public import Complexity.Foundations.NumberTheory.Harmonic public import Complexity.Foundations.NumberTheory.PrimeCounting +public import Complexity.Foundations.NumberTheory.PrimeWindow public import Complexity.Foundations.Order.BooleanSubalgebra public import Complexity.Foundations.Order.Changes public import Complexity.Foundations.Order.UpperLower @@ -2200,6 +2203,7 @@ public import Complexity.Foundations.Probability.SecondMoment public import Complexity.Foundations.Quantum.Measurement public import Complexity.Foundations.Relation.InductiveCounting public import Complexity.Foundations.Relation.RelatesWithinSteps +public import Complexity.Foundations.RingTheory.CyclicIntRing public import Complexity.Foundations.Semantics.Stateful.Amortized public import Complexity.Foundations.Semantics.Stateful.Basic public import Complexity.Foundations.Semantics.Stateful.Composition diff --git a/Complexity/Algorithms/ExactTriangle/Hashing.lean b/Complexity/Algorithms/ExactTriangle/Hashing.lean new file mode 100644 index 0000000000..6799e79fc7 --- /dev/null +++ b/Complexity/Algorithms/ExactTriangle/Hashing.lean @@ -0,0 +1,482 @@ +/- +Copyright (c) 2026 Anthropic, PBC, Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Anthropic, PBC, Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.ExactTriangle.Weights +public import Complexity.Foundations.NumberTheory.PrimeWindow +public import Complexity.Foundations.RingTheory.CyclicIntRing +public import Mathlib.Data.Matrix.Mul + +/-! +# Modular counts and prime selection for integer Exact Triangle + +Let `S(a,b,c)` denote `weightSum`, `F(p)` the false-positive count and `Z₀` the zero count. A +minimum-count prime in `[sqrt D / 2, sqrt D)` also minimises false positives. For weights +bounded by `n^κ`, with `16 ≤ D ≤ n` and `1 ≤ κ`, a universal constant bounds its false positives +by `C * κ * n^3 * log n / sqrt D`. The prime count used here is proved for every admitted `D`. + +## Main definitions + +* `countZeroMod`, `falsePositiveCount`, `zeroCount`: counts of finite triples. +* `IsSelectedPrime`: a minimum-count prime in the specified window. +* `matP`, `matQ`: cyclic polynomial matrices encoding edge residues. + +## Main results + +* `countZeroMod_eq`: modular count is the sum of false positives and true zeros. +* `coeff_matP_mul_matQ`, `falsePositiveCount_add_zeroCount_eq_sum_coeff`: coefficient + identities. +* `exists_isSelectedPrime`: the minimum exists. +* `Hashing.falsePositiveCount_le_const_mul`: the uniform false-positive bound. + +This is a mathematical reduction step. No native program or charged running-time bound is +asserted. The supplied three matrices retain the existing `HasZeroTriangle` semantics. + +## References + +Direct adaptation of [Anthropic's +Hashing.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/Theorem17/Hashing.lean), +with the minimal count/selection definitions extracted from its `PaperStatements.lean`. +Source-derived prose below is paraphrased. The source pin is +`anthropics/formal-math@16a5ce41304fcdaea4e0be7c132e1c361dc26da5`. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.ExactTriangle + +open Classical in +/-- Primes in the half-open window `[sqrt D / 2, sqrt D)`. -/ +noncomputable def primesInRange (D : ℕ) : Finset ℕ := + (Finset.range D).filter fun p => p.Prime ∧ Real.sqrt D / 2 ≤ (p : ℝ) ∧ (p : ℝ) < Real.sqrt D + +/-- Number of triangles whose sum vanishes modulo `p`. -/ +def countZeroMod {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) (p : ℕ) : ℕ := + (Finset.univ.filter fun t : Fin n × Fin n × Fin n => + weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≡ 0 [ZMOD (p : ℤ)]).card + +/-- A window prime minimises the modular count over all window primes. -/ +def IsSelectedPrime {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) (D p : ℕ) : Prop := + p ∈ primesInRange D ∧ ∀ q ∈ primesInRange D, + countZeroMod wAB wBC wAC p ≤ countZeroMod wAB wBC wAC q + +private theorem four_le_sqrt {D : ℕ} (hD : 16 ≤ D) : 4 ≤ Real.sqrt (D : ℝ) := + Real.le_sqrt_of_sq_le (by norm_num; exact_mod_cast hD) + +private theorem natCast_toNat_emod {p : ℕ} (hp : 0 < p) (u : ℤ) : + ((u % (p : ℤ)).toNat : ℤ) = u % (p : ℤ) := + Int.toNat_of_nonneg (Int.emod_nonneg _ (by omega)) + +private theorem toNat_emod_lt {p : ℕ} (hp : 0 < p) (u : ℤ) : (u % (p : ℤ)).toNat < p := by + have := Int.emod_nonneg u (by omega : (p : ℤ) ≠ 0) + have := Int.emod_lt_of_pos u (by omega : (0 : ℤ) < p) + omega + +/-! +### The primes in the range -/ + +/-- Proof of Theorem 17, a prime p ∈ [√D/2, √D): the primes in the range are the primes `p` with +`√D/2 ≤ p < √D`. (The bound `p < D` in the definition excludes none of them.) +-/ +theorem mem_primesInRange {D p : ℕ} : + p ∈ primesInRange D ↔ p.Prime ∧ Real.sqrt D / 2 ≤ (p : ℝ) ∧ (p : ℝ) < Real.sqrt D := by + unfold primesInRange + rw [Finset.mem_filter, Finset.mem_range] + refine ⟨fun h => h.2, fun h => ⟨?_, h⟩⟩ + -- From `p < √D` we get `p ≤ p² < D`. + have hsq : p ^ 2 < D := by exact_mod_cast (Real.lt_sqrt (Nat.cast_nonneg p)).1 h.2.2 + exact (Nat.le_self_pow two_ne_zero p).trans_lt hsq + +/-- Window membership is an exact integer-square test, suitable for native enumeration. -/ +theorem mem_primesInRange_iff_sq {D p : ℕ} : + p ∈ primesInRange D ↔ p.Prime ∧ D ≤ 4 * p ^ 2 ∧ p ^ 2 < D := by + rw [mem_primesInRange] + have hp0 : (0 : ℝ) ≤ p := Nat.cast_nonneg p + have hsqrt0 : 0 ≤ Real.sqrt (D : ℝ) := Real.sqrt_nonneg _ + have hsqrt : (Real.sqrt (D : ℝ)) ^ 2 = D := Real.sq_sqrt (Nat.cast_nonneg D) + constructor + · rintro ⟨hp, hlo, hhi⟩ + refine ⟨hp, ?_, ?_⟩ + · have h : (D : ℝ) ≤ 4 * (p : ℝ) ^ 2 := by nlinarith + exact_mod_cast h + · exact_mod_cast (Real.lt_sqrt hp0).1 hhi + · rintro ⟨hp, hlo, hhi⟩ + refine ⟨hp, ?_, (Real.lt_sqrt hp0).2 (by exact_mod_cast hhi)⟩ + have h : (D : ℝ) ≤ 4 * (p : ℝ) ^ 2 := by exact_mod_cast hlo + nlinarith + +/-- A computable finite enumeration of exactly the same half-open prime window. -/ +theorem primesInRange_eq_filter (D : ℕ) : + primesInRange D = (Finset.range D).filter fun p => p.Prime ∧ D ≤ 4 * p ^ 2 ∧ p ^ 2 < D := by + ext p + rw [mem_primesInRange_iff_sq, Finset.mem_filter, Finset.mem_range] + exact ⟨fun h => ⟨(Nat.le_self_pow two_ne_zero p).trans_lt h.2.2, h⟩, fun h => h.2⟩ + +/-- A number `r < p` is the sum of the residues of `u` and `v` in `{0, …, p − 1}`, reduced modulo +`p`, exactly if `u + v ≡ r`. +-/ +private theorem eq_add_toNat_emod_iff {p r : ℕ} (hp : p ≠ 0) (hr : r < p) (u v : ℤ) : + r = ((u % (p : ℤ)).toNat + (v % (p : ℤ)).toNat) % p ↔ u + v ≡ (r : ℤ) [ZMOD (p : ℤ)] := by + rw [← Int.natCast_inj, Int.ModEq, Int.emod_eq_of_lt (Int.natCast_nonneg r) (by exact_mod_cast hr), + eq_comm] + push_cast + rw [natCast_toNat_emod (Nat.pos_of_ne_zero hp), + natCast_toNat_emod (Nat.pos_of_ne_zero hp), + ← Int.add_emod] + +/-! +### Counting the triples with `S(a,b,c) ≡ 0 (mod p)` -/ + +section Counting + +variable {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) + +/-- A false positive of p is a triple (a,b,c) with S(a,b,c) ≠ 0 and p ∣ S(a,b,c). -/ +abbrev IsFalsePositive (p : ℕ) (t : Fin n × Fin n × Fin n) : Prop := + weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≠ 0 ∧ (p : ℤ) ∣ weightSum wAB wBC wAC t.1 t.2.1 t.2.2 + +/-- let falsePositiveCount(p) denote the number of false positives of p. -/ +def falsePositiveCount (p : ℕ) : ℕ := + (Finset.univ.filter fun t : Fin n × Fin n × Fin n => IsFalsePositive wAB wBC wAC p t).card + +/-- zeroCount, the number of zero triangles. -/ +def zeroCount : ℕ := + (Finset.univ.filter fun t : Fin n × Fin n × Fin n => weightSum wAB wBC wAC t.1 t.2.1 t.2.2 = + 0).card + +/-- A positive zero count is exactly the existing Exact Triangle decision predicate. -/ +theorem zeroCount_pos_iff : 0 < zeroCount wAB wBC wAC ↔ HasZeroTriangle wAB wBC wAC := by + classical + simp only [zeroCount, Finset.card_pos, HasZeroTriangle] + constructor + · rintro ⟨⟨a, b, c⟩, h⟩ + exact ⟨a, b, c, (Finset.mem_filter.mp h).2⟩ + · rintro ⟨a, b, c, h⟩ + exact ⟨(a, b, c), Finset.mem_filter.mpr ⟨Finset.mem_univ _, h⟩⟩ + +/-- This number is the number of false positives falsePositiveCount(p) plus zeroCount, the number of +zero triangles, which does not depend on p. +-/ +theorem countZeroMod_eq (p : ℕ) : (countZeroMod wAB wBC wAC) p = (falsePositiveCount wAB wBC + wAC) p + (zeroCount wAB wBC wAC) := by + classical + simp only [countZeroMod, falsePositiveCount, zeroCount, Finset.card_filter, ← + Finset.sum_add_distrib] + -- A triple with `p ∣ S(a,b,c)` is a zero triangle or a false positive, and not both. + refine Finset.sum_congr rfl fun t _ => ?_ + by_cases h0 : weightSum wAB wBC wAC t.1 t.2.1 t.2.2 = 0 <;> + simp [IsFalsePositive, Int.modEq_zero_iff_dvd, h0] + +/-- let P[a,c] := x^{w(a,c) mod p}, a matrix over `ℤ[x]/(x^p − 1)`. `w mod p` is the residue in +`{0, …, p − 1}`. +-/ +noncomputable def matP (p : ℕ) : Matrix (Fin n) (Fin n) (CyclicIntRing p) := + Matrix.of fun a c => CyclicIntRing.x p ^ (wAC a c % (p : ℤ)).toNat + +/-- and Q[c,b] := x^{w(b,c) mod p}. -/ +noncomputable def matQ (p : ℕ) : Matrix (Fin n) (Fin n) (CyclicIntRing p) := + Matrix.of fun c b => CyclicIntRing.x p ^ (wBC b c % (p : ℤ)).toNat + +/-- the coefficient of x^r in (PQ)[a,b] is the number of c ∈ C with w(a,c) + w(b,c) ≡ r (mod p). -/ +theorem coeff_matP_mul_matQ {p : ℕ} (hp : p ≠ 0) (a b : Fin n) (r : ℕ) (hr : r < p) : + CyclicIntRing.coeff hp r ((matP wAC p * matQ wBC p) a b) = + ((Finset.univ.filter fun c : Fin n => + wAC a c + wBC b c ≡ (r : ℤ) [ZMOD (p : ℤ)]).card : ℤ) := by + classical + -- `(PQ)[a,b]` is the sum over `c` of `x^(w(a,c) mod p + w(b,c) mod p)`, and the coefficient of + -- `x^r` in each term is 1 or 0. + rw [Matrix.mul_apply] + simp only [matP, matQ, Matrix.of_apply, ← pow_add] + rw [map_sum] + simp only [CyclicIntRing.coeff_x_pow] + rw [Finset.card_filter] + push_cast + refine Finset.sum_congr rfl fun c _ => ?_ + simp only [eq_add_toNat_emod_iff hp hr] + +/-- then falsePositiveCount(p) + zeroCount is the sum over the pairs (a,b) ∈ A × B of the +coefficient of x^{−w(a,b) mod p} in (PQ)[a,b]. +-/ +theorem falsePositiveCount_add_zeroCount_eq_sum_coeff {p : ℕ} (hp : p ≠ 0) : + (((falsePositiveCount wAB wBC wAC) p + (zeroCount wAB wBC wAC) : ℕ) : ℤ) = + ∑ a : Fin n, ∑ b : Fin n, + CyclicIntRing.coeff hp ((-wAB a b) % (p : ℤ)).toNat ((matP wAC p * matQ wBC p) a b) := by + classical + rw [← countZeroMod_eq] + unfold countZeroMod + rw [Finset.card_filter, Fintype.sum_prod_type] + push_cast + refine Finset.sum_congr rfl fun a _ => ?_ + rw [Fintype.sum_prod_type] + refine Finset.sum_congr rfl fun b _ => ?_ + have hcast := natCast_toNat_emod (Nat.pos_of_ne_zero hp) (-wAB a b) + rw [(coeff_matP_mul_matQ wBC wAC) hp a b _ (toNat_emod_lt (Nat.pos_of_ne_zero hp) _), + Finset.card_filter] + push_cast + refine Finset.sum_congr rfl fun c _ => ?_ + -- `S(a,b,c) ≡ 0` says the same as `w(a,c) + w(b,c) ≡ -w(a,b)`, and `-w(a,b)` is congruent to its + -- residue. + have hres : (-wAB a b) % (p : ℤ) ≡ -wAB a b [ZMOD (p : ℤ)] := Int.mod_modEq _ _ + have hiff : weightSum wAB wBC wAC a b c ≡ 0 [ZMOD (p : ℤ)] ↔ + wAC a c + wBC b c ≡ ((((-wAB a b) % (p : ℤ)).toNat : ℕ) : ℤ) [ZMOD (p : ℤ)] := by + rw [hcast, show wAC a c + wBC b c ≡ (-wAB a b) % (p : ℤ) [ZMOD (p : ℤ)] ↔ + wAC a c + wBC b c ≡ -wAB a b [ZMOD (p : ℤ)] from + ⟨fun h => h.trans hres, fun h => h.trans hres.symm⟩, + Int.modEq_iff_dvd, Int.modEq_iff_dvd, + show 0 - weightSum wAB wBC wAC a b c = -wAB a b - (wAC a c + wBC b c) by + simp only [weightSum]; ring] + simp only [hiff] + +end Counting + +/-- the fewer than √D primes in the range. -/ +theorem card_primesInRange_lt (D : ℕ) (hD : 16 ≤ D) : + ((primesInRange D).card : ℝ) < Real.sqrt D := by + have hpos : 0 < Real.sqrt D := by linarith [four_le_sqrt hD] + -- The primes in the range are among the numbers `1, …, ⌈√D⌉ - 1`. + have hsub : primesInRange D ⊆ Finset.Ico 1 ⌈Real.sqrt D⌉₊ := by + intro p hp + obtain ⟨hprime, -, hlt⟩ := mem_primesInRange.1 hp + exact Finset.mem_Ico.2 ⟨hprime.one_lt.le, Nat.lt_ceil.2 hlt⟩ + have hcard : (primesInRange D).card + 1 ≤ ⌈Real.sqrt D⌉₊ := by + have hle := Finset.card_le_card hsub + have hone : 1 ≤ ⌈Real.sqrt D⌉₊ := Nat.one_le_ceil_iff.2 hpos + rw [Nat.card_Ico] at hle + omega + have hcast : ((primesInRange D).card : ℝ) + 1 ≤ (⌈Real.sqrt D⌉₊ : ℝ) := by exact_mod_cast hcard + linarith [Nat.ceil_lt_add_one hpos.le] + +/-! +### Selecting the prime -/ + +section Selection + +variable {n D p : ℕ} {κ : ℝ} (wAB wBC wAC : Fin n → Fin n → ℤ) + +/-- There is a prime to select. -/ +theorem exists_isSelectedPrime (D : ℕ) (hD : 16 ≤ D) : ∃ p, IsSelectedPrime wAB wBC wAC D p := by + obtain ⟨q, hq⟩ := + Nat.exists_prime_half_le_and_lt (Real.sqrt D) (four_le_sqrt hD) + exact Finset.exists_min_image (primesInRange D) (countZeroMod wAB wBC wAC) ⟨q, + mem_primesInRange.2 hq⟩ + +/-- We select the prime with the smallest count, which is also the prime with the fewest false +positives. +-/ +theorem IsSelectedPrime.falsePositiveCount_le {wAB wBC wAC : Fin n → Fin n → ℤ} (hp : + IsSelectedPrime wAB wBC wAC D p) : + ∀ q ∈ primesInRange D, (falsePositiveCount wAB wBC wAC) p ≤ (falsePositiveCount wAB wBC wAC) + q := by + intro q hq + have h := hp.2 q hq + rw [countZeroMod_eq, countZeroMod_eq] at h + omega + +/-! +### The bound on the number of false positives of the selected prime -/ + +open Classical in +/-- a triple with S(a,b,c) ≠ 0 is a false positive exactly of the primes in the range that divide +S(a,b,c). Since 0 < |S(a,b,c)| ≤ 3n^ν and these primes are at least √D/2, there are at most +log_{√D/2}(3n^ν) of them. +-/ +theorem card_falsePositive_primes_le (hD : 16 ≤ D) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) + (t : Fin n × Fin n × Fin n) (ht : weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≠ 0) : + (((primesInRange D).filter fun p => IsFalsePositive wAB wBC wAC p t).card : ℝ) + ≤ Real.logb (Real.sqrt D / 2) (3 * (n : ℝ) ^ κ) := by + have hsqrt := four_le_sqrt hD + set divisors := (primesInRange D).filter fun p => IsFalsePositive wAB wBC wAC p t + have hmem : ∀ p ∈ divisors, + p.Prime ∧ Real.sqrt D / 2 ≤ (p : ℝ) ∧ (p : ℤ) ∣ weightSum wAB wBC wAC t.1 t.2.1 t.2.2 := by + intro p hp + obtain ⟨hrange, -, hdvd⟩ := Finset.mem_filter.1 hp + obtain ⟨hprime, hge, -⟩ := mem_primesInRange.1 hrange + exact ⟨hprime, hge, hdvd⟩ + -- The product of these primes divides `S(a,b,c)`, so it is at most `|S(a,b,c)| ≤ 3n^κ`. + have hdvd : ∏ p ∈ divisors, p ∣ (weightSum wAB wBC wAC t.1 t.2.1 t.2.2).natAbs := + Finset.prod_primes_dvd _ (fun p hp => (hmem p hp).1.prime) + fun p hp => Int.natCast_dvd.1 (hmem p hp).2.2 + have hprod : ∏ p ∈ divisors, (p : ℝ) ≤ 3 * (n : ℝ) ^ κ := + calc ∏ p ∈ divisors, (p : ℝ) = ((∏ p ∈ divisors, p : ℕ) : ℝ) := (Nat.cast_prod _ _).symm + _ ≤ ((weightSum wAB wBC wAC t.1 t.2.1 t.2.2).natAbs : ℝ) := by + exact_mod_cast Nat.le_of_dvd (Int.natAbs_pos.2 ht) hdvd + _ = ((|weightSum wAB wBC wAC t.1 t.2.1 t.2.2| : ℤ) : ℝ) := Nat.cast_natAbs _ + _ ≤ 3 * (n : ℝ) ^ κ := (abs_weightSum_le wAB wBC wAC) hT _ _ _ + -- Each of these primes is at least `√D/2`. + have hpow : (Real.sqrt D / 2) ^ divisors.card ≤ ∏ p ∈ divisors, (p : ℝ) := by + rw [← Finset.prod_const] + exact Finset.prod_le_prod₀ (fun _ _ => by linarith) fun p hp => (hmem p hp).2.1 + have hY : 0 < 3 * (n : ℝ) ^ κ := lt_of_lt_of_le (by positivity) (hpow.trans hprod) + rw [Real.le_logb_iff_rpow_le (by linarith) hY, Real.rpow_natCast] + exact hpow.trans hprod + +/-- Hence the numbers of false positives of all the primes in the range add up to at most n³ +log_{√D/2}(3n^ν). +-/ +theorem sum_falsePositiveCount_le (hD : 16 ≤ D) (hκ : 1 ≤ κ) (hT : WeightsBoundedBy wAB wBC wAC + ((n : ℝ) ^ κ)) : + ((∑ p ∈ primesInRange D, (falsePositiveCount wAB wBC wAC) p : ℕ) : ℝ) + ≤ (n : ℝ) ^ 3 * Real.logb (Real.sqrt D / 2) (3 * (n : ℝ) ^ κ) := by + classical + have hsqrt := four_le_sqrt hD + -- Count the pairs (prime, false positive of that prime) triple by triple. + have hswap : ∑ p ∈ primesInRange D, (falsePositiveCount wAB wBC wAC) p = ∑ t : Fin n × Fin n × + Fin n, + ((primesInRange D).filter fun p => IsFalsePositive wAB wBC wAC p t).card := by + simp only [falsePositiveCount, Finset.card_filter] + exact Finset.sum_comm + have hterm : ∀ t : Fin n × Fin n × Fin n, + (((primesInRange D).filter fun p => IsFalsePositive wAB wBC wAC p t).card : ℝ) + ≤ Real.logb (Real.sqrt D / 2) (3 * (n : ℝ) ^ κ) := by + intro t + by_cases ht : weightSum wAB wBC wAC t.1 t.2.1 t.2.2 = 0 + · -- A zero triangle is a false positive of no prime. + rw [Finset.filter_eq_empty_iff.2 fun p _ h => h.1 ht, Finset.card_empty, Nat.cast_zero] + have hn : (1 : ℝ) ≤ n := by exact_mod_cast Fin.pos t.1 + have hpow : (1 : ℝ) ≤ (n : ℝ) ^ κ := Real.one_le_rpow hn (by linarith) + exact Real.logb_nonneg (by linarith) (by linarith) + · -- Otherwise this is the bound of `card_falsePositive_primes_le` for this triple. + convert (card_falsePositive_primes_le wAB wBC wAC) hD hT t ht + rw [hswap, Nat.cast_sum] + calc ∑ t : Fin n × Fin n × Fin n, + (((primesInRange D).filter fun p => IsFalsePositive wAB wBC wAC p t).card : ℝ) + ≤ ∑ _t : Fin n × Fin n × Fin n, Real.logb (Real.sqrt D / 2) (3 * (n : ℝ) ^ κ) := + Finset.sum_le_sum fun t _ => hterm t + _ = (n : ℝ) ^ 3 * Real.logb (Real.sqrt D / 2) (3 * (n : ℝ) ^ κ) := by + rw [Finset.sum_const, Finset.card_univ, Fintype.card_prod, Fintype.card_prod, + Fintype.card_fin, nsmul_eq_mul] + push_cast + ring + +end Selection + +/-- By the prime number theorem there are Ω(√D/log D) primes in the range. + +NOTE. The modular reduction needs this for every `D ≥ 16`, not only for large `D`, and so it is +stated. It follows from a counting form of Bertrand's postulate; the prime number theorem is +not needed. +-/ +theorem exists_le_card_primesInRange : + ∃ c : ℝ, 0 < c ∧ ∀ D : ℕ, 16 ≤ D → + c * (Real.sqrt D / Real.log D) ≤ ((primesInRange D).card : ℝ) := by + obtain ⟨c, hc, hcount⟩ := Nat.exists_forall_mul_sqrt_div_log_le_card_primes + exact ⟨c, hc, fun D hD => hcount D hD _ + fun _ hp hge hlt => mem_primesInRange.2 ⟨hp, hge, hlt⟩⟩ + +/-- falsePositiveCount(p) is at most the average over them, multiplied out. -/ +theorem IsSelectedPrime.falsePositiveCount_mul_card_le_sum {n D p : ℕ} {wAB wBC wAC : Fin n → + Fin n → ℤ} + (hp : IsSelectedPrime wAB wBC wAC D p) : + (falsePositiveCount wAB wBC wAC) p * (primesInRange D).card ≤ ∑ q ∈ primesInRange D, + (falsePositiveCount wAB wBC wAC) q := by + rw [mul_comm, ← smul_eq_mul] + exact Finset.card_nsmul_le_sum _ _ _ hp.falsePositiveCount_le + +/-- log D = O(log(√D/2)) for D ≥ 16. -/ +theorem log_le_four_mul_log_sqrt_div_two (D : ℕ) (hD : 16 ≤ D) : + Real.log D ≤ 4 * Real.log (Real.sqrt D / 2) := by + have hD16 : (16 : ℝ) ≤ D := by exact_mod_cast hD + have hsqrt : 0 < Real.sqrt D := by linarith [four_le_sqrt hD] + -- `log(√D/2) = (log D)/2 - log 2`, and `4 log 2 = log 16 ≤ log D`. + have hlog16 : 4 * Real.log 2 ≤ Real.log D := by + have h := Real.log_le_log (by norm_num) (show (2 : ℝ) ^ 4 ≤ D by linarith) + rw [Real.log_pow] at h + exact_mod_cast h + rw [Real.log_div hsqrt.ne' (by norm_num), Real.log_sqrt (by linarith)] + linarith [hlog16] + +/-- falsePositiveCount(p) is at most the average over them, so [...] falsePositiveCount(p) = O(n³ +log(3n^ν)/√D): if there are at least `c √D/log D` primes in the range, then +`falsePositiveCount(p) ≤ (4/c) n³ log(3n^κ)/√D`. +-/ +theorem falsePositiveCount_le_of_le_card_primesInRange {n D p : ℕ} {κ c : ℝ} + (wAB wBC wAC : Fin n → Fin n → ℤ) (hc : 0 < c) + (hcard : c * (Real.sqrt D / Real.log D) ≤ ((primesInRange D).card : ℝ)) (hD : 16 ≤ D) + (hDn : D ≤ n) (hκ : 1 ≤ κ) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) (hp : + IsSelectedPrime wAB wBC wAC D p) : + ((falsePositiveCount wAB wBC wAC) p : ℝ) ≤ 4 / c * ((n : ℝ) ^ 3 * Real.log (3 * (n : ℝ) ^ κ) + / Real.sqrt D) := by + -- The quantities of the paper's sentence: `r = √D`, `Y = 3n^κ`, and the number of primes in the + -- range. + set r := Real.sqrt D + set Y := 3 * (n : ℝ) ^ κ with hY + set numPrimes := (primesInRange D).card + have hr : 4 ≤ r := four_le_sqrt hD + have hn : (16 : ℝ) ≤ n := by exact_mod_cast hD.trans hDn + have hlogD_pos : 0 < Real.log D := Real.log_pos (by exact_mod_cast (by omega : 1 < D)) + have hlogr_pos : 0 < Real.log (r / 2) := Real.log_pos (by linarith) + have hlogY_nonneg : 0 ≤ Real.log Y := by + have : 1 ≤ (n : ℝ) ^ κ := Real.one_le_rpow (by linarith) (by linarith) + exact Real.log_nonneg (by rw [hY]; linarith) + -- The selected prime is at most average, after multiplying by the number of primes. + -- False positives of all primes in the + -- range add up to at most `n³ log_{r/2} Y`. + have average : ((falsePositiveCount wAB wBC wAC) p : ℝ) * numPrimes ≤ (n : ℝ) ^ 3 * Real.logb + (r / 2) Y := + le_trans (by exact_mod_cast hp.falsePositiveCount_mul_card_le_sum) + ((sum_falsePositiveCount_le wAB wBC wAC) hD hκ hT) + -- The base of the logarithm changes: `log_{r/2} Y ≤ 4 log Y / log D`. + have hlogD_le : Real.log D ≤ 4 * Real.log (r / 2) := log_le_four_mul_log_sqrt_div_two D hD + have changeBase : Real.logb (r / 2) Y ≤ 4 * Real.log Y / Real.log D := by + rw [Real.logb, div_le_div_iff₀ hlogr_pos hlogD_pos] + linarith [mul_le_mul_of_nonneg_left hlogD_le hlogY_nonneg] + rw [show 4 / c * ((n : ℝ) ^ 3 * Real.log Y / r) = 4 * (n : ℝ) ^ 3 * Real.log Y / (c * r) by + field_simp, le_div_iff₀ (by positivity)] + calc ((falsePositiveCount wAB wBC wAC) p : ℝ) * (c * r) = (falsePositiveCount wAB wBC wAC) p * + (c * (r / Real.log D)) * Real.log D := by field_simp + _ ≤ (falsePositiveCount wAB wBC wAC) p * numPrimes * Real.log D := by gcongr + _ ≤ (n : ℝ) ^ 3 * Real.logb (r / 2) Y * Real.log D := by gcongr + _ ≤ (n : ℝ) ^ 3 * (4 * Real.log Y / Real.log D) * Real.log D := by gcongr + _ = 4 * (n : ℝ) ^ 3 * Real.log Y := by field_simp + +/-- falsePositiveCount(p) = O(n³ log(3n^ν)/√D) = O(ν n³ log n/√D). -/ +theorem exists_falsePositiveCount_le : + ∃ C : ℝ, 1 ≤ C ∧ ∀ {n D p : ℕ} {κ : ℝ}, 16 ≤ D → D ≤ n → 1 ≤ κ → ∀ wAB wBC wAC : Fin n → Fin + n → ℤ, + WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ) → IsSelectedPrime wAB wBC wAC D p → + ((falsePositiveCount wAB wBC wAC) p : ℝ) ≤ C * (κ * (n : ℝ) ^ 3 * Real.log n / Real.sqrt + D) := by + obtain ⟨c, hc, hcard⟩ := exists_le_card_primesInRange + refine ⟨max 1 (8 / c), le_max_left _ _, fun {n D p κ} hD hDn hκ wAB wBC wAC hT hp => ?_⟩ + have hn : (16 : ℝ) ≤ n := by exact_mod_cast hD.trans hDn + have hlogn : 0 ≤ Real.log n := Real.log_nonneg (by linarith) + -- `log(3n^κ) = log 3 + κ log n ≤ 2κ log n`, because `3 ≤ n` and `κ ≥ 1`. + have hlog : Real.log (3 * (n : ℝ) ^ κ) ≤ 2 * κ * Real.log n := by + have hlog3 : Real.log 3 ≤ κ * Real.log n := + (Real.log_le_log (by norm_num) (by linarith)).trans (le_mul_of_one_le_left hlogn hκ) + rw [Real.log_mul (by norm_num) (by positivity), Real.log_rpow (by linarith)] + linarith [hlog3] + calc ((falsePositiveCount wAB wBC wAC) p : ℝ) ≤ 4 / c * ((n : ℝ) ^ 3 * Real.log (3 * (n : ℝ) ^ + κ) / Real.sqrt D) := + (falsePositiveCount_le_of_le_card_primesInRange wAB wBC wAC) hc (hcard D hD) hD hDn hκ hT hp + _ ≤ 4 / c * ((n : ℝ) ^ 3 * (2 * κ * Real.log n) / Real.sqrt D) := by gcongr + _ = 8 / c * (κ * (n : ℝ) ^ 3 * Real.log n / Real.sqrt D) := by ring + _ ≤ max 1 (8 / c) * (κ * (n : ℝ) ^ 3 * Real.log n / Real.sqrt D) := by + gcongr + exact le_max_right _ _ + +/-- A name for the constant in falsePositiveCount(p) = [...] = O(ν n³ log n/√D), for the bounds on +running times that are built on this one. +-/ +noncomputable def Hashing.falsePositiveConst : ℝ := Classical.choose exists_falsePositiveCount_le + +/-- The constant is at least 1. -/ +theorem Hashing.one_le_falsePositiveConst : 1 ≤ Hashing.falsePositiveConst := + (Classical.choose_spec exists_falsePositiveCount_le).1 + +/-- falsePositiveCount(p) = [...] = O(ν n³ log n/√D), with `Hashing.falsePositiveConst` as the +constant. +-/ +theorem Hashing.falsePositiveCount_le_const_mul {n D p : ℕ} {κ : ℝ} (hD : 16 ≤ D) (hDn : D ≤ n) + (hκ : 1 ≤ κ) (wAB wBC wAC : Fin n → Fin n → ℤ) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) + (hp : IsSelectedPrime wAB wBC wAC D p) : + ((falsePositiveCount wAB wBC wAC) p : ℝ) ≤ Hashing.falsePositiveConst * (κ * (n : ℝ) ^ 3 * + Real.log n / Real.sqrt D) := + (Classical.choose_spec exists_falsePositiveCount_le).2 hD hDn hκ wAB wBC wAC hT hp + +end Cslib.Algorithms.ExactTriangle diff --git a/Complexity/Algorithms/ExactTriangle/Weights.lean b/Complexity/Algorithms/ExactTriangle/Weights.lean new file mode 100644 index 0000000000..d358d0cd80 --- /dev/null +++ b/Complexity/Algorithms/ExactTriangle/Weights.lean @@ -0,0 +1,61 @@ +/- +Copyright (c) 2026 Anthropic, PBC, Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Anthropic, PBC, Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.ExactTriangle.Defs +public import Mathlib.Algebra.Order.Group.Abs +public import Mathlib.Basic.Real.Basic +import Mathlib.Tactic.Linarith + +/-! +# Integer triangle sums and magnitude bounds + +## Main definitions + +* `weightSum`: the sum of the three supplied integer edge weights. +* `WeightsBoundedBy`: a shared real bound on their absolute values. + +## Main results + +* `abs_weightSum_le`: a triangle's magnitude is at most three times the edge bound. + +## References + +Adapted from the weight sum and bound definitions in [Anthropic's +PaperStatements.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/PaperStatements.lean#L893-L914) +and the triangle inequality in its +[Hashing.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/Theorem17/Hashing.lean#L246-L255). +The existing `HasZeroTriangle` remains the decision predicate. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.ExactTriangle + +variable {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) + +/-- Sum of the A–B, B–C and A–C weights at a triangle. -/ +def weightSum (a b c : Fin n) : ℤ := wAB a b + wBC b c + wAC a c + +/-- Every supplied edge weight has absolute value at most `U`. -/ +def WeightsBoundedBy (U : ℝ) : Prop := + (∀ a b, ((|wAB a b| : ℤ) : ℝ) ≤ U) ∧ + (∀ b c, ((|wBC b c| : ℤ) : ℝ) ≤ U) ∧ ∀ a c, ((|wAC a c| : ℤ) : ℝ) ≤ U + +/-- The three edge bounds imply a threefold bound on every triangle sum. -/ +theorem abs_weightSum_le {U : ℝ} (h : WeightsBoundedBy wAB wBC wAC U) (a b c : Fin n) : + ((|weightSum wAB wBC wAC a b c| : ℤ) : ℝ) ≤ 3 * U := by + obtain ⟨hAB, hBC, hAC⟩ := h + have habs : |weightSum wAB wBC wAC a b c| ≤ |wAB a b| + |wBC b c| + |wAC a c| := + abs_add_three _ _ _ + have hcast : ((|weightSum wAB wBC wAC a b c| : ℤ) : ℝ) ≤ + ((|wAB a b| : ℤ) : ℝ) + ((|wBC b c| : ℤ) : ℝ) + ((|wAC a c| : ℤ) : ℝ) := by + exact_mod_cast habs + linarith [hAB a b, hBC b c, hAC a c] + +end Cslib.Algorithms.ExactTriangle diff --git a/Complexity/Foundations/NumberTheory/PrimeWindow.lean b/Complexity/Foundations/NumberTheory/PrimeWindow.lean new file mode 100644 index 0000000000..dbeb6390e2 --- /dev/null +++ b/Complexity/Foundations/NumberTheory/PrimeWindow.lean @@ -0,0 +1,346 @@ +/- +Copyright (c) 2020 Patrick Stevens. All rights reserved. +Copyright (c) 2026 Anthropic, PBC. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Patrick Stevens, Bolton Bailey, Anthropic, PBC +-/ + +module +-- upstream: candidate + +public import Complexity.Init +public import Mathlib.Analysis.SpecialFunctions.Log.Base +public import Mathlib.Analysis.SpecialFunctions.Pow.Real +public import Mathlib.NumberTheory.PrimeCounting +import Mathlib.Analysis.Complex.ExponentialBounds +import Mathlib.Data.Nat.Prime.Factorial +import Mathlib.NumberTheory.Bertrand +import Mathlib.NumberTheory.Primorial +import Mathlib.Analysis.Convex.SpecificFunctions.Deriv + +/-! +# Explicit prime counts in a window of ratio two + +For every real `x ≥ 4`, the half-open window `[x/2, x)` contains at least `c * x / log x` primes +for one universal positive constant `c`. The proof counts prime factors of a central binomial +coefficient, so no prime-number-theorem hypothesis is assumed. + +## Main results + +* `Nat.exists_prime_half_le_and_lt`: a prime exists in every such window. +* `Nat.exists_forall_mul_div_log_le_card_primes`: the uniform counting bound. +* `Nat.exists_forall_mul_sqrt_div_log_le_card_primes`: the bound at `x = sqrt D`, for all `D ≥ + 16`. + +## References + +Adapted from [Anthropic's +PrimesInWindow.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Util/PrimesInWindow.lean). +Its central-binomial and concavity estimates adapt `Mathlib/NumberTheory/Bertrand.lean`, +originally by Patrick Stevens and Bolton Bailey; that source credit is retained. +-/ + +@[expose] public section + +open Finset + +namespace Nat + +private theorem one_le_log_of_three_le {x : ℝ} (hx : 3 ≤ x) : 1 ≤ Real.log x := + (Real.le_log_iff_exp_le (by linarith)).2 (Real.exp_one_lt_three.le.trans hx) + + +/-- The primes `p` with `n < p ≤ 2n`: the primes that Bertrand's postulate is about. -/ +private def bertrandPrimes (n : ℕ) : Finset ℕ := (Ioc n (2 * n)).filter Nat.Prime + +/-- Membership in `bertrandPrimes`, spelled out. -/ +private theorem mem_bertrandPrimes {n p : ℕ} : + p ∈ bertrandPrimes n ↔ (n < p ∧ p ≤ 2 * n) ∧ p.Prime := by + simp only [bertrandPrimes, mem_filter, mem_Ioc] + +/-! +## 1. Erdős's upper bound for the central binomial coefficient -/ + +/-- Every prime power in `C(2n, n)` is at most `2n`, so `|s|` of them multiply to at most +`(2n)^|s|`. +-/ +private theorem prod_pow_factorization_le {n : ℕ} (hn : 0 < n) (s : Finset ℕ) : + ∏ p ∈ s, p ^ n.centralBinom.factorization p ≤ (2 * n) ^ s.card := + prod_le_pow_card _ _ _ fun _ _ => Nat.pow_factorization_choose_le (by omega) + +/-- The primes in `(√(2n), n]` contribute at most `4^(2n/3)` to `C(2n, n)`: their exponent is at +most 1, it is 0 above `2n/3`, and the product of the primes up to `2n/3` is at most `4^(2n/3)`. +Adapted from Mathlib, Mathlib/NumberTheory/Bertrand.lean, Apache 2.0. +-/ +private theorem prod_middle_le {n : ℕ} (hn : 2 < n) (s : Finset ℕ) + (hs : ∀ p ∈ s, p.Prime ∧ Nat.sqrt (2 * n) < p ∧ p ≤ n) : + ∏ p ∈ s, p ^ n.centralBinom.factorization p ≤ 4 ^ (2 * n / 3) := + calc ∏ p ∈ s, p ^ n.centralBinom.factorization p + ≤ ∏ p ∈ s, (if p ≤ 2 * n / 3 then p else 1) := by + refine prod_le_prod fun p hp => ?_ + obtain ⟨hprime, hsqrt, hpn⟩ := hs p hp + split_ifs with h3 + · exact (pow_right_mono₀ hprime.one_lt.le + (Nat.factorization_choose_le_one (Nat.sqrt_lt'.mp hsqrt))).trans (pow_one p).le + · rw [Nat.factorization_centralBinom_of_two_mul_self_lt_three_mul hn hpn (by omega), + pow_zero] + _ = ∏ p ∈ s.filter (· ≤ 2 * n / 3), p := (prod_filter _ _).symm + _ ≤ primorial (2 * n / 3) := by + refine prod_le_prod_of_subset_of_one_le (fun p hp => ?_) + fun p hp _ => (mem_filter.1 hp).2.one_lt.le + obtain ⟨hps, hp3⟩ := mem_filter.1 hp + exact mem_filter.2 ⟨mem_range.2 (by omega), (hs p hps).1⟩ + _ ≤ 4 ^ (2 * n / 3) := primorial_le_four_pow _ + +/-- Erdős's bound, keeping the primes of `(n, 2n]`: `C(2n, n) ≤ (2n)^√(2n) · 4^(2n/3) · (2n)^k`, +where `k` is the number of primes in `(n, 2n]`. Mathlib proves the case `k = 0` on the way to +Bertrand's postulate (`centralBinom_le_of_no_bertrand_prime`). Adapted from Mathlib, +Mathlib/NumberTheory/Bertrand.lean, Apache 2.0. +-/ +private theorem centralBinom_le (n : ℕ) (hn : 2 < n) : + n.centralBinom ≤ + (2 * n) ^ Nat.sqrt (2 * n) * 4 ^ (2 * n / 3) * (2 * n) ^ (bertrandPrimes n).card := by + -- `C(2n, n)` is the product of its prime powers over the set `T` of the primes up to `2n`. + set T := (range (2 * n + 1)).filter Nat.Prime with hT + have hprod : ∏ p ∈ T, p ^ n.centralBinom.factorization p = n.centralBinom := + (prod_filter_of_ne fun p _ h => by + contrapose h + rw [Nat.factorization_eq_zero_of_not_prime _ h, pow_zero]).trans + n.prod_pow_factorization_centralBinom + -- Split the primes at `n`, and those up to `n` at `√(2n)`. + rw [← hprod, ← prod_filter_mul_prod_filter_not T (· ≤ n), + ← prod_filter_mul_prod_filter_not (T.filter (· ≤ n)) (· ≤ Nat.sqrt (2 * n))] + refine mul_le_mul' (mul_le_mul' ?_ ?_) ?_ + · -- There are at most `√(2n)` primes up to `√(2n)`. + have hsub : (T.filter (· ≤ n)).filter (· ≤ Nat.sqrt (2 * n)) ⊆ Icc 1 (Nat.sqrt (2 * n)) := by + intro p hp + simp only [hT, mem_filter, mem_range] at hp + obtain ⟨⟨⟨-, hprime⟩, -⟩, hsqrt⟩ := hp + exact mem_Icc.mpr ⟨hprime.one_lt.le, hsqrt⟩ + exact (prod_pow_factorization_le (by omega) _).trans (pow_right_mono₀ (by omega) + ((card_le_card hsub).trans_eq (by rw [Nat.card_Icc, Nat.add_sub_cancel]))) + · -- The primes in `(√(2n), n]`. + refine prod_middle_le hn _ fun p hp => ?_ + simp only [hT, mem_filter, mem_range, not_le] at hp + obtain ⟨⟨⟨-, hprime⟩, hpn⟩, hsqrt⟩ := hp + exact ⟨hprime, hsqrt, hpn⟩ + · -- The primes above `n` are those of `(n, 2n]`. + refine (prod_pow_factorization_le (by omega) _).trans + (pow_right_mono₀ (by omega) (card_le_card fun p hp => ?_)) + simp only [hT, mem_filter, mem_range] at hp + obtain ⟨⟨hp2n, hprime⟩, hnp⟩ := hp + exact mem_bertrandPrimes.2 ⟨⟨by omega, by omega⟩, hprime⟩ + +/-! +## 2. The elementary estimate `n^6 · (2n)^(6√(2n)) ≤ 2^n` for `n ≥ 16562` -/ + +/-- For `a ≥ 0` the function `log x + √(2x) · log (2x) - a x` is concave for `x > 1/2`. Adapted from +Mathlib, Mathlib/NumberTheory/Bertrand.lean, Apache 2.0 (there inside +`Bertrand.real_main_inequality`, for `a = log 4 / 3`). +-/ +private theorem concaveOn_log_add_sqrt_mul_log_sub (a : ℝ) (ha : 0 ≤ a) : + ConcaveOn ℝ (Set.Ioi (1 / 2 : ℝ)) + (fun x : ℝ => Real.log x + √(2 * x) * Real.log (2 * x) - a * x) := by + apply ConcaveOn.sub + · apply ConcaveOn.add + · exact strictConcaveOn_log_Ioi.concaveOn.subset + (Set.Ioi_subset_Ioi (by norm_num)) (convex_Ioi _) + · have h := strictConcaveOn_sqrt_mul_log_Ioi.concaveOn.comp_linearMap + ((2 : ℝ) • (LinearMap.id : ℝ →ₗ[ℝ] ℝ)) + have hset : (⇑((2 : ℝ) • (LinearMap.id : ℝ →ₗ[ℝ] ℝ))) ⁻¹' Set.Ioi 1 = Set.Ioi (1 / 2) := by + ext x + simp only [Set.mem_Ioi, Set.mem_preimage, LinearMap.smul_apply, LinearMap.id_coe, id_eq, + smul_eq_mul] + constructor <;> intro h <;> linarith + rw [hset] at h + exact h + · exact ConvexOn.smul ha (convexOn_id (convex_Ioi _)) + +/-- `log x + √(2x) · log (2x) ≤ (log 2 / 6) · x` for `x ≥ 16562`. The difference of the two sides +is concave, it is nonnegative at `x = 2` and nonpositive at `x = 16562 = 2 · 91²` (where `√(2x) += 182`; the only numerical input is `91^13 ≤ 2^85`), so it is nonpositive from there on. Adapted +from Mathlib, Mathlib/NumberTheory/Bertrand.lean, Apache 2.0 (the same argument proves +`Bertrand.real_main_inequality` there, with other numbers). +-/ +private theorem log_add_sqrt_mul_log_le {x : ℝ} (hx : 16562 ≤ x) : + Real.log x + √(2 * x) * Real.log (2 * x) ≤ Real.log 2 / 6 * x := by + have hlog2 : 0 < Real.log 2 := Real.log_pos (by norm_num) + have hat2 : (0 : ℝ) ≤ + Real.log 2 + √(2 * 2) * Real.log (2 * 2) - Real.log 2 / 6 * 2 := by + rw [Real.sqrt_mul_self (by norm_num), Real.log_mul (by norm_num) (by norm_num)] + linarith + have hat16562 : Real.log 16562 + √(2 * 16562) * Real.log (2 * 16562) + - Real.log 2 / 6 * 16562 ≤ 0 := by + have hlog91 : 13 * Real.log 91 ≤ 85 * Real.log 2 := by + have h := Real.log_le_log (by positivity) (show (91 : ℝ) ^ 13 ≤ 2 ^ 85 by norm_num) + rw [Real.log_pow, Real.log_pow] at h + exact_mod_cast h + have hlogx : Real.log 16562 = Real.log 2 + 2 * Real.log 91 := by + rw [show (16562 : ℝ) = 2 * 91 ^ 2 by norm_num, Real.log_mul (by norm_num) (by norm_num), + Real.log_pow] + norm_num + have hlog2x : Real.log (182 * 182) = 2 * Real.log 2 + 2 * Real.log 91 := by + rw [show (182 * 182 : ℝ) = 2 ^ 2 * 91 ^ 2 by norm_num, + Real.log_mul (by norm_num) (by norm_num), Real.log_pow, Real.log_pow] + norm_num + rw [show (2 * 16562 : ℝ) = 182 * 182 by norm_num, Real.sqrt_mul_self (by norm_num), hlogx, + hlog2x] + -- `366 log 91 ≤ 366 · (85/13) log 2 < (16562/6 - 365) log 2` + linarith only [hlog91, hlog2] + have hconcave := concaveOn_log_add_sqrt_mul_log_sub (Real.log 2 / 6) (by positivity) + have hatx := hconcave.right_le_of_le_left'' + (x := 2) (y := 16562) (z := x) (by norm_num [Set.mem_Ioi]) (by rw [Set.mem_Ioi]; linarith) + (by norm_num) hx (hat16562.trans hat2) + linarith only [hatx, hat16562] + +/-- The same estimate for natural numbers, without logarithms. -/ +private theorem pow_mul_pow_sqrt_le_two_pow {n : ℕ} (hn : 16562 ≤ n) : + n ^ 6 * (2 * n) ^ (6 * Nat.sqrt (2 * n)) ≤ 2 ^ n := by + have hx : (16562 : ℝ) ≤ n := by exact_mod_cast hn + have hreal := log_add_sqrt_mul_log_le hx + have hsqrt : (Nat.sqrt (2 * n) : ℝ) ≤ √(2 * n) := by + exact_mod_cast Real.nat_sqrt_le_real_sqrt (a := 2 * n) + have hlog : 0 ≤ Real.log (2 * n) := Real.log_nonneg (by linarith) + have hpos : (0 : ℝ) < n := by linarith + rw [← Nat.cast_le (α := ℝ)] + push_cast + rw [← Real.log_le_log_iff (by positivity) (by positivity), + Real.log_mul (by positivity) (by positivity), Real.log_pow, Real.log_pow, Real.log_pow] + push_cast + -- `6 log n + 6 ⌊√(2n)⌋ log(2n) ≤ 6 (log n + √(2n) log(2n)) ≤ n log 2` + linarith [hreal, mul_le_mul_of_nonneg_right hsqrt hlog] + +/-! +## 3. The count for `n ≥ 16562` -/ + +/-- `2^n ≤ (2n)^(2k)` for `n ≥ 16562`. If it failed, the sixth power of `4^n < n · (2n)^√(2n) · +4^(2n/3) · (2n)^k` would give `2^(12n) < (n^6 (2n)^(6√(2n))) · 2^(8n) · 2^(3n) ≤ 2^(12n)`. +-/ +private theorem two_pow_le_of_large {n : ℕ} (hn : 16562 ≤ n) : + 2 ^ n ≤ (2 * n) ^ (2 * (bertrandPrimes n).card) := by + by_contra hcon + have hfew : (2 * n) ^ (2 * (bertrandPrimes n).card) ≤ 2 ^ n := (not_le.mp hcon).le + have herdos : 4 ^ n < n * ((2 * n) ^ Nat.sqrt (2 * n) * 4 ^ (2 * n / 3) + * (2 * n) ^ (bertrandPrimes n).card) := + (Nat.four_pow_lt_mul_centralBinom n (by omega)).trans_le + (Nat.mul_le_mul_left _ (centralBinom_le n (by omega))) + have hsmall := pow_mul_pow_sqrt_le_two_pow hn + have hmiddle : 4 ^ (6 * (2 * n / 3)) ≤ 4 ^ (4 * n) := + Nat.pow_le_pow_right (by norm_num) (by omega) + refine lt_irrefl ((4 ^ n) ^ 6) ?_ + calc (4 ^ n) ^ 6 + < (n * ((2 * n) ^ Nat.sqrt (2 * n) * 4 ^ (2 * n / 3) + * (2 * n) ^ (bertrandPrimes n).card)) ^ 6 := Nat.pow_lt_pow_left herdos (by norm_num) + _ = (n ^ 6 * (2 * n) ^ (6 * Nat.sqrt (2 * n))) * 4 ^ (6 * (2 * n / 3)) + * ((2 * n) ^ (2 * (bertrandPrimes n).card)) ^ 3 := by ring + _ ≤ 2 ^ n * 4 ^ (4 * n) * (2 ^ n) ^ 3 := by gcongr + _ = (4 ^ n) ^ 6 := by + rw [mul_comm 4 n, pow_mul 4 n 4, show (4 : ℕ) ^ n = 2 ^ n * 2 ^ n by rw [← mul_pow]; rfl] + ring + +/-- The count of step 3 with natural logarithms: `n log 2 ≤ 2 k log (2n)`. -/ +private theorem mul_log_two_le_card_bertrandPrimes_mul_log {n : ℕ} (hn : 16562 ≤ n) : + (n : ℝ) * Real.log 2 ≤ 2 * (bertrandPrimes n).card * Real.log (2 * n) := by + have hpow : (2 : ℝ) ^ n ≤ (2 * (n : ℝ)) ^ (2 * (bertrandPrimes n).card) := by + exact_mod_cast two_pow_le_of_large hn + have hlog := Real.log_le_log (by positivity) hpow + rw [Real.log_pow, Real.log_pow] at hlog + exact_mod_cast hlog + +/-! +## 4. Half-open windows `[x/2, x)` with real endpoints -/ + +/-- For a real `x ≥ 4` the natural number `n = ⌈x/2⌉ - 1`, written below as `(⌈x⌉ - 1) / 2`, +satisfies `n ≥ 1`, `n ≥ x/4` and `2n < x`, and every natural number in `(n, 2n]` lies in `[x/2, +x)`. +-/ +private theorem exists_nat_window {x : ℝ} (hx : 4 ≤ x) : + ∃ n : ℕ, 1 ≤ n ∧ x / 4 ≤ (n : ℝ) ∧ 2 * (n : ℝ) < x ∧ + ∀ p : ℕ, n < p → p ≤ 2 * n → x / 2 ≤ (p : ℝ) ∧ (p : ℝ) < x := by + have hle : x ≤ (⌈x⌉₊ : ℝ) := Nat.le_ceil x + have hlt : (⌈x⌉₊ : ℝ) < x + 1 := Nat.ceil_lt_add_one (by linarith) + have hceil : 4 ≤ ⌈x⌉₊ := by exact_mod_cast hx.trans hle + obtain ⟨n, hn⟩ : ∃ n, n = (⌈x⌉₊ - 1) / 2 := ⟨_, rfl⟩ + have hlow : (⌈x⌉₊ : ℝ) ≤ 2 * n + 2 := by exact_mod_cast (by omega : ⌈x⌉₊ ≤ 2 * n + 2) + have hhigh : (2 * n + 1 : ℝ) ≤ (⌈x⌉₊ : ℝ) := by exact_mod_cast (by omega : 2 * n + 1 ≤ ⌈x⌉₊) + refine ⟨n, by omega, by linarith, by linarith, fun p hnp hp2n => ?_⟩ + have hnp' : (n : ℝ) + 1 ≤ p := by exact_mod_cast hnp + have hp2n' : (p : ℝ) ≤ 2 * n := by exact_mod_cast hp2n + exact ⟨by linarith, by linarith⟩ + +/-- For every real `x ≥ 4` there is a prime `p` with `x/2 ≤ p < x`, by Bertrand's postulate. -/ +theorem exists_prime_half_le_and_lt (x : ℝ) (hx : 4 ≤ x) : + ∃ p : ℕ, p.Prime ∧ x / 2 ≤ (p : ℝ) ∧ (p : ℝ) < x := by + obtain ⟨n, hn1, -, -, hwindow⟩ := exists_nat_window hx + obtain ⟨p, hprime, hnp, hp2n⟩ := Nat.exists_prime_lt_and_le_two_mul n (by omega) + exact ⟨p, hprime, hwindow p hnp hp2n⟩ + +/-- If `(n, 2n]` lies in the window and `n ≥ x/4` is large, the window has at least `x / (12 log x)` +primes: `(x/4) log 2 ≤ n log 2 ≤ 2 k log (2n) ≤ 2 k log x`, and `log 2 ≥ 2/3`. +-/ +private theorem div_log_le_card_of_large {x : ℝ} {n : ℕ} {S : Finset ℕ} (hx : 4 ≤ x) + (hn : 16562 ≤ n) (hn4 : x / 4 ≤ (n : ℝ)) (h2n : 2 * (n : ℝ) < x) (hsub : bertrandPrimes n ⊆ S) : + x / (12 * Real.log x) ≤ (S.card : ℝ) := by + have hcard : ((bertrandPrimes n).card : ℝ) ≤ S.card := by exact_mod_cast card_le_card hsub + have hn1 : (1 : ℝ) ≤ n := by exact_mod_cast (by omega : 1 ≤ n) + have hlog : Real.log (2 * n) ≤ Real.log x := Real.log_le_log (by linarith) h2n.le + have hlog0 : 0 ≤ Real.log (2 * n) := Real.log_nonneg (by linarith) + have hlog2 : 2 / 3 ≤ Real.log 2 := by linarith [Real.log_two_gt_d9] + rw [div_le_iff₀ (mul_pos (by norm_num) (Real.log_pos (by linarith)))] + calc x = 6 * (x / 4 * (2 / 3)) := by ring + _ ≤ 6 * (n * Real.log 2) := by gcongr + _ ≤ 6 * (2 * (bertrandPrimes n).card * Real.log (2 * n)) := + mul_le_mul_of_nonneg_left (mul_log_two_le_card_bertrandPrimes_mul_log hn) (by norm_num) + _ ≤ 6 * (2 * S.card * Real.log x) := by gcongr + _ = S.card * (12 * Real.log x) := by ring + +/-- There is a `c > 0` such that, for every real `x ≥ 4`, every finite set of natural numbers that +contains all primes `p` with `x/2 ≤ p < x` has at least `c · x / log x` elements. (This form +applies to the set of those primes however it is written down.) +-/ +theorem exists_forall_mul_div_log_le_card_primes : + ∃ c : ℝ, 0 < c ∧ ∀ x : ℝ, 4 ≤ x → ∀ S : Finset ℕ, + (∀ p : ℕ, p.Prime → x / 2 ≤ (p : ℝ) → (p : ℝ) < x → p ∈ S) → + c * (x / Real.log x) ≤ (S.card : ℝ) := by + -- `66248 = 4 · 16562`: for `x` below it one prime is enough + refine ⟨1 / 66248, by norm_num, fun x hx S hS => ?_⟩ + obtain ⟨n, hn1, hn4, h2n, hwindow⟩ := exists_nat_window hx + have hlogx : 1 ≤ Real.log x := one_le_log_of_three_le (by linarith) + rcases le_or_gt 16562 n with hlarge | hsmall + · -- The primes of `(n, 2n]` lie in the window. + have hsub : bertrandPrimes n ⊆ S := fun p hp => by + obtain ⟨⟨hnp, hp2n⟩, hprime⟩ := mem_bertrandPrimes.1 hp + exact hS p hprime (hwindow p hnp hp2n).1 (hwindow p hnp hp2n).2 + calc 1 / 66248 * (x / Real.log x) ≤ 1 / 12 * (x / Real.log x) := by gcongr; norm_num + _ = x / (12 * Real.log x) := by field_simp + _ ≤ (S.card : ℝ) := div_log_le_card_of_large hx hlarge hn4 h2n hsub + · -- Here `x ≤ 4n < 66248`, and one prime is enough. + obtain ⟨p, hprime, hlow, hhigh⟩ := exists_prime_half_le_and_lt x hx + have hone : (1 : ℝ) ≤ S.card := by + exact_mod_cast Finset.card_pos.mpr ⟨p, hS p hprime hlow hhigh⟩ + have hn : (n : ℝ) ≤ 16562 := by exact_mod_cast hsmall.le + calc 1 / 66248 * (x / Real.log x) ≤ 1 / 66248 * (x / 1) := by gcongr + _ ≤ 1 := by linarith + _ ≤ (S.card : ℝ) := hone + +/-! +## 5. The window `[√D/2, √D)` used in the modular reduction -/ + +/-- There is a `c > 0` such that, for every natural number `D ≥ 16`, every finite set of natural +numbers that contains all primes `p` with `√D/2 ≤ p < √D` has at least `c · √D / log D` +elements. +-/ +theorem exists_forall_mul_sqrt_div_log_le_card_primes : + ∃ c : ℝ, 0 < c ∧ ∀ D : ℕ, 16 ≤ D → ∀ S : Finset ℕ, + (∀ p : ℕ, p.Prime → √D / 2 ≤ (p : ℝ) → (p : ℝ) < √D → p ∈ S) → + c * (√D / Real.log D) ≤ (S.card : ℝ) := by + obtain ⟨c, hc, hcount⟩ := exists_forall_mul_div_log_le_card_primes + refine ⟨c, hc, fun D hD S hS => ?_⟩ + have hlog : 0 < Real.log D := Real.log_pos (by exact_mod_cast (by omega : 1 < D)) + calc c * (√D / Real.log D) ≤ c * (√D / (Real.log D / 2)) := by + gcongr + linarith + _ = c * (√D / Real.log (√D)) := by rw [Real.log_sqrt (Nat.cast_nonneg D)] + _ ≤ (S.card : ℝ) := hcount _ ((Real.le_sqrt_of_sq_le (by norm_num; exact_mod_cast hD))) S hS + +end Nat diff --git a/Complexity/Foundations/RingTheory/CyclicIntRing.lean b/Complexity/Foundations/RingTheory/CyclicIntRing.lean new file mode 100644 index 0000000000..94ae6d738f --- /dev/null +++ b/Complexity/Foundations/RingTheory/CyclicIntRing.lean @@ -0,0 +1,84 @@ +/- +Copyright (c) 2026 Anthropic, PBC, Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Anthropic, PBC, Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Init +public import Mathlib.RingTheory.AdjoinRoot + +/-! +# Integer polynomials modulo `X^p - 1` + +## Main definitions + +* `CyclicIntRing`: the integer polynomial quotient. +* `CyclicIntRing.coeff`: coefficients in the unique representative of degree less than `p`. + +## Main results + +* `CyclicIntRing.coeff_x_pow`: powers have a single coefficient at the residue of their + exponent. + +## References + +Adapted from [Anthropic's +Hashing.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/Theorem17/Hashing.lean#L57-L101). +-/ + +@[expose] public section + +namespace Cslib + +open Polynomial in +/-- The quotient ring `ℤ[X]/(X^p - 1)`. -/ +abbrev CyclicIntRing (p : ℕ) : Type := AdjoinRoot ((X : ℤ[X]) ^ p - C 1) + +namespace CyclicIntRing + +open Polynomial + +variable {p : ℕ} + +/-- The element `x` of `ℤ[x]/(x^p − 1)`. -/ +noncomputable def x (p : ℕ) : CyclicIntRing p := AdjoinRoot.root _ + +/-- The coefficient of degree `r` in an element `z` of the quotient, as a linear map. It is the +coefficient of `x^r` in the unique polynomial of degree less than `p` that represents `z` +(Mathlib's `AdjoinRoot.modByMonicHom` gives that polynomial; it asks for the fact that `x^p − 1` +is monic, which holds as `p ≠ 0`). +-/ +noncomputable def coeff (hp : p ≠ 0) (r : ℕ) : CyclicIntRing p →ₗ[ℤ] ℤ := + lcoeff ℤ r ∘ₗ AdjoinRoot.modByMonicHom (monic_X_pow_sub_C (1 : ℤ) hp) + +/-- The coefficient, in terms of the representing polynomial. -/ +theorem coeff_apply (hp : p ≠ 0) (r : ℕ) (z : CyclicIntRing p) : + coeff hp r z = (AdjoinRoot.modByMonicHom (monic_X_pow_sub_C (1 : ℤ) hp) z).coeff r := + rfl + +/-- In `ℤ[x]/(x^p − 1)` we have `x^p = 1`. -/ +private theorem x_pow_self (p : ℕ) : x p ^ p = 1 := by + have h : AdjoinRoot.mk ((X : ℤ[X]) ^ p - C 1) ((X : ℤ[X]) ^ p - C 1) = 0 := AdjoinRoot.mk_self + rw [map_sub, map_pow, AdjoinRoot.mk_X, C_1, map_one] at h + exact sub_eq_zero.1 h + +/-- The coefficient of `x^r` in the power `x^k` of `ℤ[x]/(x^p − 1)` is 1 if `r = k mod p` and 0 +otherwise. +-/ +theorem coeff_x_pow (hp : p ≠ 0) (k r : ℕ) : + coeff hp r (x p ^ k) = if r = k % p then 1 else 0 := by + have hred : x p ^ k = x p ^ (k % p) := by + conv_lhs => rw [← Nat.div_add_mod k p, pow_add, pow_mul, x_pow_self, one_pow, one_mul] + have hdeg : ((X : ℤ[X]) ^ (k % p)).degree < ((X : ℤ[X]) ^ p - C 1).degree := by + rw [degree_X_pow, degree_X_pow_sub_C (Nat.pos_of_ne_zero hp)] + exact_mod_cast Nat.mod_lt k (Nat.pos_of_ne_zero hp) + rw [coeff_apply, hred, x, ← AdjoinRoot.mk_X, ← map_pow, AdjoinRoot.modByMonicHom_mk, + (modByMonic_eq_self_iff (monic_X_pow_sub_C (1 : ℤ) hp)).2 hdeg, coeff_X_pow] + +end CyclicIntRing + + +end Cslib diff --git a/ComplexityTests/ExactTriangleHashing.lean b/ComplexityTests/ExactTriangleHashing.lean new file mode 100644 index 0000000000..41df09a874 --- /dev/null +++ b/ComplexityTests/ExactTriangleHashing.lean @@ -0,0 +1,115 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +import Complexity.Algorithms.ExactTriangle.Hashing + +/-! +# Modular Exact Triangle counting and prime-selection regressions + +Concrete zeros, false positives, coefficient wraparound, minimum selection, the integer window +boundary and an admitted instance of the universal bound. These assert no native running time. +-/ + +namespace ComplexityTests.ExactTriangleHashing + +open Cslib Cslib.Algorithms.ExactTriangle + +private def ab (a _ : Fin 2) : ℤ := if a = 0 then 0 else 2 +private def z (_ _ : Fin 2) : ℤ := 0 + +private theorem window16 : primesInRange 16 = {2, 3} := by + rw [primesInRange_eq_filter] + decide + +-- Both endpoints matter: prime 2 is admitted at the lower endpoint for D=16; +-- prime 5 is excluded at the upper endpoint for D=25, and prime 3 is admitted there. +example : 2 ∈ primesInRange 16 ∧ 5 ∉ primesInRange 25 ∧ 3 ∈ primesInRange 25 := by + simp only [mem_primesInRange_iff_sq] + decide + +example : 2 ∉ primesInRange 17 := by rw [mem_primesInRange_iff_sq]; decide + +example : countZeroMod ab z z 2 = 8 ∧ countZeroMod ab z z 3 = 4 := by + unfold countZeroMod weightSum ab z + decide + +example : zeroCount ab z z = 4 ∧ falsePositiveCount ab z z 2 = 4 ∧ + falsePositiveCount ab z z 3 = 0 := by + unfold zeroCount falsePositiveCount IsFalsePositive weightSum ab z + decide + +example : countZeroMod ab z z 0 = 4 ∧ falsePositiveCount ab z z 0 = 0 := by + unfold countZeroMod falsePositiveCount IsFalsePositive weightSum ab z + decide + +-- Counting distinguishes the two candidate primes; the minimum is 3. +example : IsSelectedPrime ab z z 16 3 ∧ ¬IsSelectedPrime ab z z 16 2 := by + simp only [IsSelectedPrime, window16] + unfold countZeroMod weightSum ab z + decide + +-- A tie permits either minimum prime; the predicate does not impose a tie-break. +example : IsSelectedPrime (fun _ _ : Fin 1 => 0) (fun _ _ => 0) (fun _ _ => 0) 16 2 ∧ + IsSelectedPrime (fun _ _ : Fin 1 => 0) (fun _ _ => 0) (fun _ _ => 0) 16 3 := by + constructor <;> refine ⟨by rw [window16]; decide, ?_⟩ <;> intro q _ <;> + simp [countZeroMod, weightSum] + +example : HasZeroTriangle ab z z := ⟨0, 0, 0, rfl⟩ +example : 0 < zeroCount ab z z := (zeroCount_pos_iff ab z z).mpr ⟨0, 0, 0, rfl⟩ + +private def one (_ _ : Fin 2) : ℤ := 1 + +example : ¬HasZeroTriangle one one one := by + rintro ⟨a, b, c, h⟩ + norm_num [one] at h + +example : zeroCount one one one = 0 ∧ countZeroMod one one one 3 = 8 ∧ + falsePositiveCount one one one 3 = 8 := by + unfold zeroCount countZeroMod falsePositiveCount IsFalsePositive weightSum one + decide + +-- A negative residue and a positive residue wrap to exponent zero modulo 3. +example : CyclicIntRing.coeff (by decide : 3 ≠ 0) 0 + ((matP (fun _ _ : Fin 2 => (-1 : ℤ)) 3 * + matQ (fun _ _ : Fin 2 => (1 : ℤ)) 3) 0 0) = 2 := by + rw [coeff_matP_mul_matQ (fun _ _ : Fin 2 => (1 : ℤ)) + (fun _ _ : Fin 2 => (-1 : ℤ)) (by decide) 0 0 0 (by decide)] + decide + +-- These matrices are asymmetric: transposing either input changes the selected count. +private def ac (a c : Fin 2) : ℤ := if a = 0 then if c = 0 then -1 else 2 else 0 +private def bc (b c : Fin 2) : ℤ := if b = 0 ∧ c = 1 then 1 else 0 + +example : CyclicIntRing.coeff (by decide : 3 ≠ 0) 0 ((matP ac 3 * matQ bc 3) 0 0) = 1 := by + rw [coeff_matP_mul_matQ bc ac (by decide) 0 0 0 (by decide)] + unfold ac bc + decide + +example : CyclicIntRing.coeff (by decide : 2 ≠ 0) 1 (CyclicIntRing.x 2 ^ 5) = 1 := by + rw [CyclicIntRing.coeff_x_pow] + norm_num + +example : CyclicIntRing.coeff (by decide : 2 ≠ 0) 3 (CyclicIntRing.x 2 ^ 5) = 0 := by + rw [CyclicIntRing.coeff_x_pow] + norm_num + +example : zeroCount (fun _ _ : Fin 0 => 0) (fun _ _ => 0) (fun _ _ => 0) = 0 := by + simp [zeroCount] + +example : ¬WeightsBoundedBy (fun _ _ : Fin 1 => (0 : ℤ)) (fun _ _ => 0) (fun _ _ => 0) (-1) := by + norm_num [WeightsBoundedBy] + +-- The full bound's n,D,kappa hypotheses and input bound have a concrete witness. +example : (falsePositiveCount (fun _ _ : Fin 16 => 0) (fun _ _ => 0) (fun _ _ => 0) 2 : ℝ) ≤ + Hashing.falsePositiveConst * (1 * (16 : ℝ) ^ 3 * Real.log 16 / Real.sqrt 16) := by + apply Hashing.falsePositiveCount_le_const_mul (by decide : 16 ≤ 16) (by decide : 16 ≤ 16) + (by norm_num : (1 : ℝ) ≤ 1) + · norm_num [WeightsBoundedBy] + · refine ⟨by rw [window16]; decide, ?_⟩ + intro q _ + simp [countZeroMod, weightSum] + +end ComplexityTests.ExactTriangleHashing diff --git a/NOTICE b/NOTICE index 3ab01a2df3..95b9c8d57a 100644 --- a/NOTICE +++ b/NOTICE @@ -1,10 +1,11 @@ The 3SUM/Exact Triangle reductions, their finite counting lemmas and native matrix builder in -Complexity/Algorithms/ExactTriangle/{Defs,Lopsided}.lean, +Complexity/Algorithms/ExactTriangle/{Defs,Weights,Lopsided,Hashing}.lean, Complexity/Algorithms/ThreeSum/, Complexity/Algorithms/WordRAM/ThreeSum/, Complexity/Foundations/Data/Finset/Average.lean, Complexity/Foundations/Data/Matrix/IntegerBits.lean and -Complexity/Foundations/NumberTheory/PrimeCounting.lean are adapted from +Complexity/Foundations/NumberTheory/{PrimeCounting,PrimeWindow}.lean and +Complexity/Foundations/RingTheory/CyclicIntRing.lean are adapted from anthropics/formal-math, 3sum-apsp, revision 16a5ce41304fcdaea4e0be7c132e1c361dc26da5. @@ -13,3 +14,7 @@ Sparse Lopsided Graphs" Copyright 2026 Anthropic, PBC This product is licensed under the Apache License, Version 2.0 (see LICENSE). + +The prime-window proof also adapts Mathlib/NumberTheory/Bertrand.lean, +Copyright (c) 2020 Patrick Stevens; authors Patrick Stevens, Bolton Bailey, +under the Apache License, Version 2.0. diff --git a/ROADMAP.md b/ROADMAP.md index 5b3b230013..d31b9ede12 100644 --- a/ROADMAP.md +++ b/ROADMAP.md @@ -422,8 +422,9 @@ each statement against the numbered result. Organize modules by subject, not cha has the Convolution-3SUM/Exact Triangle correctness reduction, a charged native single-instance matrix builder and the complete finite mathematical Chan–He reduction ([#1435]), including ordinary distinct-position inputs and explicit size bounds. The finite lopsided matrix/counting - interface ([#1438]) adds padding, detection and signed digit reconstruction; it still needs - modular Exact Triangle hashing/false-positive splitting, tensor/sparse-query programs, native + interface ([#1438]) adds padding, detection and signed digit reconstruction. Modular counting + and the selected-prime false-positive bound ([#1442]) are proved mathematically; it still needs + residue/chunk/piece splitting and witness recovery, tensor/sparse-query programs, native Chan–He preprocessing, the reduction host, the fast Exact Triangle solver and full input/time/width/ asymptotic composition, then the improved Section 4 route to exponent 1.9992; size-based width selection with an input-fit guarantee and explicit overhead; expected-time rules for unbounded @@ -550,4 +551,5 @@ not planned, not as proved. [#1394]: https://github.com/cslib-community/complexitylib/issues/1394 [#1435]: https://github.com/cslib-community/complexitylib/issues/1435 [#1438]: https://github.com/cslib-community/complexitylib/issues/1438 +[#1442]: https://github.com/cslib-community/complexitylib/issues/1442 [#1449]: https://github.com/cslib-community/complexitylib/issues/1449 diff --git a/UPSTREAM.md b/UPSTREAM.md index 87d9cb5288..2bc29e783b 100644 --- a/UPSTREAM.md +++ b/UPSTREAM.md @@ -174,8 +174,16 @@ Do not add a ledger row for an ordinary candidate or repeat its source marker he `Set.EqOn`; no additional weighted-triangle predicate or source machine is introduced. This is the mathematical interface milestone [#1438](https://github.com/cslib-community/complexitylib/issues/1438), with no sparsity assumption, - preprocessing/time claim or solver. Modular hashing/false-positive splitting and tensor/sparse - query programs remain in #1394, alongside the native and asymptotic composition below. + preprocessing/time claim or solver. `Algorithms/ExactTriangle/{Weights,Hashing}.lean` now adapts + the modular counting/coefficient identities, minimum-count prime and uniform false-positive + bound from `Sec3/Theorem17/Hashing.lean`, using supplied matrices and the existing decision + predicate. Its integer-square window characterisation supports later native enumeration. + The quotient/coefficient API lives in `Foundations/RingTheory/CyclicIntRing.lean`; the uniform + all-window prime-count bound lives in `Foundations/NumberTheory/PrimeWindow.lean`, retaining + the source credit to Mathlib, Patrick Stevens and Bolton Bailey. These are mathematical + statements, not native instructions or an assumed solver. See + [#1442](https://github.com/cslib-community/complexitylib/issues/1442). Residue/chunk/piece splitting, + witness recovery and tensor/sparse query programs remain in #1394, alongside native composition. The existing signed and unsigned APIs were reused; no source machine or separate run semantics were copied. The full reduction host, native Chan–He preprocessing, fast Exact Triangle and full input/time/width/asymptotic composition remain in From 0d8da3e78edceaa2e11fa3c31a1b44fce10dbdbd Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 06:09:24 -0400 Subject: [PATCH 03/11] doc(threeSum): cite the lopsided triangle paper and clarify adaptation scope --- Complexity/Algorithms/ExactTriangle/Hashing.lean | 1 + Complexity/Algorithms/ExactTriangle/Lopsided.lean | 1 + UPSTREAM.md | 13 +++++++------ references.bib | 12 ++++++++++++ 4 files changed, 21 insertions(+), 6 deletions(-) diff --git a/Complexity/Algorithms/ExactTriangle/Hashing.lean b/Complexity/Algorithms/ExactTriangle/Hashing.lean index 6799e79fc7..affaadd2e2 100644 --- a/Complexity/Algorithms/ExactTriangle/Hashing.lean +++ b/Complexity/Algorithms/ExactTriangle/Hashing.lean @@ -39,6 +39,7 @@ asserted. The supplied three matrices retain the existing `HasZeroTriangle` sema ## References +The modular reduction is part of [AlmanWilliams2026, Section 3]. Direct adaptation of [Anthropic's Hashing.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/Theorem17/Hashing.lean), with the minimal count/selection definitions extracted from its `PaperStatements.lean`. diff --git a/Complexity/Algorithms/ExactTriangle/Lopsided.lean b/Complexity/Algorithms/ExactTriangle/Lopsided.lean index d52e779437..84694440fe 100644 --- a/Complexity/Algorithms/ExactTriangle/Lopsided.lean +++ b/Complexity/Algorithms/ExactTriangle/Lopsided.lean @@ -36,6 +36,7 @@ fast solver or subquadratic running time. ## References +The lopsided matrix interface is used in [AlmanWilliams2026, Section 3]. Adapted from [Anthropic's MatrixLanguage.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/MatrixLanguage.lean) and the lopsided definitions in [PaperStatements.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/PaperStatements.lean#L798-L852). The adaptation uses `Set.EqOn`, pointwise signed parts, and direct finite cardinalities. diff --git a/UPSTREAM.md b/UPSTREAM.md index 2bc29e783b..8b49547dae 100644 --- a/UPSTREAM.md +++ b/UPSTREAM.md @@ -158,6 +158,13 @@ Do not add a ledger row for an ordinary candidate or repeat its source marker he Given its arguments in registers, it writes one instance in at most `131*t² + 5` instructions at width `w + 2`, excluding the halt; its theorem identifies the three encoded matrices and preserves the input. Layout, signed index/filler bounds and target address bounds are explicit. + The existing signed and unsigned APIs were reused; no source machine or separate run semantics + were copied. The full reduction host, native Chan–He preprocessing, fast Exact Triangle and full + input/time/width/asymptotic composition remain in + [#1394](https://github.com/cslib-community/complexitylib/issues/1394), including the improved + Section 4 route needed for the published exponent 1.9992. Anthropic's copyright and attribution + are retained; the native proofs and representation choices are local. + `Algorithms/ThreeSum/ThreeSets.lean` and `ArrayPacking.lean` provide the generic distinct-position frontend and three-array packing. `Algorithms/ThreeSum/ChanHe/` now adapts the core mathematical reduction: minimum-prime searches, collision/heavy-element bounds, light-element arrays, bounded @@ -184,12 +191,6 @@ Do not add a ledger row for an ordinary candidate or repeat its source marker he statements, not native instructions or an assumed solver. See [#1442](https://github.com/cslib-community/complexitylib/issues/1442). Residue/chunk/piece splitting, witness recovery and tensor/sparse query programs remain in #1394, alongside native composition. - The existing signed and unsigned APIs were reused; no source machine or separate run semantics - were copied. The full reduction host, native Chan–He preprocessing, fast Exact Triangle and full - input/time/width/asymptotic composition remain in - [#1394](https://github.com/cslib-community/complexitylib/issues/1394), including the improved - Section 4 route needed for the published exponent 1.9992. Anthropic's copyright and attribution - are retained; the native proofs and representation choices are local. * **Typed RAM upstream core.** [CSLib #1158](https://github.com/leanprover/cslib/pull/1158) at `beea60a15eddc91d93e92fa5af9fb965c351d0d5` proposes `Cslib/Computability/Machines/RAM/Basic.lean`, in namespace `Cslib.RAM.Typed`. diff --git a/references.bib b/references.bib index 8cf6daceb4..f42e7901a7 100644 --- a/references.bib +++ b/references.bib @@ -103,6 +103,18 @@ @misc{Algolean2026 url = {https://github.com/Shreyas4991/Algolean} } +@misc{AlmanWilliams2026, + author = {Alman, Josh and Vassilevska Williams, Virginia}, + title = {Truly Subquadratic {3SUM} and Truly Subcubic {APSP} via Triangles in + Sparse Lopsided Graphs}, + year = {2026}, + eprint = {2610.06783}, + archivePrefix = {arXiv}, + primaryClass = {cs.DS}, + doi = {10.48550/arXiv.2610.06783}, + url = {https://arxiv.org/abs/2610.06783v1} +} + @article{AllenderBuhrmanKouckyVanMelkebeekRonneburger2006, author = {Allender, Eric and Buhrman, Harry and Kouck{\'y}, Michal and van Melkebeek, Dieter and Ronneburger, Detlef}, title = {Power from Random Strings}, From 8d80ecce623c8ab4a08742ed6944892ebe9f9b21 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 06:09:49 -0400 Subject: [PATCH 04/11] doc(bibliography): sort the new paper entry --- references.bib | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/references.bib b/references.bib index f42e7901a7..af7ca6f2c2 100644 --- a/references.bib +++ b/references.bib @@ -103,18 +103,6 @@ @misc{Algolean2026 url = {https://github.com/Shreyas4991/Algolean} } -@misc{AlmanWilliams2026, - author = {Alman, Josh and Vassilevska Williams, Virginia}, - title = {Truly Subquadratic {3SUM} and Truly Subcubic {APSP} via Triangles in - Sparse Lopsided Graphs}, - year = {2026}, - eprint = {2610.06783}, - archivePrefix = {arXiv}, - primaryClass = {cs.DS}, - doi = {10.48550/arXiv.2610.06783}, - url = {https://arxiv.org/abs/2610.06783v1} -} - @article{AllenderBuhrmanKouckyVanMelkebeekRonneburger2006, author = {Allender, Eric and Buhrman, Harry and Kouck{\'y}, Michal and van Melkebeek, Dieter and Ronneburger, Detlef}, title = {Power from Random Strings}, @@ -137,6 +125,18 @@ @article{AllenderRubinstein1988 doi = {10.1137/0217075} } +@misc{AlmanWilliams2026, + author = {Alman, Josh and Vassilevska Williams, Virginia}, + title = {Truly Subquadratic {3SUM} and Truly Subcubic {APSP} via Triangles in + Sparse Lopsided Graphs}, + year = {2026}, + eprint = {2610.06783}, + archivePrefix = {arXiv}, + primaryClass = {cs.DS}, + doi = {10.48550/arXiv.2610.06783}, + url = {https://arxiv.org/abs/2610.06783v1} +} + @article{AlonBoppana1987, author = {Alon, Noga and Boppana, Ravi B.}, title = {The Monotone Circuit Complexity of {B}oolean Functions}, From fd69be7bc459d83cab3863f6f6e2f59239b21595 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 09:02:59 -0400 Subject: [PATCH 05/11] fix(exactTriangle): repair the prime-window header and hashing documentation --- .../Algorithms/ExactTriangle/Hashing.lean | 115 +++++++++++------- .../Foundations/NumberTheory/PrimeWindow.lean | 4 +- .../Foundations/RingTheory/CyclicIntRing.lean | 1 - 3 files changed, 71 insertions(+), 49 deletions(-) diff --git a/Complexity/Algorithms/ExactTriangle/Hashing.lean b/Complexity/Algorithms/ExactTriangle/Hashing.lean index affaadd2e2..fbb2834b6e 100644 --- a/Complexity/Algorithms/ExactTriangle/Hashing.lean +++ b/Complexity/Algorithms/ExactTriangle/Hashing.lean @@ -43,7 +43,8 @@ The modular reduction is part of [AlmanWilliams2026, Section 3]. Direct adaptation of [Anthropic's Hashing.lean](https://github.com/anthropics/formal-math/blob/16a5ce41304fcdaea4e0be7c132e1c361dc26da5/3sum-apsp/ThreeSumApsp/Sec3/Theorem17/Hashing.lean), with the minimal count/selection definitions extracted from its `PaperStatements.lean`. -Source-derived prose below is paraphrased. The source pin is +The declaration docstrings describe the Lean definitions and statements in local notation. +The source pin is `anthropics/formal-math@16a5ce41304fcdaea4e0be7c132e1c361dc26da5`. -/ @@ -52,16 +53,20 @@ Source-derived prose below is paraphrased. The source pin is namespace Cslib.Algorithms.ExactTriangle open Classical in -/-- Primes in the half-open window `[sqrt D / 2, sqrt D)`. -/ +/-- Prime numbers in the half-open real interval `[sqrt D / 2, sqrt D)`. +-/ noncomputable def primesInRange (D : ℕ) : Finset ℕ := (Finset.range D).filter fun p => p.Prime ∧ Real.sqrt D / 2 ≤ (p : ℝ) ∧ (p : ℝ) < Real.sqrt D -/-- Number of triangles whose sum vanishes modulo `p`. -/ +/-- Count vertex triples whose three edge weights have sum congruent to zero modulo `p`. +-/ def countZeroMod {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) (p : ℕ) : ℕ := (Finset.univ.filter fun t : Fin n × Fin n × Fin n => weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≡ 0 [ZMOD (p : ℤ)]).card -/-- A window prime minimises the modular count over all window primes. -/ +/-- A prime in `primesInRange D` whose modular-zero count is minimal among all primes in that +window. Tied minima all satisfy this predicate. +-/ def IsSelectedPrime {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) (D p : ℕ) : Prop := p ∈ primesInRange D ∧ ∀ q ∈ primesInRange D, countZeroMod wAB wBC wAC p ≤ countZeroMod wAB wBC wAC q @@ -81,8 +86,8 @@ private theorem toNat_emod_lt {p : ℕ} (hp : 0 < p) (u : ℤ) : (u % (p : ℤ)) /-! ### The primes in the range -/ -/-- Proof of Theorem 17, a prime p ∈ [√D/2, √D): the primes in the range are the primes `p` with -`√D/2 ≤ p < √D`. (The bound `p < D` in the definition excludes none of them.) +/-- Window membership is equivalent to primality and the two real endpoint inequalities. The +finite enumeration bound in `primesInRange` excludes no member of the interval. -/ theorem mem_primesInRange {D p : ℕ} : p ∈ primesInRange D ↔ p.Prime ∧ Real.sqrt D / 2 ≤ (p : ℝ) ∧ (p : ℝ) < Real.sqrt D := by @@ -93,7 +98,9 @@ theorem mem_primesInRange {D p : ℕ} : have hsq : p ^ 2 < D := by exact_mod_cast (Real.lt_sqrt (Nat.cast_nonneg p)).1 h.2.2 exact (Nat.le_self_pow two_ne_zero p).trans_lt hsq -/-- Window membership is an exact integer-square test, suitable for native enumeration. -/ +/-- Characterize window membership using primality and integer-square inequalities, including the +closed lower and open upper endpoints. +-/ theorem mem_primesInRange_iff_sq {D p : ℕ} : p ∈ primesInRange D ↔ p.Prime ∧ D ≤ 4 * p ^ 2 ∧ p ^ 2 < D := by rw [mem_primesInRange] @@ -111,7 +118,8 @@ theorem mem_primesInRange_iff_sq {D p : ℕ} : have h : (D : ℝ) ≤ 4 * (p : ℝ) ^ 2 := by exact_mod_cast hlo nlinarith -/-- A computable finite enumeration of exactly the same half-open prime window. -/ +/-- Enumerate the real prime window using a decidable filter of the integers below `D`. +-/ theorem primesInRange_eq_filter (D : ℕ) : primesInRange D = (Finset.range D).filter fun p => p.Prime ∧ D ≤ 4 * p ^ 2 ∧ p ^ 2 < D := by ext p @@ -137,20 +145,26 @@ section Counting variable {n : ℕ} (wAB wBC wAC : Fin n → Fin n → ℤ) -/-- A false positive of p is a triple (a,b,c) with S(a,b,c) ≠ 0 and p ∣ S(a,b,c). -/ +/-- A vertex triple has a nonzero integer weight sum that is congruent to zero modulo `p`. +-/ abbrev IsFalsePositive (p : ℕ) (t : Fin n × Fin n × Fin n) : Prop := weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≠ 0 ∧ (p : ℤ) ∣ weightSum wAB wBC wAC t.1 t.2.1 t.2.2 -/-- let falsePositiveCount(p) denote the number of false positives of p. -/ +/-- Count vertex triples satisfying `IsFalsePositive` for the supplied three weight matrices and +modulus. +-/ def falsePositiveCount (p : ℕ) : ℕ := (Finset.univ.filter fun t : Fin n × Fin n × Fin n => IsFalsePositive wAB wBC wAC p t).card -/-- zeroCount, the number of zero triangles. -/ +/-- Count vertex triples whose three supplied edge weights sum to zero as an integer. +-/ def zeroCount : ℕ := (Finset.univ.filter fun t : Fin n × Fin n × Fin n => weightSum wAB wBC wAC t.1 t.2.1 t.2.2 = 0).card -/-- A positive zero count is exactly the existing Exact Triangle decision predicate. -/ +/-- The zero-triple count is positive exactly when the supplied matrices satisfy +`HasZeroTriangle`. +-/ theorem zeroCount_pos_iff : 0 < zeroCount wAB wBC wAC ↔ HasZeroTriangle wAB wBC wAC := by classical simp only [zeroCount, Finset.card_pos, HasZeroTriangle] @@ -160,8 +174,8 @@ theorem zeroCount_pos_iff : 0 < zeroCount wAB wBC wAC ↔ HasZeroTriangle wAB wB · rintro ⟨a, b, c, h⟩ exact ⟨(a, b, c), Finset.mem_filter.mpr ⟨Finset.mem_univ _, h⟩⟩ -/-- This number is the number of false positives falsePositiveCount(p) plus zeroCount, the number of -zero triangles, which does not depend on p. +/-- Partition the modular-zero triples into false positives and exact zero-sum triples. The +exact-zero count is independent of the modulus. -/ theorem countZeroMod_eq (p : ℕ) : (countZeroMod wAB wBC wAC) p = (falsePositiveCount wAB wBC wAC) p + (zeroCount wAB wBC wAC) := by @@ -173,17 +187,21 @@ theorem countZeroMod_eq (p : ℕ) : (countZeroMod wAB wBC wAC) p = (falsePositiv by_cases h0 : weightSum wAB wBC wAC t.1 t.2.1 t.2.2 = 0 <;> simp [IsFalsePositive, Int.modEq_zero_iff_dvd, h0] -/-- let P[a,c] := x^{w(a,c) mod p}, a matrix over `ℤ[x]/(x^p − 1)`. `w mod p` is the residue in -`{0, …, p − 1}`. +/-- Encode the A–C weight at `(a,c)` as a power of the cyclic polynomial generator, with exponent +equal to its nonnegative residue modulo `p`. -/ noncomputable def matP (p : ℕ) : Matrix (Fin n) (Fin n) (CyclicIntRing p) := Matrix.of fun a c => CyclicIntRing.x p ^ (wAC a c % (p : ℤ)).toNat -/-- and Q[c,b] := x^{w(b,c) mod p}. -/ +/-- Encode the B–C weight at `(b,c)` in matrix position `(c,b)`. This orientation makes the +product sum over the middle index `c`. +-/ noncomputable def matQ (p : ℕ) : Matrix (Fin n) (Fin n) (CyclicIntRing p) := Matrix.of fun c b => CyclicIntRing.x p ^ (wBC b c % (p : ℤ)).toNat -/-- the coefficient of x^r in (PQ)[a,b] is the number of c ∈ C with w(a,c) + w(b,c) ≡ r (mod p). -/ +/-- For a nonzero modulus and a residue `r < p`, the coefficient of degree `r` in the cyclic +matrix product counts middle vertices whose two incident weights sum to `r` modulo `p`. +-/ theorem coeff_matP_mul_matQ {p : ℕ} (hp : p ≠ 0) (a b : Fin n) (r : ℕ) (hr : r < p) : CyclicIntRing.coeff hp r ((matP wAC p * matQ wBC p) a b) = ((Finset.univ.filter fun c : Fin n => @@ -200,8 +218,8 @@ theorem coeff_matP_mul_matQ {p : ℕ} (hp : p ≠ 0) (a b : Fin n) (r : ℕ) (hr refine Finset.sum_congr rfl fun c _ => ?_ simp only [eq_add_toNat_emod_iff hp hr] -/-- then falsePositiveCount(p) + zeroCount is the sum over the pairs (a,b) ∈ A × B of the -coefficient of x^{−w(a,b) mod p} in (PQ)[a,b]. +/-- Recover the total modular-zero count from the cyclic matrix product by summing, over `(a,b)`, +the coefficient at the residue of `-wAB a b`. -/ theorem falsePositiveCount_add_zeroCount_eq_sum_coeff {p : ℕ} (hp : p ≠ 0) : (((falsePositiveCount wAB wBC wAC) p + (zeroCount wAB wBC wAC) : ℕ) : ℤ) = @@ -235,7 +253,8 @@ theorem falsePositiveCount_add_zeroCount_eq_sum_coeff {p : ℕ} (hp : p ≠ 0) : end Counting -/-- the fewer than √D primes in the range. -/ +/-- For `D ≥ 16`, the number of primes in `primesInRange D` is strictly less than `sqrt D`. +-/ theorem card_primesInRange_lt (D : ℕ) (hD : 16 ≤ D) : ((primesInRange D).card : ℝ) < Real.sqrt D := by have hpos : 0 < Real.sqrt D := by linarith [four_le_sqrt hD] @@ -259,15 +278,17 @@ section Selection variable {n D p : ℕ} {κ : ℝ} (wAB wBC wAC : Fin n → Fin n → ℤ) -/-- There is a prime to select. -/ +/-- For every `D ≥ 16` and every three supplied weight matrices, the prime window contains a +minimum of the modular-zero count. +-/ theorem exists_isSelectedPrime (D : ℕ) (hD : 16 ≤ D) : ∃ p, IsSelectedPrime wAB wBC wAC D p := by obtain ⟨q, hq⟩ := Nat.exists_prime_half_le_and_lt (Real.sqrt D) (four_le_sqrt hD) exact Finset.exists_min_image (primesInRange D) (countZeroMod wAB wBC wAC) ⟨q, mem_primesInRange.2 hq⟩ -/-- We select the prime with the smallest count, which is also the prime with the fewest false -positives. +/-- A selected minimum of the modular-zero count also minimizes the false-positive count over the +same prime window. -/ theorem IsSelectedPrime.falsePositiveCount_le {wAB wBC wAC : Fin n → Fin n → ℤ} (hp : IsSelectedPrime wAB wBC wAC D p) : @@ -282,9 +303,8 @@ theorem IsSelectedPrime.falsePositiveCount_le {wAB wBC wAC : Fin n → Fin n → ### The bound on the number of false positives of the selected prime -/ open Classical in -/-- a triple with S(a,b,c) ≠ 0 is a false positive exactly of the primes in the range that divide -S(a,b,c). Since 0 < |S(a,b,c)| ≤ 3n^ν and these primes are at least √D/2, there are at most -log_{√D/2}(3n^ν) of them. +/-- For a fixed nonzero-sum triple and edge weights bounded by `n^κ`, bound the number of window +primes accepting that triple by the base-`sqrt D / 2` logarithm of `3*n^κ`. -/ theorem card_falsePositive_primes_le (hD : 16 ≤ D) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) (t : Fin n × Fin n × Fin n) (ht : weightSum wAB wBC wAC t.1 t.2.1 t.2.2 ≠ 0) : @@ -316,8 +336,9 @@ theorem card_falsePositive_primes_le (hD : 16 ≤ D) (hT : WeightsBoundedBy wAB rw [Real.le_logb_iff_rpow_le (by linarith) hY, Real.rpow_natCast] exact hpow.trans hprod -/-- Hence the numbers of false positives of all the primes in the range add up to at most n³ -log_{√D/2}(3n^ν). +/-- Sum the per-triple prime-divisibility bound over all `n^3` triples. For `D ≥ 16`, `κ ≥ 1` and +edge weights bounded by `n^κ`, the sum of false-positive counts is at most `n^3 * logb (sqrt D / +2) (3*n^κ)`. -/ theorem sum_falsePositiveCount_le (hD : 16 ≤ D) (hκ : 1 ≤ κ) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) : @@ -356,11 +377,9 @@ theorem sum_falsePositiveCount_le (hD : 16 ≤ D) (hκ : 1 ≤ κ) (hT : Weights end Selection -/-- By the prime number theorem there are Ω(√D/log D) primes in the range. - -NOTE. The modular reduction needs this for every `D ≥ 16`, not only for large `D`, and so it is -stated. It follows from a counting form of Bertrand's postulate; the prime number theorem is -not needed. +/-- There is a universal positive `c` such that every window with `D ≥ 16` contains at least `c * +sqrt D / log D` primes. The proof uses the all-window counting form of Bertrand’s postulate from +`PrimeWindow`, with no prime-number-theorem assumption. -/ theorem exists_le_card_primesInRange : ∃ c : ℝ, 0 < c ∧ ∀ D : ℕ, 16 ≤ D → @@ -369,7 +388,9 @@ theorem exists_le_card_primesInRange : exact ⟨c, hc, fun D hD => hcount D hD _ fun _ hp hge hlt => mem_primesInRange.2 ⟨hp, hge, hlt⟩⟩ -/-- falsePositiveCount(p) is at most the average over them, multiplied out. -/ +/-- Multiply the selected prime’s minimum false-positive count by the number of window primes to +obtain a bound by the sum of their false-positive counts. +-/ theorem IsSelectedPrime.falsePositiveCount_mul_card_le_sum {n D p : ℕ} {wAB wBC wAC : Fin n → Fin n → ℤ} (hp : IsSelectedPrime wAB wBC wAC D p) : @@ -378,7 +399,8 @@ theorem IsSelectedPrime.falsePositiveCount_mul_card_le_sum {n D p : ℕ} {wAB wB rw [mul_comm, ← smul_eq_mul] exact Finset.card_nsmul_le_sum _ _ _ hp.falsePositiveCount_le -/-- log D = O(log(√D/2)) for D ≥ 16. -/ +/-- For `D ≥ 16`, compare `log D` with four times `log (sqrt D / 2)`. +-/ theorem log_le_four_mul_log_sqrt_div_two (D : ℕ) (hD : 16 ≤ D) : Real.log D ≤ 4 * Real.log (Real.sqrt D / 2) := by have hD16 : (16 : ℝ) ≤ D := by exact_mod_cast hD @@ -391,9 +413,9 @@ theorem log_le_four_mul_log_sqrt_div_two (D : ℕ) (hD : 16 ≤ D) : rw [Real.log_div hsqrt.ne' (by norm_num), Real.log_sqrt (by linarith)] linarith [hlog16] -/-- falsePositiveCount(p) is at most the average over them, so [...] falsePositiveCount(p) = O(n³ -log(3n^ν)/√D): if there are at least `c √D/log D` primes in the range, then -`falsePositiveCount(p) ≤ (4/c) n³ log(3n^κ)/√D`. +/-- For `16 ≤ D ≤ n`, `κ ≥ 1` and edge weights bounded by `n^κ`, a positive pointwise prime-count +lower bound `c * sqrt D / log D` gives every selected prime at most `(4/c) * n^3 * log (3*n^κ) / +sqrt D` false positives. -/ theorem falsePositiveCount_le_of_le_card_primesInRange {n D p : ℕ} {κ c : ℝ} (wAB wBC wAC : Fin n → Fin n → ℤ) (hc : 0 < c) @@ -435,7 +457,9 @@ theorem falsePositiveCount_le_of_le_card_primesInRange {n D p : ℕ} {κ c : ℝ _ ≤ (n : ℝ) ^ 3 * (4 * Real.log Y / Real.log D) * Real.log D := by gcongr _ = 4 * (n : ℝ) ^ 3 * Real.log Y := by field_simp -/-- falsePositiveCount(p) = O(n³ log(3n^ν)/√D) = O(ν n³ log n/√D). -/ +/-- There is a universal `C ≥ 1` such that, for `16 ≤ D ≤ n`, `κ ≥ 1` and edge weights bounded by +`n^κ`, every selected prime has at most `C * κ * n^3 * log n / sqrt D` false positives. +-/ theorem exists_falsePositiveCount_le : ∃ C : ℝ, 1 ≤ C ∧ ∀ {n D p : ℕ} {κ : ℝ}, 16 ≤ D → D ≤ n → 1 ≤ κ → ∀ wAB wBC wAC : Fin n → Fin n → ℤ, @@ -461,17 +485,18 @@ theorem exists_falsePositiveCount_le : gcongr exact le_max_right _ _ -/-- A name for the constant in falsePositiveCount(p) = [...] = O(ν n³ log n/√D), for the bounds on -running times that are built on this one. +/-- A universal real constant witnessing `exists_falsePositiveCount_le`. It depends on none of the +input dimensions, weight matrices, exponent or selected prime, and is at least one. -/ noncomputable def Hashing.falsePositiveConst : ℝ := Classical.choose exists_falsePositiveCount_le -/-- The constant is at least 1. -/ +/-- The universal false-positive constant `Hashing.falsePositiveConst` is at least one. +-/ theorem Hashing.one_le_falsePositiveConst : 1 ≤ Hashing.falsePositiveConst := (Classical.choose_spec exists_falsePositiveCount_le).1 -/-- falsePositiveCount(p) = [...] = O(ν n³ log n/√D), with `Hashing.falsePositiveConst` as the -constant. +/-- For `16 ≤ D ≤ n`, `κ ≥ 1` and edge weights bounded by `n^κ`, bound every selected prime’s +false-positive count by `Hashing.falsePositiveConst * κ * n^3 * log n / sqrt D`. -/ theorem Hashing.falsePositiveCount_le_const_mul {n D p : ℕ} {κ : ℝ} (hD : 16 ≤ D) (hDn : D ≤ n) (hκ : 1 ≤ κ) (wAB wBC wAC : Fin n → Fin n → ℤ) (hT : WeightsBoundedBy wAB wBC wAC ((n : ℝ) ^ κ)) diff --git a/Complexity/Foundations/NumberTheory/PrimeWindow.lean b/Complexity/Foundations/NumberTheory/PrimeWindow.lean index dbeb6390e2..9693a53028 100644 --- a/Complexity/Foundations/NumberTheory/PrimeWindow.lean +++ b/Complexity/Foundations/NumberTheory/PrimeWindow.lean @@ -1,6 +1,5 @@ /- -Copyright (c) 2020 Patrick Stevens. All rights reserved. -Copyright (c) 2026 Anthropic, PBC. All rights reserved. +Copyright (c) 2020 Patrick Stevens, 2026 Anthropic, PBC. All rights reserved. Released under Apache 2.0 license as described in the file LICENSE. Authors: Patrick Stevens, Bolton Bailey, Anthropic, PBC -/ @@ -49,7 +48,6 @@ namespace Nat private theorem one_le_log_of_three_le {x : ℝ} (hx : 3 ≤ x) : 1 ≤ Real.log x := (Real.le_log_iff_exp_le (by linarith)).2 (Real.exp_one_lt_three.le.trans hx) - /-- The primes `p` with `n < p ≤ 2n`: the primes that Bertrand's postulate is about. -/ private def bertrandPrimes (n : ℕ) : Finset ℕ := (Ioc n (2 * n)).filter Nat.Prime diff --git a/Complexity/Foundations/RingTheory/CyclicIntRing.lean b/Complexity/Foundations/RingTheory/CyclicIntRing.lean index 94ae6d738f..d15f95b193 100644 --- a/Complexity/Foundations/RingTheory/CyclicIntRing.lean +++ b/Complexity/Foundations/RingTheory/CyclicIntRing.lean @@ -80,5 +80,4 @@ theorem coeff_x_pow (hp : p ≠ 0) (k r : ℕ) : end CyclicIntRing - end Cslib From 4fc6a70f0da95833bedf57afcacd117b8700c22c Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 07:59:17 -0400 Subject: [PATCH 06/11] feat(wordram): add native radix sorting and multiplicity arrays MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Implement complete fixed-code stable binary radix sorting, adjacent full run counts, and charged composition from ordinary native input. Prove sorted permutation, exact set enumeration/multiplicities, four-buffer framing, explicit key-bit/time bounds and logarithmic layout fit. The 81-instruction program uses 16 registers. Scope: #1452 under #1446/#1394; full Chan–He preprocessing remains open. No Histogram or ThreeSum/ChanHe changes. Original implementation/proofs; no source machine import. Curated abstracts and rendering: no impact. Independent review remains requested. Validation on the final staged sources above 8e182b091: all six scripts/check.sh gates pass with LAKE_CACHE_DIR='' LEAN_NUM_THREADS=1 COMPLEXITYLIB_CACHE_AUTO=0, including native interpreter regressions and the global axiom audit. PYTHONUTF8=1 hygiene --check-simp, changed-file lints --base 8e182b091 and git diff --cached --check pass. Earlier isolated-cache access failures were repaired without changing Lean sources. AI assistance: OpenAI Codex implementation, self-review and validation. --- Complexity.lean | 9 + .../Algorithms/WordRAM/Deduplicate/Count.lean | 148 +++++++++ .../Algorithms/WordRAM/Deduplicate/Model.lean | 146 +++++++++ .../WordRAM/Deduplicate/Program.lean | 284 ++++++++++++++++++ Complexity/Algorithms/WordRAM/RadixSort.lean | 106 +++++++ .../WordRAM/RadixSort/Deduplicate.lean | 135 +++++++++ .../Algorithms/WordRAM/RadixSort/Model.lean | 134 +++++++++ .../Algorithms/WordRAM/RadixSort/Program.lean | 204 +++++++++++++ .../Algorithms/WordRAM/RadixSort/Scan.lean | 254 ++++++++++++++++ .../Computability/WordRAM/ListArray.lean | 94 ++++++ ComplexityTests/WordRAMRadixSort.lean | 67 +++++ 11 files changed, 1581 insertions(+) create mode 100755 Complexity/Algorithms/WordRAM/Deduplicate/Count.lean create mode 100755 Complexity/Algorithms/WordRAM/Deduplicate/Model.lean create mode 100755 Complexity/Algorithms/WordRAM/Deduplicate/Program.lean create mode 100755 Complexity/Algorithms/WordRAM/RadixSort.lean create mode 100755 Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean create mode 100755 Complexity/Algorithms/WordRAM/RadixSort/Model.lean create mode 100755 Complexity/Algorithms/WordRAM/RadixSort/Program.lean create mode 100755 Complexity/Algorithms/WordRAM/RadixSort/Scan.lean create mode 100755 Complexity/Computability/WordRAM/ListArray.lean create mode 100755 ComplexityTests/WordRAMRadixSort.lean diff --git a/Complexity.lean b/Complexity.lean index e6ac945545..f7f030e1ba 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -49,6 +49,9 @@ 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.Deduplicate.Count +public import Complexity.Algorithms.WordRAM.Deduplicate.Model +public import Complexity.Algorithms.WordRAM.Deduplicate.Program public import Complexity.Algorithms.WordRAM.Dijkstra.BinaryHeap public import Complexity.Algorithms.WordRAM.Dijkstra.ExactNotation public import Complexity.Algorithms.WordRAM.Dijkstra.HeapProgram @@ -97,6 +100,11 @@ public import Complexity.Algorithms.WordRAM.Queue.Dynamic public import Complexity.Algorithms.WordRAM.Queue.Model public import Complexity.Algorithms.WordRAM.Quickselect public import Complexity.Algorithms.WordRAM.Quicksort +public import Complexity.Algorithms.WordRAM.RadixSort +public import Complexity.Algorithms.WordRAM.RadixSort.Deduplicate +public import Complexity.Algorithms.WordRAM.RadixSort.Model +public import Complexity.Algorithms.WordRAM.RadixSort.Program +public import Complexity.Algorithms.WordRAM.RadixSort.Scan public import Complexity.Algorithms.WordRAM.Stack public import Complexity.Algorithms.WordRAM.Stack.Controller public import Complexity.Algorithms.WordRAM.ThreeSum.Cell @@ -1603,6 +1611,7 @@ public import Complexity.Computability.WordRAM.Dynamic public import Complexity.Computability.WordRAM.Function public import Complexity.Computability.WordRAM.Generic public import Complexity.Computability.WordRAM.Heap +public import Complexity.Computability.WordRAM.ListArray public import Complexity.Computability.WordRAM.MapOp public import Complexity.Computability.WordRAM.Multiword public import Complexity.Computability.WordRAM.Multiword.Arithmetic diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Count.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Count.lean new file mode 100755 index 0000000000..7ecf8a3c0a --- /dev/null +++ b/Complexity/Algorithms/WordRAM/Deduplicate/Count.lean @@ -0,0 +1,148 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.Deduplicate.Model +public import Complexity.Computability.WordRAM.ListArray + +/-! +# Native counting of one adjacent run + +The counter tests the array bound before every read and stops at the first different word. +No sorting or counting oracle is used. The mathematical end index is a proof parameter and +does not occur in the code or a register. + +## Main definitions + +* `Count.body`, `Count.code`: the native run-counting loop. + +## Main results + +* `Count.exec`: every matching input position advances the cursor, with linear charged time. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.Deduplicate + +open Cslib.WordRAM + +/-- Input length. -/ +abbrev len : Fin 16 := 0 +/-- Sorted input array. -/ +abbrev src : Fin 16 := 1 +/-- Unique-key output array. -/ +abbrev keys : Fin 16 := 2 +/-- Multiplicity output array. -/ +abbrev counts : Fin 16 := 3 +/-- Input cursor. -/ +abbrev idx : Fin 16 := 4 +/-- Number of completed output runs. -/ +abbrev out : Fin 16 := 5 +/-- Current run key. -/ +abbrev key : Fin 16 := 6 +/-- Current run multiplicity. -/ +abbrev count : Fin 16 := 7 +/-- Loaded input word. -/ +abbrev value : Fin 16 := 8 +/-- Address temporary. -/ +abbrev addr : Fin 16 := 9 +/-- Constant one. -/ +abbrev one : Fin 16 := 10 +/-- Run-counting loop flag. -/ +abbrev active : Fin 16 := 11 + +namespace Count + +/-- Read the next word only if the cursor is in bounds; count it or end the run. -/ +def body : Code .standard 16 := ram[ + if idx < len then [ + addr ← src + idx; value ← *addr; + if value = key then [ count ← count + one; idx ← idx + one ] + else [ active ← #0 ] + ] else [ active ← #0 ] ] + +/-- Count from the current cursor, starting with zero count and flag one. -/ +def code : Code .standard 16 := .loop (.eq active one) body + +/-- Count a maximal adjacent run in `8c+10` instructions, including the terminating test. +`stop` is a proof parameter: execution learns the boundary solely through charged comparisons +and loads. The source memory is unchanged. Empty runs are also admitted by this kernel. -/ +theorem exec {w s n start stop k : ℕ} (x : Data 16) + (hn : x.regs len = n) (hs : x.regs src = s) (hi : x.regs idx = start) + (hk : x.regs key = k) (hc : x.regs count = 0) + (ha : x.regs active = 1) (ho : x.regs one = 1) + (hle : start ≤ stop) (hstop : stop ≤ n) + (hmem : ∀ j, start ≤ j → j < stop → x.mem (s + j) = k) + (hnext : stop < n → x.mem (s + stop) ≠ k) + (hfit : Fits w (s + n)) : + ∃ y, ExecWithin w code x y (8 * (stop - start) + 10) ∧ + y.regs idx = stop ∧ y.regs count = stop - start ∧ y.mem = x.mem := by + let I (y : Data 16) : Prop := + y.regs len = n ∧ y.regs src = s ∧ y.regs key = k ∧ y.regs one = 1 ∧ + start ≤ y.regs idx ∧ y.regs idx ≤ stop ∧ y.regs count = y.regs idx - start ∧ + y.regs active ≤ 1 ∧ (y.regs active = 0 → y.regs idx = stop) ∧ y.mem = x.mem + have hbody (y : Data 16) (hy : I y) (hcond : (Cond.eq active one).Holds y) : + ∃ z, ExecWithin w body y z 6 ∧ I z ∧ + stop - z.regs idx + z.regs active < stop - y.regs idx + y.regs active := by + obtain ⟨hyn, hys, hyk, hyo, hlo, hhi, hyc, hya, hyend, hym⟩ := hy + have hyflag : y.regs active = 1 := by simpa only [Cond.holds_eq, hyo] using hcond + simp only [len, src, key, one, idx, count, active] at * + have hlenLink : y.regs len = y.regs 0 := rfl + have hidxLink : y.regs idx = y.regs 4 := rfl + have hcountLink : y.regs count = y.regs 7 := rfl + have hactiveLink : y.regs active = y.regs 11 := rfl + by_cases hlt : y.regs idx < stop + · have hidx : y.regs idx < y.regs len := by omega + have hread : y.mem (y.regs src + y.regs idx) = k := by + rw [hym, hys]; exact hmem _ hlo hlt + have had : Fits w (y.regs src + y.regs idx) := by + rw [hys]; exact hfit.of_le (by omega) + have hic : Fits w (y.regs idx + 1) := hfit.of_le (by omega) + have hcc : Fits w (y.regs count + 1) := hfit.of_le (by omega) + refine ⟨?_, ?run, ?_, ?_⟩ + case run => ram_eval [body, hidx, hread, hyk, hyo, had, hic, hcc] + · simp +decide only [I, Data.set_regs, Data.set_mem, Function.update_apply, + idx, count, value, addr, len, src, key, one, active, ↓reduceIte] + exact ⟨hyn, hys, hyk, hyo, by omega, by omega, by omega, by omega, by omega, hym⟩ + · simp +decide only [Data.set_regs, Function.update_apply, idx, count, value, addr, + ↓reduceIte] + omega + · have hend : y.regs idx = stop := by omega + by_cases hsn : stop < n + · have hidx : y.regs idx < y.regs len := by omega + have hread : y.mem (y.regs src + y.regs idx) ≠ y.regs key := by + simpa only [hym, hys, hend, hyk] using hnext hsn + have had : Fits w (y.regs src + y.regs idx) := by + rw [hys]; exact hfit.of_le (by omega) + refine ⟨?_, ?run, ?_, ?_⟩ + case run => ram_eval [body, hidx, hread, had] + · simp +decide only [I, Data.set_regs, Data.set_mem, Function.update_apply, + idx, count, value, addr, len, src, key, one, active, ↓reduceIte] + exact ⟨hyn, hys, hyk, hyo, hlo, hhi, hyc, by omega, fun _ => hend, hym⟩ + · simp +decide [Data.set_regs, hyflag, hend] + · have hidx : ¬y.regs idx < y.regs len := by omega + refine ⟨?_, ?run, ?_, ?_⟩ + case run => ram_eval [body, hidx] + · simp +decide only [I, Data.set_regs, Data.set_mem, Function.update_apply, + idx, count, len, src, key, one, active, ↓reduceIte] + exact ⟨hyn, hys, hyk, hyo, hlo, hhi, hyc, by omega, fun _ => hend, hym⟩ + · simp +decide [Data.set_regs, hyflag, hend] + have hx : I x := ⟨hn, hs, hk, ho, by omega, by omega, by omega, by omega, by omega, rfl⟩ + obtain ⟨y, he, hy, hnot⟩ := ExecWithin.loop_invariant I + (fun y => stop - y.regs idx + y.regs active) hbody hx + obtain ⟨_, _, _, hyo, _, _, hyc, hya, hyend, hym⟩ := hy + have hz : y.regs active = 0 := by + simp only [Cond.holds_eq, hyo] at hnot + omega + have hend := hyend hz + refine ⟨y, ?_, hend, by simpa only [hend] using hyc, hym⟩ + simpa only [code, hi, ha, Nat.add_mul, Nat.one_mul, Nat.mul_comm, Nat.add_assoc] using he + +end Count +end Cslib.Algorithms.WordRAM.Deduplicate diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Model.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Model.lean new file mode 100755 index 0000000000..fb6feb225b --- /dev/null +++ b/Complexity/Algorithms/WordRAM/Deduplicate/Model.lean @@ -0,0 +1,146 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Init +public import Mathlib.Data.List.Sort +public import Mathlib.Data.List.Nodup + +/-! +# Adjacent deduplication with full multiplicities + +The output is a list of key/count pairs with positive counts and strictly increasing keys. +Expanding those runs recovers the sorted input exactly. This specifies the whole finite set and +all multiplicities, rather than a one-way membership condition. These proofs are original. + +## Main definitions + +* `expand`: expand key/count pairs into their repeated keys. +* `RunEncoding`: positive, strictly ordered runs reconstructing the input. + +## Main results + +* `RunEncoding.keys_nodup`: each key occurs once. +* `RunEncoding.mem_iff`: the keys enumerate exactly the input set. +* `RunEncoding.multiplicity`: each count is the full input multiplicity. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.Deduplicate + +/-- Expand a list of key/count pairs into repeated keys. -/ +def expand (rs : List (ℕ × ℕ)) : List ℕ := rs.flatMap (fun p => List.replicate p.2 p.1) + +/-- Runs with positive counts and strictly increasing keys reconstruct exactly `xs`. -/ +structure RunEncoding (xs : List ℕ) (rs : List (ℕ × ℕ)) : Prop where + /-- Each emitted run is nonempty. -/ + positive : ∀ p ∈ rs, 0 < p.2 + /-- A key never reappears after a different key. -/ + ordered : rs.Pairwise (fun p q => p.1 < q.1) + /-- All input positions are accounted for, with their exact order and multiplicity. -/ + reconstruct : expand rs = xs + +/-- Expansion distributes over concatenation. -/ +theorem expand_append (rs ss : List (ℕ × ℕ)) : expand (rs ++ ss) = expand rs ++ expand ss := by + simp [expand] + +/-- Expansion of one run. -/ +theorem expand_singleton (k c : ℕ) : expand [(k, c)] = List.replicate c k := by + simp [expand] + +/-- Every key in an expansion occurs in some nonempty run. -/ +theorem mem_expand {rs : List (ℕ × ℕ)} {k : ℕ} : + k ∈ expand rs ↔ ∃ p ∈ rs, p.1 = k ∧ 0 < p.2 := by + simp only [expand, List.mem_flatMap, List.mem_replicate] + constructor + · rintro ⟨p, hp, hc, hk⟩ + exact ⟨p, hp, hk.symm, by omega⟩ + · rintro ⟨p, hp, hk, hc⟩ + exact ⟨p, hp, by omega, hk.symm⟩ + +/-- Positive runs contain at least one element each. -/ +theorem length_le_expand {rs : List (ℕ × ℕ)} (hp : ∀ p ∈ rs, 0 < p.2) : + rs.length ≤ (expand rs).length := by + induction rs with + | nil => simp [expand] + | cons p rs ih => + have hc := hp p (by simp) + have ht := ih (fun q hq => hp q (by simp [hq])) + simp only [expand, List.flatMap_cons, List.length_append, List.length_replicate, + List.length_cons] at * + omega + +private theorem count_expand {rs : List (ℕ × ℕ)} + (ho : rs.Pairwise (fun p q => p.1 < q.1)) {p : ℕ × ℕ} (hp : p ∈ rs) : + (expand rs).count p.1 = p.2 := by + induction rs with + | nil => simp at hp + | cons q rs ih => + obtain ⟨hq, ht⟩ := List.pairwise_cons.mp ho + rcases List.mem_cons.mp hp with rfl | hp + · have hnot : p.1 ∉ expand rs := by + rintro hk + obtain ⟨p, hp, he, _⟩ := mem_expand.mp hk + have := hq p hp + omega + change (List.replicate p.2 p.1 ++ expand rs).count p.1 = p.2 + rw [List.count_append, List.count_replicate_self, List.count_eq_zero_of_not_mem hnot, + Nat.add_zero] + · have hlt := hq p hp + have hc := ih ht hp + simp only [expand, List.flatMap_cons, List.count_append] + rw [show (List.replicate q.2 q.1).count p.1 = 0 by + simp [List.count_replicate, Nat.ne_of_lt hlt]] + simpa only [Nat.zero_add, expand] using hc + +/-- The key array contains no repeated key. -/ +theorem RunEncoding.keys_nodup {xs : List ℕ} {rs : List (ℕ × ℕ)} + (h : RunEncoding xs rs) : (rs.map Prod.fst).Nodup := by + rw [List.nodup_iff_pairwise_ne, List.pairwise_map] + exact h.ordered.imp (fun hlt => Nat.ne_of_lt hlt) + +/-- The key array enumerates precisely the input set. -/ +theorem RunEncoding.mem_iff {xs : List ℕ} {rs : List (ℕ × ℕ)} + (h : RunEncoding xs rs) (k : ℕ) : k ∈ rs.map Prod.fst ↔ k ∈ xs := by + rw [← h.reconstruct, mem_expand] + rw [List.mem_map] + constructor + · rintro ⟨p, hp, hk⟩ + exact ⟨p, hp, hk, h.positive p hp⟩ + · rintro ⟨p, hp, hk, _⟩ + exact ⟨p, hp, hk⟩ + +/-- Each emitted count equals the full multiplicity in the input list. -/ +theorem RunEncoding.multiplicity {xs : List ℕ} {rs : List (ℕ × ℕ)} + (h : RunEncoding xs rs) {p : ℕ × ℕ} (hp : p ∈ rs) : xs.count p.1 = p.2 := by + rw [← h.reconstruct] + exact count_expand h.ordered hp + +/-- A nonempty sorted list starts with a positive run of its head, followed only by larger +keys. This is a mathematical decomposition used to verify the adjacent scan, not an oracle +operation of the program. -/ +theorem first_run (k : ℕ) (xs : List ℕ) (hs : (k :: xs).Pairwise (· ≤ ·)) : + ∃ c rest, 0 < c ∧ k :: xs = List.replicate c k ++ rest ∧ ∀ v ∈ rest, k < v := by + induction xs generalizing k with + | nil => exact ⟨1, [], by decide, rfl, by simp⟩ + | cons j xs ih => + obtain ⟨hk, ht⟩ := List.pairwise_cons.mp hs + by_cases heq : k = j + · subst j + obtain ⟨c, rest, hc, he, hr⟩ := ih k ht + refine ⟨c + 1, rest, by omega, ?_, hr⟩ + simpa only [List.replicate_succ, List.cons_append] using congrArg (List.cons k) he + · have hkj : k < j := lt_of_le_of_ne (hk j (by simp)) heq + refine ⟨1, j :: xs, by decide, rfl, ?_⟩ + intro v hv + rcases List.mem_cons.mp hv with rfl | hv + · exact hkj + · exact lt_of_lt_of_le hkj ((List.pairwise_cons.mp ht).1 v hv) + +end Cslib.Algorithms.WordRAM.Deduplicate diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Program.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Program.lean new file mode 100755 index 0000000000..df4da01c3d --- /dev/null +++ b/Complexity/Algorithms/WordRAM/Deduplicate/Program.lean @@ -0,0 +1,284 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.Deduplicate.Count + +/-! +# Complete native adjacent deduplication and multiplicity arrays + +The input is a sorted array. The program counts one maximal adjacent run, emits its key and +full count, and repeats. The output consists of two adjacent arrays, each with capacity equal +to the input length. All run tests, loads, stores, address calculations and control instructions +are charged. The total time is linear, including repeated keys and empty inputs. + +## Main definitions + +* `flush`: store one key/count pair and advance the output cursor. +* `body`: count and emit one complete run. +* `code`: the full adjacent deduplication program. + +## Main results + +* `exec`: full unique-key and multiplicity output with a linear native instruction bound. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.Deduplicate + +open Cslib.WordRAM + +/-- Store one key/count pair and advance the output cursor. -/ +def flush : Code .standard 16 := ram[ + addr ← keys + out; *addr ← key; + addr ← counts + out; *addr ← count; out ← out + one ] + +/-- Count and emit the next maximal run. -/ +def body : Code .standard 16 := + .seq (ram[ addr ← src + idx; key ← *addr; count ← #0; active ← #1 ]) + (.seq Count.code flush) + +/-- Initialize the two cursors and count all adjacent runs. -/ +def code : Code .standard 16 := + .seq (ram[ one ← #1; idx ← #0; out ← #0 ]) (.loop (.lt idx len) body) + +private theorem flush_exec {w : ℕ} (x : Data 16) (ho : x.regs one = 1) + (hk : Fits w (x.regs keys + x.regs out)) + (hc : Fits w (x.regs counts + x.regs out)) (hi : Fits w (x.regs out + 1)) : + ∃ y, ExecWithin w flush x y 5 ∧ y.regs out = x.regs out + 1 ∧ + y.mem = Function.update (Function.update x.mem (x.regs keys + x.regs out) (x.regs key)) + (x.regs counts + x.regs out) (x.regs count) := by + refine ⟨?_, ?run, ?_, ?_⟩ + case run => ram_eval [flush, ho, hk, hc, hi] + all_goals simp +decide [Data.set_regs, Data.write_mem] + +/-- Emit each distinct key once and its full multiplicity in at most `29n+5` instructions. +The input precedes the two output blocks; no input cells or cells outside the output blocks +change. Counts are ordinary unsigned words bounded by `n`. No restriction on key signs or +canonical encodings is needed: in particular canonical `Signed.encode` inputs are admitted. -/ +theorem exec {w s base : ℕ} (xs : List ℕ) (x : Data 16) + (hn : x.regs len = xs.length) (hs : x.regs src = s) + (hk : x.regs keys = base) (hc : x.regs counts = base + xs.length) + (hmem : Block.HoldsList s xs x.mem) (hsorted : xs.Pairwise (· ≤ ·)) + (hbelow : s + xs.length ≤ base) (hfit : Fits w (base + 2 * xs.length)) + (h1 : Fits w 1) : + ∃ y rs, ExecWithin w code x y (29 * xs.length + 5) ∧ + RunEncoding xs rs ∧ y.regs out = rs.length ∧ + Block.HoldsList base (rs.map Prod.fst) y.mem ∧ + Block.HoldsList (base + xs.length) (rs.map Prod.snd) y.mem ∧ + Block.HoldsList s xs y.mem ∧ + Modifies (Block.mk base (2 * xs.length)).footprint x.mem y.mem := by + let n := xs.length + let I (y : Data 16) : Prop := + y.regs len = n ∧ y.regs src = s ∧ y.regs keys = base ∧ + y.regs counts = base + n ∧ y.regs one = 1 ∧ y.regs idx ≤ n ∧ + Block.HoldsList s xs y.mem ∧ + Modifies (Block.mk base (2 * n)).footprint x.mem y.mem ∧ + ∃ rs, RunEncoding (xs.take (y.regs idx)) rs ∧ y.regs out = rs.length ∧ + Block.HoldsList base (rs.map Prod.fst) y.mem ∧ + Block.HoldsList (base + n) (rs.map Prod.snd) y.mem ∧ + ∀ p ∈ rs, ∀ v ∈ xs.drop (y.regs idx), p.1 < v + have hbody (y : Data 16) (hy : I y) (hcond : (Cond.lt idx len).Holds y) : + ∃ z T, ExecWithin w body y z T ∧ I z ∧ + T + 2 + 29 * (n - z.regs idx) ≤ 29 * (n - y.regs idx) := by + obtain ⟨hyn, hys, hyk, hyc, hyo, hiy, hym, hyf, rs, hrs, hyout, hyrk, hyrc, hrlt⟩ := hy + have hi : y.regs idx < xs.length := by simpa only [Cond.holds_lt, hyn] using hcond + have htail := hsorted.sublist (List.drop_sublist (y.regs idx) xs) + rw [List.drop_eq_getElem_cons hi] at htail + obtain ⟨c, rest, hcp, heq, hrest⟩ := first_run xs[y.regs idx] (xs.drop (y.regs idx + 1)) htail + have hdrop : xs.drop (y.regs idx) = List.replicate c xs[y.regs idx] ++ rest := by + rw [List.drop_eq_getElem_cons hi, heq] + have hlen := congrArg List.length hdrop + simp only [List.length_drop, List.length_append, List.length_replicate] at hlen + have hstop : y.regs idx + c ≤ n := by omega + have hdrop' : xs.drop (y.regs idx + c) = rest := by + rw [← List.drop_drop, hdrop] + exact List.drop_left' (by simp) + have htake : xs.take (y.regs idx + c) = + xs.take (y.regs idx) ++ List.replicate c xs[y.regs idx] := by + rw [List.take_add, hdrop, List.take_left' (by simp)] + have hrun (j : ℕ) (hlo : y.regs idx ≤ j) (hhi : j < y.regs idx + c) : + y.mem (s + j) = xs[y.regs idx] := by + rw [hym.get j (by omega)] + have hh : (xs.drop (y.regs idx))[j - y.regs idx]'(by simp only [List.length_drop]; omega) = + xs[y.regs idx] := by + simp [hdrop, show j - y.regs idx < c from by omega] + simpa only [List.getElem_drop, show y.regs idx + (j - y.regs idx) = j from by omega] using hh + have hnext (hlt : y.regs idx + c < n) : + y.mem (s + (y.regs idx + c)) ≠ xs[y.regs idx] := by + rw [hym.get _ hlt] + have hmember : xs[y.regs idx + c] ∈ rest := by + rw [← hdrop'] + have hh := List.getElem_mem (l := xs.drop (y.regs idx + c)) + (n := 0) (by simp only [List.length_drop]; omega) + simpa using hh + have := hrest _ hmember + omega + let y₀ := (((y.set addr (s + y.regs idx)).set key xs[y.regs idx]).set count 0).set active 1 + have hinit : ExecWithin (B := .standard) w + (ram[ addr ← src + idx; key ← *addr; count ← #0; active ← #1 ]) y y₀ 4 := by + have had : Fits w (y.regs src + y.regs idx) := by + rw [hys]; exact hfit.of_le (by omega) + have hread : y.mem (y.regs src + y.regs idx) = xs[y.regs idx] := by + rw [hys]; exact hym.get _ hi + have had' : Fits w (s + y.regs idx) := hfit.of_le (by omega) + have hread' : y.mem (s + y.regs idx) = xs[y.regs idx] := hym.get _ hi + ram_eval [y₀, had, had', hread, hread', hys, h1] + obtain ⟨u, hu, hui, huc, hum⟩ := Count.exec (k := xs[y.regs idx]) y₀ + (by simpa +decide [y₀] using hyn) (by simpa +decide [y₀] using hys) + (by simp +decide [y₀]) (by simp +decide [y₀]) + (by simp +decide [y₀]) (by simp +decide [y₀]) + (by simpa +decide [y₀] using hyo) (show y.regs idx ≤ y.regs idx + c by omega) + hstop (by simpa only [y₀, Data.set_mem] using hrun) + (by simpa only [y₀, Data.set_mem] using hnext) (hfit.of_le (by omega)) + have huc' : u.regs count = c := by simpa using huc + have hur (j : Fin 16) (hj : Count.code.Preserves j) : u.regs j = y₀.regs j := hu.regs_eq hj + have hun : u.regs len = n := by + simpa +decide [y₀] + using (hur len (by decide)).trans + (by simpa +decide [y₀] using hyn) + have hus : u.regs src = s := by + simpa +decide [y₀] + using (hur src (by decide)).trans + (by simpa +decide [y₀] using hys) + have huk : u.regs keys = base := by + simpa +decide [y₀] + using (hur keys (by decide)).trans + (by simpa +decide [y₀] using hyk) + have hucb : u.regs counts = base + n := by + simpa +decide [y₀] + using (hur counts (by decide)).trans + (by simpa +decide [y₀] using hyc) + have huo : u.regs one = 1 := by + simpa +decide [y₀] + using (hur one (by decide)).trans + (by simpa +decide [y₀] using hyo) + have huout : u.regs out = rs.length := by + simpa +decide [y₀] + using (hur out (by decide)).trans + (by simpa +decide [y₀] using hyout) + have hukey : u.regs key = xs[y.regs idx] := by simpa +decide [y₀] using hur key (by decide) + have hrl : rs.length ≤ y.regs idx := by + have hh := length_le_expand hrs.positive + rw [hrs.reconstruct, List.length_take_of_le hiy] at hh + exact hh + obtain ⟨z, hz, hzout, hzm⟩ := flush_exec u huo + (by rw [huk, huout]; exact hfit.of_le (by omega)) + (by rw [hucb, huout]; exact hfit.of_le (by omega)) + (by rw [huout]; exact hfit.of_le (by omega)) + have hzmem : z.mem = Function.update + (Function.update y.mem (base + rs.length) xs[y.regs idx]) (base + n + rs.length) c := by + simpa only [huk, hucb, huout, hukey, huc', hum, y₀, Data.set_mem] using hzm + have hzi : z.regs idx = y.regs idx + c := (hz.regs_eq (by decide)).trans hui + have hzf : Modifies (Block.mk base (2 * n)).footprint y.mem z.mem := by + intro a hnot + have ha : ¬(base ≤ a ∧ a < base + 2 * n) := by + simpa only [Block.mem_footprint, Block.Owns] using hnot + rw [hzmem, Function.update_of_ne (by omega), Function.update_of_ne (by omega)] + have hkeyrep : Block.HoldsList base + ((rs ++ [(xs[y.regs idx], c)]).map Prod.fst) z.mem := by + rw [List.map_append] + simp only [List.map_cons, List.map_nil] + rw [hzmem] + have hkr : Block.HoldsList base (rs.map Prod.fst ++ [xs[y.regs idx]]) + (Function.update y.mem (base + rs.length) xs[y.regs idx]) := by + simpa only [List.length_map] using hyrk.append_one xs[y.regs idx] + apply hkr.frame + intro a hlo hhi + rw [Function.update_of_ne (by + simp only [List.length_append, List.length_map, List.length_singleton] at hhi + omega)] + have hcountrep : Block.HoldsList (base + n) + ((rs ++ [(xs[y.regs idx], c)]).map Prod.snd) z.mem := by + rw [List.map_append] + simp only [List.map_cons, List.map_nil] + rw [hzmem] + have hcr : Block.HoldsList (base + n) (rs.map Prod.snd) + (Function.update y.mem (base + rs.length) xs[y.regs idx]) := by + apply hyrc.frame + intro a hlo hhi + rw [Function.update_of_ne (by omega)] + simpa only [List.length_map] using hcr.append_one c + have hrs' : RunEncoding (xs.take (y.regs idx + c)) (rs ++ [(xs[y.regs idx], c)]) := by + refine ⟨?_, ?_, ?_⟩ + · intro p hp + rcases List.mem_append.mp hp with hp | hp + · exact hrs.positive p hp + · have hp' := List.mem_singleton.mp hp + subst p + exact hcp + · apply List.pairwise_append.mpr + refine ⟨hrs.ordered, by simp, ?_⟩ + intro p hp q hq + have hq' : q = (xs[y.regs idx], c) := List.mem_singleton.mp hq + subst q + apply hrlt p hp + rw [List.drop_eq_getElem_cons hi] + exact List.mem_cons_self + · rw [expand_append, expand_singleton, hrs.reconstruct, htake] + refine ⟨z, 8 * c + 19, ?_, ?_, ?_⟩ + · have hall := hinit.seq (hu.seq hz) + exact hall.mono (by omega) + · refine ⟨(hz.regs_eq (by decide)).trans hun, (hz.regs_eq (by decide)).trans hus, + (hz.regs_eq (by decide)).trans huk, (hz.regs_eq (by decide)).trans hucb, + (hz.regs_eq (by decide)).trans huo, by omega, ?_, ?_, + rs ++ [(xs[y.regs idx], c)], by simpa only [hzi] using hrs', ?_, hkeyrep, hcountrep, ?_⟩ + · apply hym.frame + intro a hlo hhi + apply hzf + rw [Block.mem_footprint] + change ¬(base ≤ a ∧ a < base + 2 * n) + omega + · intro a ha + exact (hzf a ha).trans (hyf a ha) + · rw [hzout, huout, List.length_append, List.length_singleton] + · intro p hp v hv + rw [hzi, hdrop'] at hv + rcases List.mem_append.mp hp with hp | hp + · apply hrlt p hp + rw [hdrop] + exact List.mem_append_right _ hv + · have hp' : p = (xs[y.regs idx], c) := List.mem_singleton.mp hp + subst p + exact hrest v hv + · rw [hzi] + omega + let x₀ := ((x.set one 1).set idx 0).set out 0 + have hx₀ : I x₀ := by + refine ⟨?_, ?_, ?_, ?_, ?_, ?_, hmem, Modifies.refl _ _, [], ?_, ?_, ?_, ?_, ?_⟩ + · simpa +decide [x₀] using hn + · simpa +decide [x₀] using hs + · simpa +decide [x₀] using hk + · simpa +decide [x₀] using hc + · simp +decide [x₀] + · simp +decide [x₀] + · refine ⟨by simp, by simp, ?_⟩ + simp +decide [x₀, expand] + · simp [x₀] + · exact Block.holdsList_nil _ _ + · exact Block.holdsList_nil _ _ + · simp + obtain ⟨y, he, hy, hcond⟩ := ExecWithin.loop_potential I + (fun y => 29 * (n - y.regs idx)) hbody hx₀ + obtain ⟨hyn, _, _, _, _, hiy, hym, hyf, rs, hrs, hyout, hyrk, hyrc, _⟩ := hy + have hend : y.regs idx = xs.length := by + simp only [Cond.holds_lt, hyn] at hcond + omega + have hinit : ExecWithin (B := .standard) w (ram[ one ← #1; idx ← #0; out ← #0 ]) x x₀ 3 := by + ram_eval [x₀, h1] + refine ⟨y, rs, ?_, by simpa only [hend, List.take_length] using hrs, + hyout, hyrk, hyrc, hym, hyf⟩ + have hall := hinit.seq he + convert hall using 1 + · rfl + · simp +decide only [x₀, Data.set_regs, Function.update_apply, idx, out, one, + ↓reduceIte, Nat.sub_zero, n] + omega + +end Cslib.Algorithms.WordRAM.Deduplicate diff --git a/Complexity/Algorithms/WordRAM/RadixSort.lean b/Complexity/Algorithms/WordRAM/RadixSort.lean new file mode 100755 index 0000000000..d3a1a6584c --- /dev/null +++ b/Complexity/Algorithms/WordRAM/RadixSort.lean @@ -0,0 +1,106 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.RadixSort.Deduplicate +public import Complexity.Computability.WordRAM.Time + +/-! +# Fixed-code native sorting and multiplicity construction from machine input + +The input is the ordinary native list of unsigned words. The code reads its length and the +physical word width, and uses that width as an explicit key-bit bound. The lower-level +`deduplicateCode` accepts a smaller caller-supplied bit bound when known. + +The key array starts at `1+2n`, the count array at `1+3n`; the return value is their common +logical length. The two sorting buffers occupy `[1,1+2n)`. The implementation accepts repeated +keys, zero keys and canonical signed integer encodings. Its order is unsigned encoding order. + +## Main definitions + +* `inputCode`, `program`: the complete structured and flat native programs. +* `layoutWidth`: a logarithmic width sufficient for all four array blocks. + +## Main results + +* `input_exec`: charged execution from ordinary native input, with exact key/count output. +* `layout_fits`: logarithmic-width layout fit, including empty input. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.RadixSort + +open Cslib.WordRAM + +/-- Read the native input header and width, then sort and construct both output arrays. -/ +def inputCode : Code .standard 16 := + .seq (ram[ addr ← #0; len ← *addr; src ← #1; dst ← src + len; bits ← width% ]) + deduplicateCode + +/-- The fixed native instruction list, returning the number of unique keys. -/ +def program : Program .standard 16 := inputCode.program Deduplicate.out + +/-- A logarithmic layout-width schedule; the key magnitude bound is independent. -/ +def layoutWidth (n : ℕ) : ℕ := logWidth 3 (n + 1) + +/-- Four array blocks and both initialization constants fit at the layout-width threshold. +This does not by itself bound input key magnitudes. -/ +theorem layout_fits (n : ℕ) : Fits (layoutWidth n) (1 + 4 * n) ∧ Fits (layoutWidth n) 2 := by + have hs : 4 ≤ (n + 2) ^ 2 := Nat.pow_le_pow_left (show 2 ≤ n + 2 by omega) 2 + have hm := Nat.mul_le_mul_right (n + 2) hs + have hp : 4 * (n + 2) ≤ (n + 2) ^ 3 := by simpa only [Nat.pow_succ] using hm + constructor + · apply fits_of_lt_pow + have he : n + 1 + 1 = n + 2 := by omega + rw [he] + omega + · apply fits_of_lt_pow + have he : n + 1 + 1 = n + 2 := by omega + rw [he] + omega + +/-- Execute the complete finite native program on the standard input convention. The bound +includes the final halt and every setup/address/control instruction. Input fit and the explicit +`4n+1` layout fit are independent requirements. Empty and singleton inputs need no exceptions. +The output words enumerate exactly the input set and each stored count is its full multiplicity. -/ +theorem input_exec {w : ℕ} (xs : List ℕ) (hinput : InputFits w xs) + (hfit : Fits w (1 + 4 * xs.length)) (h2 : Fits w 2) : + ∃ (y : Data 16) (rs : List (ℕ × ℕ)) (steps : ℕ), + Exec w inputCode (Data.init xs) y steps ∧ + Executes w program (State.init xs) rs.length (steps + 1) ∧ + steps + 1 ≤ deduplicateTime xs.length w + 6 ∧ + Block.HoldsList (1 + 2 * xs.length) (rs.map Prod.fst) y.mem ∧ + Block.HoldsList (1 + 3 * xs.length) (rs.map Prod.snd) y.mem ∧ + (rs.map Prod.fst).Nodup ∧ + (∀ k, k ∈ rs.map Prod.fst ↔ k ∈ xs) ∧ + (∀ p ∈ rs, 0 < p.2 ∧ xs.count p.1 = p.2) ∧ + Modifies (Block.mk 1 (4 * xs.length)).footprint (Data.init xs : Data 16).mem y.mem := by + let x := Data.init xs (r := 16) + have h1 := h2.of_le (show 1 ≤ 2 by omega) + have had : Fits w (1 + xs.length) := hfit.of_le (by omega) + let u := ((((x.set addr 0).set len xs.length).set src 1).set dst (1 + xs.length)).set bits w + have hi : ExecWithin w + (ram[ addr ← #0; len ← *addr; src ← #1; dst ← src + len; bits ← width% ]) x u 5 := by + ram_eval [u, x, Data.init_mem_zero, h1, had] + have hn : u.regs len = xs.length := by simp +decide [u] + have hs : u.regs src = 1 := by simp +decide [u] + have hd : u.regs dst = 1 + xs.length := by simp +decide [u] + have hb : u.regs bits = w := by simp [u] + have hm : Block.HoldsList 1 xs u.mem := by + intro i + change (Data.init xs : Data 16).mem (1 + i.val) = xs[i.val] + rw [Nat.add_comm, Data.init_mem_succ] + simp only [List.getD, List.getElem?_eq_getElem i.isLt, Option.getD_some] + obtain ⟨y, rs, he, hyout, hyrk, hyrc, hnodup, hset, hcount, hframe⟩ := + deduplicate_exec xs u hn hs hd hb hm hinput.2 hinput.2 hfit (fits_self w) h2 + obtain ⟨steps, hsteps, hexec⟩ := hi.seq he + refine ⟨y, rs, steps, hexec, ?_, by omega, hyrk, hyrc, hnodup, hset, hcount, hframe⟩ + simpa only [program, inputCode, hyout, State.init, x] using hexec.executes Deduplicate.out + +end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean b/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean new file mode 100755 index 0000000000..eed159b431 --- /dev/null +++ b/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean @@ -0,0 +1,135 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.RadixSort.Program +public import Complexity.Algorithms.WordRAM.Deduplicate.Program +public import Complexity.Computability.WordRAM.Signed.Encoding + +/-! +# Native radix sorting followed by unique-key/multiplicity array construction + +This composes the complete native sorter and adjacent-run builder. Two length-`n` sorting +buffers are followed by length-`n` key and count output buffers. All four blocks lie in one +explicit frame. The returned keys are in unsigned encoding order; decoding signed keys does +not turn that order into signed integer order. + +## Main definitions + +* `deduplicateCode`: fixed 16-register sorting and multiplicity program. +* `deduplicateTime`: complete structured instruction bound. + +## Main results + +* `deduplicate_exec`: exact set enumeration and full multiplicities, with no sorting premise. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.RadixSort + +open Cslib.WordRAM + +/-- Original base retained while radix sorting swaps its two pointers. -/ +abbrev origin : Fin 16 := 14 + +/-- Install the two output pointers after the two sorting buffers. -/ +def outputLayout : Code .standard 16 := ram[ + Deduplicate.keys ← len + len; Deduplicate.keys ← Deduplicate.keys + origin; + Deduplicate.counts ← Deduplicate.keys + len ] + +/-- Complete sorting and deduplication; `Deduplicate.out` returns the number of unique words. -/ +def deduplicateCode : Code .standard 16 := + .seq (ram[ origin ← src ]) (.seq code (.seq outputLayout Deduplicate.code)) + +/-- All instructions, including pointer construction and both complete algorithms, excluding +only the final halt. For positive `n` and key-bit bound this is `O(n * keyBits)`. -/ +def deduplicateTime (n keyBits : ℕ) : ℕ := keyBits * (22 * n + 15) + 29 * n + 14 + +/-- Complete fixed-code radix sorting and adjacent deduplication, with full multiplicities. +Keys below `2^keyBits` may repeat arbitrarily. Only the four length-`n` buffers change. The +conditions constrain arithmetic/address fit, not sorting or deduplication behavior. -/ +theorem deduplicate_exec {w base keyBits : ℕ} (xs : List ℕ) (x : Data 16) + (hn : x.regs len = xs.length) (hs : x.regs src = base) + (hd : x.regs dst = base + xs.length) (hb : x.regs bits = keyBits) + (hmem : Block.HoldsList base xs x.mem) + (hkeys : ∀ v ∈ xs, v < 2 ^ keyBits) (hwords : ∀ v ∈ xs, Fits w v) + (hfit : Fits w (base + 4 * xs.length)) (hbits : Fits w keyBits) (h2 : Fits w 2) : + ∃ (y : Data 16) (rs : List (ℕ × ℕ)), + ExecWithin w deduplicateCode x y (deduplicateTime xs.length keyBits) ∧ + y.regs Deduplicate.out = rs.length ∧ + Block.HoldsList (base + 2 * xs.length) (rs.map Prod.fst) y.mem ∧ + Block.HoldsList (base + 3 * xs.length) (rs.map Prod.snd) y.mem ∧ + (rs.map Prod.fst).Nodup ∧ + (∀ k, k ∈ rs.map Prod.fst ↔ k ∈ xs) ∧ + (∀ p ∈ rs, 0 < p.2 ∧ xs.count p.1 = p.2) ∧ + Modifies (Block.mk base (4 * xs.length)).footprint x.mem y.mem := by + let x₀ := x.set origin base + have hinit : ExecWithin (B := .standard) w (ram[ origin ← src ]) x x₀ 1 := by ram_eval [x₀, hs] + obtain ⟨u, hu, hum, hp, hsorted, hptr, huf⟩ := exec xs x₀ + (by simpa +decide [x₀] using hn) (by simpa +decide [x₀] using hs) + (by simpa +decide [x₀] using hd) (by simpa +decide [x₀] using hb) + hmem hkeys hwords (hfit.of_le (by omega)) hbits h2 + have hun : u.regs len = xs.length := by + simpa +decide [x₀] using (hu.regs_eq (i := len) (by decide)).trans + (show x₀.regs len = xs.length by simpa +decide [x₀] using hn) + have huo : u.regs origin = base := by + simpa only [x₀, Data.set_regs, Function.update_self] using hu.regs_eq (i := origin) (by decide) + have htwice : Fits w (u.regs len + u.regs len) := hfit.of_le (by rw [hun]; omega) + have hkeyaddr : Fits w (u.regs len + u.regs len + base) := hfit.of_le (by rw [hun]; omega) + have hcntaddr : Fits w (u.regs len + u.regs len + base + u.regs len) := + hfit.of_le (by rw [hun]; omega) + let v := ((u.set Deduplicate.keys (u.regs len + u.regs len)).set + Deduplicate.keys (u.regs len + u.regs len + base)).set + Deduplicate.counts (u.regs len + u.regs len + base + u.regs len) + have hl : ExecWithin w outputLayout u v 3 := by + ram_eval [v, outputLayout, huo, htwice, hkeyaddr, hcntaddr] + have hvl : v.regs Deduplicate.len = (passes keyBits xs).length := by + change u.regs len = (passes keyBits xs).length + rw [hun, hp.length_eq] + have hvs : v.regs Deduplicate.src = u.regs src := rfl + have hvk : v.regs Deduplicate.keys = base + 2 * xs.length := by + change u.regs len + u.regs len + base = base + 2 * xs.length + omega + have hvc : v.regs Deduplicate.counts = base + 3 * xs.length := by + change u.regs len + u.regs len + base + u.regs len = base + 3 * xs.length + omega + have hlen : (passes keyBits xs).length = xs.length := hp.length_eq + obtain ⟨y, rs, hy, hrs, hyout, hyrk, hyrc, _, hyf⟩ := Deduplicate.exec (passes keyBits xs) v + hvl hvs hvk (by rw [hvc, hlen]; omega) hum hsorted + (by rw [hlen]; rcases hptr with hptr | hptr <;> omega) + (by rw [hlen]; exact hfit.of_le (by omega)) (h2.of_le (by omega)) + have heq : base + 2 * xs.length + xs.length = base + 3 * xs.length := by omega + refine ⟨y, rs, ?_, hyout, hyrk, by simpa only [hlen, heq] using hyrc, + hrs.keys_nodup, ?_, ?_, ?_⟩ + · have hall := hinit.seq (hu.seq (hl.seq hy)) + apply hall.mono + simp only [deduplicateTime, timeBound, hlen] + omega + · intro k + exact (hrs.mem_iff k).trans hp.mem_iff + · intro p hp' + exact ⟨hrs.positive p hp', (hp.count_eq p.1).symm.trans (hrs.multiplicity hp')⟩ + · intro a hnot + have hnotSort : a ∉ (Block.mk base (2 * xs.length)).footprint := by + intro ha + apply hnot + rw [Block.mem_footprint] at ha ⊢ + dsimp only [Block.Owns] at ha ⊢ + omega + have hnotOut : a ∉ (Block.mk (base + 2 * xs.length) + (2 * (passes keyBits xs).length)).footprint := by + intro ha + apply hnot + rw [Block.mem_footprint] at ha ⊢ + dsimp only [Block.Owns] at ha ⊢ + rw [hlen] at ha + omega + exact (hyf a hnotOut).trans (huf a hnotSort) + +end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Model.lean b/Complexity/Algorithms/WordRAM/RadixSort/Model.lean new file mode 100755 index 0000000000..c92aeef96c --- /dev/null +++ b/Complexity/Algorithms/WordRAM/RadixSort/Model.lean @@ -0,0 +1,134 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Init +public import Mathlib.Data.List.Sort + +/-! +# Stable binary radix sorting: mathematical specification + +The order is the unsigned order on encoded words. In particular this does not claim that +sign-and-magnitude integer encodings are ordered as integers. The construction and proofs here +are original; the executable native program is supplied separately. + +## Main definitions + +* `digit`: one binary digit, represented as zero or one. +* `select`: the stable subsequence with a specified digit. +* `partition`: zeros followed by ones, preserving the order within each part. +* `passes`: least-significant-digit radix passes. + +## Main results + +* `partition_perm`: a pass preserves the multiset. +* `passes_correct`: sufficiently many passes give a sorted permutation. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.RadixSort + +/-- Binary digit `b` of unsigned word `v`. -/ +def digit (b v : ℕ) : ℕ := v / 2 ^ b % 2 + +/-- The stable subsequence of words whose digit `b` equals `t`. -/ +def select (b t : ℕ) (xs : List ℕ) : List ℕ := + xs.filter (fun v => decide (digit b v = t)) + +/-- Stable zero/one partition at digit `b`. -/ +def partition (b : ℕ) (xs : List ℕ) : List ℕ := select b 0 xs ++ select b 1 xs + +/-- Apply the first `bits` binary radix passes, starting with the least significant bit. -/ +def passes : ℕ → List ℕ → List ℕ + | 0, xs => xs + | b + 1, xs => partition b (passes b xs) + +/-- Every binary digit is below two. -/ +theorem digit_lt_two (b v : ℕ) : digit b v < 2 := Nat.mod_lt _ (by decide) + +/-- The native shift/remainder implementation reads precisely `digit`. -/ +theorem digit_eq_shift (b v : ℕ) : digit b v = (v >>> b) % 2 := by + rw [Nat.shiftRight_eq_div_pow, digit] + +/-- A selection never increases the length. -/ +theorem length_select_le (b t : ℕ) (xs : List ℕ) : (select b t xs).length ≤ xs.length := + List.length_filter_le _ _ + +/-- A selected word occurred in the original list and has the requested digit. -/ +theorem mem_select {b t v : ℕ} {xs : List ℕ} : + v ∈ select b t xs ↔ v ∈ xs ∧ digit b v = t := by + simp [select] + +/-- Appending an input word appends it to precisely its selected subsequence. -/ +theorem select_append (b t : ℕ) (xs ys : List ℕ) : + select b t (xs ++ ys) = select b t xs ++ select b t ys := by + simp [select] + +/-- A stable bit partition is a permutation of the input, including repeated words. -/ +theorem partition_perm (b : ℕ) (xs : List ℕ) : (partition b xs).Perm xs := by + have h (v : ℕ) : decide (digit b v = 1) = !decide (digit b v = 0) := by + have := digit_lt_two b v + by_cases hz : digit b v = 0 + · simp [hz] + · have ho : digit b v = 1 := by omega + simp [ho] + simpa only [partition, select, h] using + List.filter_append_perm (fun v => decide (digit b v = 0)) xs + +/-- Each pass keeps exactly the same number of words. -/ +theorem length_partition (b : ℕ) (xs : List ℕ) : (partition b xs).length = xs.length := + (partition_perm b xs).length_eq + +/-- All radix passes preserve the input multiset. -/ +theorem passes_perm (bits : ℕ) (xs : List ℕ) : (passes bits xs).Perm xs := by + induction bits with + | zero => exact .refl _ + | succ b ih => exact (partition_perm b _).trans ih + +/-- A pass extends sortedness modulo `2^b` to sortedness modulo `2^(b+1)`. -/ +theorem partition_pairwise {b : ℕ} {xs : List ℕ} + (hs : xs.Pairwise (fun u v => u % 2 ^ b ≤ v % 2 ^ b)) : + (partition b xs).Pairwise (fun u v => u % 2 ^ (b + 1) ≤ v % 2 ^ (b + 1)) := by + have he (v : ℕ) : v % 2 ^ (b + 1) = v % 2 ^ b + 2 ^ b * digit b v := by + rw [Nat.pow_succ, Nat.mod_mul, digit] + have hp (t : ℕ) : + (select b t xs).Pairwise (fun u v => u % 2 ^ (b + 1) ≤ v % 2 ^ (b + 1)) := by + apply (hs.filter _).imp_of_mem + intro u v hu hv huv + change u ∈ select b t xs at hu + change v ∈ select b t xs at hv + rw [mem_select] at hu hv + rw [he, he, hu.2, hv.2] + omega + apply List.pairwise_append.mpr + refine ⟨hp 0, hp 1, ?_⟩ + intro u hu v hv + rw [mem_select] at hu hv + rw [he, he, hu.2, hv.2] + have hmod := Nat.mod_lt u (Nat.two_pow_pos b) + omega + +/-- The output is sorted on the low `bits` bits. -/ +theorem passes_pairwise (bits : ℕ) (xs : List ℕ) : + (passes bits xs).Pairwise (fun u v => u % 2 ^ bits ≤ v % 2 ^ bits) := by + induction bits with + | zero => simp [passes, List.pairwise_iff_getElem, Nat.mod_one] + | succ b ih => exact partition_pairwise ih + +/-- Keys strictly below `2^bits` are sorted in unsigned order after `bits` passes. -/ +theorem passes_correct (bits : ℕ) (xs : List ℕ) + (hbound : ∀ v ∈ xs, v < 2 ^ bits) : + (passes bits xs).Perm xs ∧ (passes bits xs).Pairwise (· ≤ ·) := by + refine ⟨passes_perm bits xs, (passes_pairwise bits xs).imp_of_mem ?_⟩ + intro u v hu hv huv + have hu' := hbound u ((passes_perm bits xs).mem_iff.mp hu) + have hv' := hbound v ((passes_perm bits xs).mem_iff.mp hv) + simpa only [Nat.mod_eq_of_lt hu', Nat.mod_eq_of_lt hv'] using huv + +end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Program.lean b/Complexity/Algorithms/WordRAM/RadixSort/Program.lean new file mode 100755 index 0000000000..a4d8a6f851 --- /dev/null +++ b/Complexity/Algorithms/WordRAM/RadixSort/Program.lean @@ -0,0 +1,204 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.RadixSort.Scan +public import Complexity.Computability.WordRAM.Compile + +/-! +# Fixed-code native unsigned radix sorter + +Two stable scans per bit partition the source into zeros followed by ones in scratch. Swapping +the pointers avoids a copy-back scan. Code and registers are independent of length and key width. +The source pointer on exit identifies the sorted array; the other block is scratch. + +## Main definitions + +* `pass`: two digit scans and pointer swap. +* `code`: the complete radix sorter with a register-supplied bit bound. +* `timeBound`: explicit instruction bound, excluding halt. + +## Main results + +* `exec`: complete charged sorting, with permutation, sortedness and a two-array frame. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.RadixSort + +open Cslib.WordRAM + +/-- Exchange source and scratch, then advance to the next bit. -/ +def advance : Code .standard 16 := ram[ + swap ← src; src ← dst; dst ← swap; bit ← bit + one ] + +/-- Exactly two stable scans per bit. -/ +def pass : Code .standard 16 := + .seq (ram[ out ← #0; tag ← #0 ]) + (ram[ ⟪Scan.code⟫; tag ← #1; ⟪Scan.code⟫; ⟪advance⟫ ]) + +/-- Complete unsigned radix sorter. Parameters are length, two disjoint array bases, and bits. +The returned source register identifies the sorted block. -/ +def code : Code .standard 16 := + .seq (ram[ one ← #1; two ← #2; bit ← #0 ]) (.loop (.lt bit bits) pass) + +/-- All instructions of radix sorting, excluding a final halt. -/ +def timeBound (n keyBits : ℕ) : ℕ := keyBits * (22 * n + 15) + 5 + +private theorem pass_exec {w b s d : ℕ} (xs : List ℕ) (x : Data 16) + (hn : x.regs len = xs.length) (hs : x.regs src = s) (hd : x.regs dst = d) + (hb : x.regs bit = b) (ho : x.regs one = 1) (h2 : x.regs two = 2) + (hmem : Block.HoldsList s xs x.mem) (hwords : ∀ v ∈ xs, Fits w v) + (hdisj : s + xs.length ≤ d ∨ d + xs.length ≤ s) + (hsfit : Fits w (s + xs.length)) (hdfit : Fits w (d + xs.length)) + (hbit : Fits w (b + 1)) (h1 : Fits w 1) : + ∃ y, ExecWithin w pass x y (22 * xs.length + 13) ∧ + Block.HoldsList d (partition b xs) y.mem ∧ + y.regs src = d ∧ y.regs dst = s ∧ y.regs bit = b + 1 ∧ + Modifies (Block.mk d xs.length).footprint x.mem y.mem := by + let x₀ := (x.set out 0).set tag 0 + have hinit : ExecWithin (B := .standard) w (ram[ out ← #0; tag ← #0 ]) x x₀ 2 := by + ram_eval [x₀] + obtain ⟨y, hy, hym, hysrc, hyout, hyf⟩ := Scan.exec (t := 0) xs [] x₀ + (by simpa +decide [x₀] using hn) (by simpa +decide [x₀] using hs) + (by simpa +decide [x₀] using hd) (by simpa +decide [x₀] using hb) + (by simp [x₀]) (by simpa +decide [x₀] using ho) + (by simpa +decide [x₀] using h2) (by simp +decide [x₀]) + hmem (Block.holdsList_nil _ _) hwords + (by simpa using length_select_le b 0 xs) hdisj hsfit hdfit + have hyr (j : Fin 16) (hj : Scan.code.Preserves j) : y.regs j = x₀.regs j := + hy.regs_eq hj + let y₀ := y.set tag 1 + have hn₀ : y₀.regs len = xs.length := by + simpa +decide [y₀, hyr len (by decide), x₀] using hn + have hs₀ : y₀.regs src = s := by + simpa +decide [y₀, hyr src (by decide), x₀] using hs + have hd₀ : y₀.regs dst = d := by + simpa +decide [y₀, hyr dst (by decide), x₀] using hd + have hb₀ : y₀.regs bit = b := by + simpa +decide [y₀, hyr bit (by decide), x₀] using hb + have ho₀ : y₀.regs one = 1 := by + simpa +decide [y₀, hyr one (by decide), x₀] using ho + have h2₀ : y₀.regs two = 2 := by + simpa +decide [y₀, hyr two (by decide), x₀] using h2 + obtain ⟨z, hz, hzm, _, _, hzf⟩ := Scan.exec (t := 1) xs (select b 0 xs) y₀ + hn₀ hs₀ hd₀ hb₀ (by simp [y₀]) ho₀ h2₀ + (by simpa +decide [y₀] using hyout) hysrc + (by simpa only [y₀, Data.set_mem, List.nil_append] using hym) hwords + (by rw [← partition, length_partition]) hdisj hsfit hdfit + have hzs : z.regs src = s := (hz.regs_eq (by decide)).trans hs₀ + have hzd : z.regs dst = d := (hz.regs_eq (by decide)).trans hd₀ + have hzb : z.regs bit = b := (hz.regs_eq (by decide)).trans hb₀ + have hzo : z.regs one = 1 := (hz.regs_eq (by decide)).trans ho₀ + refine ⟨?u, ?run, ?_, ?_, ?_, ?_, ?_⟩ + case run => + have ha : ExecWithin w advance z ?u 4 := by + ram_eval [advance, hzs, hzd, hzb, hzo, hbit] + have hall := hinit.seq (hy.seq ((ExecWithin.const (d := tag) (x := y) h1).seq (hz.seq ha))) + exact hall.mono (by omega) + · exact hzm + · simp +decide [Data.set_regs] + · simp +decide [Data.set_regs] + · simp +decide [Data.set_regs] + · intro a hnot + exact (hzf a hnot).trans (hyf a hnot) + +/-- Sort all keys below `2^keyBits` in `keyBits * (22n+15)+5` instructions. The two arrays +have length `n`, are adjacent, and start at caller-supplied `base`. Only these `2n` cells may +change. The sorted block is returned in `src`; it is one of the two original blocks. This +includes empty lists and zero bit bounds. The key-bit bound and both one-past-end addresses +must fit; constant two must fit. No ordering contract is assumed. -/ +theorem exec {w base keyBits : ℕ} (xs : List ℕ) (x : Data 16) + (hn : x.regs len = xs.length) (hs : x.regs src = base) + (hd : x.regs dst = base + xs.length) (hb : x.regs bits = keyBits) + (hmem : Block.HoldsList base xs x.mem) + (hkeys : ∀ v ∈ xs, v < 2 ^ keyBits) + (hwords : ∀ v ∈ xs, Fits w v) + (hfit : Fits w (base + 2 * xs.length)) (hbits : Fits w keyBits) (h2 : Fits w 2) : + ∃ y, ExecWithin w code x y (timeBound xs.length keyBits) ∧ + Block.HoldsList (y.regs src) (passes keyBits xs) y.mem ∧ + (passes keyBits xs).Perm xs ∧ (passes keyBits xs).Pairwise (· ≤ ·) ∧ + (y.regs src = base ∨ y.regs src = base + xs.length) ∧ + Modifies (Block.mk base (2 * xs.length)).footprint x.mem y.mem := by + let n := xs.length + let I (y : Data 16) : Prop := + y.regs len = n ∧ y.regs bits = keyBits ∧ y.regs one = 1 ∧ y.regs two = 2 ∧ + y.regs bit ≤ keyBits ∧ + ((y.regs src = base ∧ y.regs dst = base + n) ∨ + (y.regs src = base + n ∧ y.regs dst = base)) ∧ + Block.HoldsList (y.regs src) (passes (y.regs bit) xs) y.mem ∧ + Modifies (Block.mk base (2 * n)).footprint x.mem y.mem + have hbody (y : Data 16) (hy : I y) (hc : (Cond.lt bit bits).Holds y) : + ∃ z, ExecWithin w pass y z (22 * n + 13) ∧ I z ∧ + keyBits - z.regs bit < keyBits - y.regs bit := by + obtain ⟨hyn, hyB, hyo, hy2, hyb, hyp, hym, hyf⟩ := hy + have hlt : y.regs bit < keyBits := by simpa only [Cond.holds_lt, hyB] using hc + have hl : (passes (y.regs bit) xs).length = n := (passes_perm _ _).length_eq + have hdis : y.regs src + n ≤ y.regs dst ∨ y.regs dst + n ≤ y.regs src := by + rcases hyp with ⟨ha, hb⟩ | ⟨ha, hb⟩ <;> omega + have hsf : Fits w (y.regs src + n) := hfit.of_le (by + rcases hyp with ⟨ha, hb⟩ | ⟨ha, hb⟩ <;> omega) + have hdf : Fits w (y.regs dst + n) := hfit.of_le (by + rcases hyp with ⟨ha, hb⟩ | ⟨ha, hb⟩ <;> omega) + have hwds : ∀ v ∈ passes (y.regs bit) xs, Fits w v := fun v hv => + hwords v ((passes_perm _ _).mem_iff.mp hv) + obtain ⟨z, he, hzm, hzs, hzd, hzb, hzf⟩ := pass_exec (passes (y.regs bit) xs) y + (hyn.trans hl.symm) rfl rfl rfl hyo hy2 hym hwds + (by simpa only [hl] using hdis) (by simpa only [hl] using hsf) + (by simpa only [hl] using hdf) (hbits.of_le (by omega)) (h2.of_le (by omega)) + refine ⟨z, by simpa only [hl] using he, ?_, by omega⟩ + refine ⟨(he.regs_eq (by decide)).trans hyn, (he.regs_eq (by decide)).trans hyB, + (he.regs_eq (by decide)).trans hyo, (he.regs_eq (by decide)).trans hy2, + by omega, ?_, ?_, ?_⟩ + · rw [hzs, hzd] + rcases hyp with ⟨ha, hb⟩ | ⟨ha, hb⟩ + · exact Or.inr ⟨hb, ha⟩ + · exact Or.inl ⟨hb, ha⟩ + · rw [hzs, hzb, passes] + exact hzm + · intro a hnot + have hnot' : a ∉ (Block.mk (y.regs dst) (passes (y.regs bit) xs).length).footprint := by + intro ha + apply hnot + rw [Block.mem_footprint] at ha ⊢ + change y.regs dst ≤ a ∧ + a < y.regs dst + (passes (y.regs bit) xs).length at ha + change base ≤ a ∧ a < base + 2 * n + rw [hl] at ha + rcases hyp with ⟨ha', hb'⟩ | ⟨ha', hb'⟩ <;> omega + exact (hzf a hnot').trans (hyf a hnot) + let x₀ := ((x.set one 1).set two 2).set bit 0 + have hx₀ : I x₀ := by + refine ⟨?_, ?_, ?_, ?_, ?_, ?_, ?_, Modifies.refl _ _⟩ + · simpa +decide [x₀] using hn + · simpa +decide [x₀] using hb + · simp +decide [x₀] + · simp +decide [x₀] + · simp [x₀] + · left; simpa +decide [x₀, n] using And.intro hs hd + · simpa +decide [x₀, passes, hs] using hmem + obtain ⟨y, he, hy, hc⟩ := ExecWithin.loop_invariant I + (fun y => keyBits - y.regs bit) hbody hx₀ + obtain ⟨hyn, hyB, hyo, hy2, hyb, hyp, hym, hyf⟩ := hy + have hend : y.regs bit = keyBits := by + simp only [Cond.holds_lt, hyB] at hc + omega + have hinit : ExecWithin (B := .standard) w (ram[ one ← #1; two ← #2; bit ← #0 ]) x x₀ 3 := by + ram_eval [x₀, h2, h2.of_le (show 1 ≤ 2 by omega)] + refine ⟨y, ?_, by simpa only [hend] using hym, + (passes_correct keyBits xs hkeys).1, (passes_correct keyBits xs hkeys).2, + hyp.imp And.left And.left, hyf⟩ + have hall := hinit.seq he + convert hall using 1 + · rfl + · simp +decide only [timeBound, n, x₀, Data.set_regs, Function.update_apply, + ↓reduceIte, Nat.sub_zero, Nat.add_assoc, Nat.reduceAdd] + omega + +end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Scan.lean b/Complexity/Algorithms/WordRAM/RadixSort/Scan.lean new file mode 100755 index 0000000000..8ebb617df9 --- /dev/null +++ b/Complexity/Algorithms/WordRAM/RadixSort/Scan.lean @@ -0,0 +1,254 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Algorithms.WordRAM.RadixSort.Model +public import Complexity.Computability.WordRAM.ListArray + +/-! +# Native stable digit-selection scan + +The scan appends a stable digit-selected subsequence to a disjoint destination array. Both +orders of the two disjoint blocks are allowed, so consecutive radix passes can swap pointers. +All loads, address calculations, stores and loop control are native instructions. + +## Main definitions + +* `Scan.body`: read, test, conditionally append, and advance. +* `Scan.code`: reset the input cursor and scan the array. + +## Main results + +* `Scan.exec`: stable selection with a linear instruction bound and memory frame. +-/ + +@[expose] public section + +namespace Cslib.Algorithms.WordRAM.RadixSort + +open Cslib.WordRAM + +/-- Array length. -/ +abbrev len : Fin 16 := 0 +/-- Source pointer; the sorter returns the final sorted pointer here. -/ +abbrev src : Fin 16 := 1 +/-- Scratch pointer. -/ +abbrev dst : Fin 16 := 2 +/-- Current bit index. -/ +abbrev bit : Fin 16 := 3 +/-- Digit selected by the current scan. -/ +abbrev tag : Fin 16 := 4 +/-- Constant one. -/ +abbrev one : Fin 16 := 5 +/-- Constant two. -/ +abbrev two : Fin 16 := 6 +/-- Input cursor. -/ +abbrev idx : Fin 16 := 7 +/-- Output cursor. -/ +abbrev out : Fin 16 := 8 +/-- Address temporary. -/ +abbrev addr : Fin 16 := 9 +/-- Loaded key. -/ +abbrev value : Fin 16 := 10 +/-- Current key digit. -/ +abbrev keyDigit : Fin 16 := 11 +/-- Number of key bits, supplied by the caller. -/ +abbrev bits : Fin 16 := 12 +/-- Temporary used to swap the two pointers. -/ +abbrev swap : Fin 16 := 13 + +namespace Scan + +/-- One stable-selection iteration, with every address and counter operation charged. -/ +def body : Code .standard 16 := ram[ + addr ← src + idx; value ← *addr; + keyDigit ← value >>> bit; keyDigit ← keyDigit % two; + if keyDigit = tag then [ addr ← dst + out; *addr ← value; out ← out + one ] + else [ addr ← addr ]; + idx ← idx + one ] + +/-- Reset the input cursor, then perform one stable-selection scan. -/ +def code : Code .standard 16 := ram[ idx ← #0; while idx < len do [ ⟪body⟫ ] ] + +private theorem preserves_low (j : Fin 16) (hj : j.val ≤ 6) : body.Preserves j := by + have hne (k : Fin 16) (hk : 7 ≤ k.val) : (some k : Option (Fin 16)) ≠ some j := by + intro h + have hv := congrArg Fin.val (Option.some.inj h) + omega + change (some addr ≠ some j) ∧ (some value ≠ some j) ∧ + (some keyDigit ≠ some j) ∧ (some keyDigit ≠ some j) ∧ + (((some addr ≠ some j) ∧ (none ≠ some j) ∧ (some out ≠ some j)) ∧ + (some addr ≠ some j)) ∧ (some idx ≠ some j) + have hnone : (none : Option (Fin 16)) ≠ some j := by intro h; cases h + exact ⟨hne addr (by decide), hne value (by decide), hne keyDigit (by decide), + hne keyDigit (by decide), + ⟨⟨hne addr (by decide), hnone, hne out (by decide)⟩, + hne addr (by decide)⟩, + hne idx (by decide)⟩ + +private theorem step {w : ℕ} (x : Data 16) + (ho : x.regs one = 1) (ht : x.regs two = 2) + (hs : Fits w (x.regs src + x.regs idx)) + (hd : Fits w (x.regs dst + x.regs out)) + (hi : Fits w (x.regs idx + 1)) + (hu : digit (x.regs bit) (x.mem (x.regs src + x.regs idx)) = x.regs tag → + Fits w (x.regs out + 1)) + (hv : Fits w (x.mem (x.regs src + x.regs idx))) : + ∃ y, ExecWithin w body x y 9 ∧ + y.regs idx = x.regs idx + 1 ∧ + y.regs out = x.regs out + + (if digit (x.regs bit) (x.mem (x.regs src + x.regs idx)) = x.regs tag then 1 else 0) ∧ + y.mem = (if digit (x.regs bit) (x.mem (x.regs src + x.regs idx)) = x.regs tag then + Function.update x.mem (x.regs dst + x.regs out) (x.mem (x.regs src + x.regs idx)) + else x.mem) := by + have hshift := hv.of_le (Nat.shiftRight_le _ (x.regs bit)) + have hmod := hshift.of_le (Nat.mod_le _ 2) + by_cases h : digit (x.regs bit) (x.mem (x.regs src + x.regs idx)) = x.regs tag + · have heq : (x.mem (x.regs src + x.regs idx) >>> x.regs bit) % 2 = x.regs tag := + (digit_eq_shift _ _).symm.trans h + have hinc := hu h + have htagfit : Fits w (x.regs tag) := heq ▸ hmod + refine ⟨?_, ?run, ?_, ?_, ?_⟩ + case run => ram_eval [body, ho, ht, hs, hd, hi, hinc, hshift, hmod, heq, htagfit] + all_goals simp +decide [Data.set_regs, h] + · have hne : (x.mem (x.regs src + x.regs idx) >>> x.regs bit) % 2 ≠ x.regs tag := by + simpa only [← digit_eq_shift] using h + refine ⟨?_, ?run, ?_, ?_, ?_⟩ + case run => ram_eval [body, ho, ht, hs, hd, hi, hshift, hmod, hne] + all_goals simp +decide [Data.set_regs, h] + +/-- Append a stable digit-selected subsequence in `11n+3` instructions. The destination has +capacity `n`; source and destination blocks are disjoint, in either order. The caller's existing +destination prefix is retained. Only destination cells may change. -/ +theorem exec {w s d b t : ℕ} (xs acc : List ℕ) (x : Data 16) + (hn : x.regs len = xs.length) (hs : x.regs src = s) (hd : x.regs dst = d) + (hb : x.regs bit = b) (ht : x.regs tag = t) + (ho : x.regs one = 1) (h2 : x.regs two = 2) (hout : x.regs out = acc.length) + (hsrc : Block.HoldsList s xs x.mem) (hacc : Block.HoldsList d acc x.mem) + (hwords : ∀ v ∈ xs, Fits w v) + (hcap : (acc ++ select b t xs).length ≤ xs.length) + (hdisj : s + xs.length ≤ d ∨ d + xs.length ≤ s) + (hsfit : Fits w (s + xs.length)) (hdfit : Fits w (d + xs.length)) : + ∃ y, ExecWithin w code x y (11 * xs.length + 3) ∧ + Block.HoldsList d (acc ++ select b t xs) y.mem ∧ + Block.HoldsList s xs y.mem ∧ + y.regs out = (acc ++ select b t xs).length ∧ + Modifies (Block.mk d xs.length).footprint x.mem y.mem := by + let selected (i : ℕ) := acc ++ select b t (xs.take i) + have hselected (i : ℕ) : (selected i).length ≤ xs.length := by + have hsub : (select b t (xs.take i)).Sublist (select b t xs) := + (List.take_sublist i xs).filter _ + simp only [selected, List.length_append] at hcap ⊢ + have := hsub.length_le + omega + let I (y : Data 16) : Prop := + (∀ j : Fin 16, j.val ≤ 6 → y.regs j = x.regs j) ∧ + y.regs idx ≤ xs.length ∧ y.regs out = (selected (y.regs idx)).length ∧ + Block.HoldsList s xs y.mem ∧ Block.HoldsList d (selected (y.regs idx)) y.mem ∧ + Modifies (Block.mk d xs.length).footprint x.mem y.mem + have hbody (y : Data 16) (hy : I y) (hc : (Cond.lt idx len).Holds y) : + ∃ z, ExecWithin w body y z 9 ∧ I z ∧ + xs.length - z.regs idx < xs.length - y.regs idx := by + obtain ⟨hr, hle, hu, hm, ha, hf⟩ := hy + have hyl : y.regs len = xs.length := (hr len (by decide)).trans hn + have hys : y.regs src = s := (hr src (by decide)).trans hs + have hyd : y.regs dst = d := (hr dst (by decide)).trans hd + have hyb : y.regs bit = b := (hr bit (by decide)).trans hb + have hyt : y.regs tag = t := (hr tag (by decide)).trans ht + have hyo : y.regs one = 1 := (hr one (by decide)).trans ho + have hy2 : y.regs two = 2 := (hr two (by decide)).trans h2 + have hi : y.regs idx < xs.length := by simpa only [Cond.holds_lt, hyl] using hc + have hread : y.mem (y.regs src + y.regs idx) = xs[y.regs idx] := by + rw [hys]; exact hm.get _ hi + have hv : Fits w (y.mem (y.regs src + y.regs idx)) := by + rw [hread]; exact hwords _ (List.getElem_mem hi) + have hsa : Fits w (y.regs src + y.regs idx) := by + rw [hys]; exact hsfit.of_le (by omega) + have hda : Fits w (y.regs dst + y.regs out) := by + rw [hyd, hu]; exact hdfit.of_le (by have := hselected (y.regs idx); omega) + have hic : Fits w (y.regs idx + 1) := hsfit.of_le (by omega) + have hsuc : selected (y.regs idx + 1) = selected (y.regs idx) ++ + (if digit b xs[y.regs idx] = t then [xs[y.regs idx]] else []) := by + dsimp only [selected] + rw [List.take_succ_eq_append_getElem hi, select_append] + simp only [select, List.filter_cons, List.filter_nil, decide_eq_true_eq] + split_ifs <;> simp only [List.append_assoc, List.append_nil] + by_cases heq : digit b xs[y.regs idx] = t + · have hlt : y.regs out < xs.length := by + have hlim := hselected (y.regs idx + 1) + rw [hsuc, ite_eq_left heq, List.length_append, List.length_singleton] at hlim + omega + obtain ⟨z, he, hzi, hzu, hzm⟩ := step y hyo hy2 hsa hda hic + (fun _ => hdfit.of_le (by omega)) hv + rw [hyb, hread, hyt, ite_eq_left heq] at hzu hzm + refine ⟨z, he, ?_, by omega⟩ + refine ⟨?_, by omega, ?_, ?_, ?_, ?_⟩ + · intro j hj + exact (he.regs_eq (preserves_low j hj)).trans (hr j hj) + · rw [hzi, hsuc, ite_eq_left heq, List.length_append, List.length_singleton] + omega + · apply hm.frame + intro a hlo hhi + rw [hzm, hyd, Function.update_of_ne (by rcases hdisj with h | h <;> omega)] + · rw [hzi, hsuc, ite_eq_left heq, hzm, hyd, hu] + exact ha.append_one _ + · intro a hnot + have hne : a ≠ d + y.regs out := by + intro heqa + apply hnot + apply Block.mem_footprint.mpr + change d ≤ a ∧ a < d + xs.length + omega + rw [hzm, hyd, Function.update_of_ne hne] + exact hf a hnot + · obtain ⟨z, he, hzi, hzu, hzm⟩ := step y hyo hy2 hsa hda hic + (by simpa only [hyb, hread, hyt] using fun h => (heq h).elim) hv + rw [hyb, hread, hyt, ite_eq_right heq] at hzu hzm + have hsel : selected (y.regs idx + 1) = selected (y.regs idx) := by + rw [hsuc, ite_eq_right heq, List.append_nil] + refine ⟨z, he, ?_, by omega⟩ + refine ⟨?_, by omega, ?_, ?_, ?_, ?_⟩ + · intro j hj + exact (he.regs_eq (preserves_low j hj)).trans (hr j hj) + · rw [hzi, hsel, hzu, Nat.add_zero, hu] + · simpa only [hzm] using hm + · simpa only [hzi, hsel, hzm] using ha + · simpa only [hzm] using hf + let x₀ := x.set idx 0 + have hx₀ : I x₀ := by + refine ⟨?_, ?_, ?_, ?_, ?_, ?_⟩ + · intro j hj + have hne : j ≠ idx := by + intro h + have hv := congrArg Fin.val h + change j.val = 7 at hv + omega + exact Function.update_of_ne hne _ _ + · simp [x₀] + · simpa +decide [x₀, selected, select] using hout + · exact hsrc + · simpa [x₀, selected, select] using hacc + · exact Modifies.refl _ _ + obtain ⟨y, he, hy, hc⟩ := ExecWithin.loop_invariant I + (fun y => xs.length - y.regs idx) hbody hx₀ + obtain ⟨hr, hle, hu, hm, ha, hf⟩ := hy + have hend : y.regs idx = xs.length := by + have hyl : y.regs len = xs.length := (hr len (by decide)).trans hn + simp only [Cond.holds_lt, hyl] at hc + omega + have hi := (ExecWithin.const (d := idx) (x := x) (fits_zero w)).seq he + refine ⟨y, ?_, ?_, hm, ?_, hf⟩ + · convert hi using 1 + · rfl + · simp only [x₀, Data.set_regs, Function.update_self, Nat.sub_zero] + omega + · simpa only [hend, selected, List.take_length] using ha + · simpa only [hend, selected, List.take_length] using hu + +end Scan +end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Computability/WordRAM/ListArray.lean b/Complexity/Computability/WordRAM/ListArray.lean new file mode 100755 index 0000000000..5d90770831 --- /dev/null +++ b/Complexity/Computability/WordRAM/ListArray.lean @@ -0,0 +1,94 @@ +/- +Copyright (c) 2026 Erik Demaine. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: Erik Demaine +-/ + +module +-- upstream: candidate + +public import Complexity.Computability.WordRAM.Heap + +/-! +# List-valued block assertions + +These are specializations of `Block.Realizes`, not an additional memory model. + +## Main definitions + +* `Block.HoldsList`: a contiguous block realizes a list. + +## Main results + +* `Block.HoldsList.append_one`: a store immediately after a list appends one element. +* `Block.HoldsList.frame`: writes outside the block preserve its list. +-/ + +@[expose] public section + +namespace Cslib.WordRAM.Block + +/-- The block starting at `base` realizes exactly the entries of `xs`. -/ +abbrev HoldsList (base : ℕ) (xs : List ℕ) (mem : ℕ → ℕ) : Prop := + (Block.mk base xs.length).Realizes mem (fun i => xs[i]) + +/-- Read an entry of a represented list. -/ +theorem HoldsList.get {base : ℕ} {xs : List ℕ} {mem : ℕ → ℕ} + (h : HoldsList base xs mem) (i : ℕ) (hi : i < xs.length) : + mem (base + i) = xs[i] := h ⟨i, hi⟩ + +/-- An empty list imposes no memory constraints. -/ +theorem holdsList_nil (base : ℕ) (mem : ℕ → ℕ) : HoldsList base [] mem := by + intro i + exact Fin.elim0 i + +/-- The prefix of a represented list is represented in the same memory. -/ +theorem HoldsList.take {base : ℕ} {xs : List ℕ} {mem : ℕ → ℕ} + (h : HoldsList base xs mem) (k : ℕ) : HoldsList base (xs.take k) mem := by + intro i + have hi : i.val < xs.length := lt_of_lt_of_le i.isLt (List.take_sublist k xs).length_le + change mem (base + i.val) = (xs.take k)[i.val] + rw [List.getElem_take] + exact h.get i.val hi + +/-- A store immediately following a represented list appends an element. -/ +theorem HoldsList.append_one {base : ℕ} {xs : List ℕ} {mem : ℕ → ℕ} + (h : HoldsList base xs mem) (v : ℕ) : + HoldsList base (xs ++ [v]) (Function.update mem (base + xs.length) v) := by + intro i + change Function.update mem (base + xs.length) v (base + i.val) = (xs ++ [v])[i.val] + by_cases hi : i.val < xs.length + · rw [Function.update_of_ne (by omega), List.getElem_append_left hi] + exact h.get i.val hi + · have hil := i.isLt + simp only [List.length_append, List.length_singleton] at hil + have he : i.val = xs.length := by omega + simp [he] + +/-- Writes outside the represented block preserve the list. -/ +theorem HoldsList.frame {base : ℕ} {xs : List ℕ} {mem mem' : ℕ → ℕ} + (h : HoldsList base xs mem) + (heq : ∀ a, base ≤ a → a < base + xs.length → mem' a = mem a) : + HoldsList base xs mem' := by + intro i + have hi : i.val < xs.length := i.isLt + change mem' (base + i.val) = xs[i.val] + rw [heq _ (by omega) (by omega)] + exact h.get i.val i.isLt + +/-- Append two adjacent represented blocks. -/ +theorem HoldsList.append {base : ℕ} {xs ys : List ℕ} {mem : ℕ → ℕ} + (h : HoldsList base xs mem) (h' : HoldsList (base + xs.length) ys mem) : + HoldsList base (xs ++ ys) mem := by + intro i + change mem (base + i.val) = (xs ++ ys)[i.val] + by_cases hi : i.val < xs.length + · rw [List.getElem_append_left hi] + exact h.get i.val hi + · rw [List.getElem_append_right (by omega)] + have hil := i.isLt + simp only [List.length_append] at hil + have he : base + xs.length + (i.val - xs.length) = base + i.val := by omega + simpa only [he] using h'.get (i.val - xs.length) (by omega) + +end Cslib.WordRAM.Block diff --git a/ComplexityTests/WordRAMRadixSort.lean b/ComplexityTests/WordRAMRadixSort.lean new file mode 100755 index 0000000000..55ddb07ef7 --- /dev/null +++ b/ComplexityTests/WordRAMRadixSort.lean @@ -0,0 +1,67 @@ +import Complexity.Algorithms.WordRAM.RadixSort + +open Cslib.WordRAM +open Cslib.Algorithms.WordRAM + +private def observe (w fuel : ℕ) (xs : List ℕ) : Option (List ℕ × List ℕ × ℕ) := + (Code.eval w fuel RadixSort.inputCode (Data.init xs)).map fun (y, _) => + let m := y.regs Deduplicate.out + ((List.range m).map (fun i => y.mem (1 + 2 * xs.length + i)), + (List.range m).map (fun i => y.mem (1 + 3 * xs.length + i)), m) + +-- Empty input still executes setup and exits both array scans safely. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 3 4 [] = some ([], [], 0)) + +-- The largest fitting four-bit canonical signed encoding occurs at the key boundary. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 4 5 [Signed.encode (-7)] = some ([15], [1], 1)) + +-- Zero is canonical; a repeated zero is one key with full multiplicity two. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 4 5 [0, 0] = some ([0], [2], 1)) + +-- Counts are full multiplicities, including values above three. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 5 6 (List.replicate 4 (Signed.encode (-1))) = some ([3], [4], 1)) + +-- Positive and negative input integers, repeated keys, and zero in unsigned encoding order. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 5 6 (([2, -1, 2, 0, -1] : List ℤ).map Signed.encode) = + some ([0, 3, 4], [1, 2, 2], 3)) + +-- Encoding order deliberately differs from signed integer order. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 4 5 (([-2, 1, -1] : List ℤ).map Signed.encode) = + some ([2, 3, 5], [1, 1, 1], 3)) + +-- Two physical bits cannot address this singleton's count output: layout fit matters. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 2 4 [2] = none) + +-- A raw unsigned input beyond the declared physical word bound faults as well. +/-- info: true -/ +#guard_msgs in +#eval decide (observe 2 4 [4] = none) + +-- The assembled program is executed using the ordinary native initial state. +/-- info: true -/ +#guard_msgs in +#eval decide ((run 4 RadixSort.program 500 (State.init [4, 3, 3])).map Prod.fst = some 2) + +-- The array layout threshold includes zero and small input lengths. +example : RadixSort.layoutWidth 0 = 3 ∧ RadixSort.layoutWidth 2 = 6 := by decide +kernel + +-- A bit bound of zero is meaningful exactly when all unsigned keys are zero. +example : RadixSort.passes 0 [0, 0] = [0, 0] ∧ + ¬([2, 0] : List ℕ).Pairwise (· ≤ ·) := by decide +kernel + +-- One finite instruction list is used at every input length and word width. +example : RadixSort.program.length = 81 := by decide +kernel From 4edf0b8ac9bba87b0f282fec989670c49b8bb641 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 08:00:31 -0400 Subject: [PATCH 07/11] chore(wordram): keep Lean sources non-executable Normalize the new file modes. Lean source bytes are unchanged from the fully validated 25474f432. --- Complexity/Algorithms/WordRAM/Deduplicate/Count.lean | 0 Complexity/Algorithms/WordRAM/Deduplicate/Model.lean | 0 Complexity/Algorithms/WordRAM/Deduplicate/Program.lean | 0 Complexity/Algorithms/WordRAM/RadixSort.lean | 0 Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean | 0 Complexity/Algorithms/WordRAM/RadixSort/Model.lean | 0 Complexity/Algorithms/WordRAM/RadixSort/Program.lean | 0 Complexity/Algorithms/WordRAM/RadixSort/Scan.lean | 0 Complexity/Computability/WordRAM/ListArray.lean | 0 ComplexityTests/WordRAMRadixSort.lean | 0 10 files changed, 0 insertions(+), 0 deletions(-) mode change 100755 => 100644 Complexity/Algorithms/WordRAM/Deduplicate/Count.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/Deduplicate/Model.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/Deduplicate/Program.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/RadixSort.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/RadixSort/Model.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/RadixSort/Program.lean mode change 100755 => 100644 Complexity/Algorithms/WordRAM/RadixSort/Scan.lean mode change 100755 => 100644 Complexity/Computability/WordRAM/ListArray.lean mode change 100755 => 100644 ComplexityTests/WordRAMRadixSort.lean diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Count.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Count.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Model.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Model.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/Deduplicate/Program.lean b/Complexity/Algorithms/WordRAM/Deduplicate/Program.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/RadixSort.lean b/Complexity/Algorithms/WordRAM/RadixSort.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean b/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Model.lean b/Complexity/Algorithms/WordRAM/RadixSort/Model.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Program.lean b/Complexity/Algorithms/WordRAM/RadixSort/Program.lean old mode 100755 new mode 100644 diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Scan.lean b/Complexity/Algorithms/WordRAM/RadixSort/Scan.lean old mode 100755 new mode 100644 diff --git a/Complexity/Computability/WordRAM/ListArray.lean b/Complexity/Computability/WordRAM/ListArray.lean old mode 100755 new mode 100644 diff --git a/ComplexityTests/WordRAMRadixSort.lean b/ComplexityTests/WordRAMRadixSort.lean old mode 100755 new mode 100644 From b2554d58b9eb90291e7d9314128433b62016e203 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 12:03:10 -0400 Subject: [PATCH 08/11] fix(wordram): export sorted radix keys from composed contracts Return strict unsigned key order already proved by RunEncoding from both composed APIs. Add a generic client, isolate key overflow with fitting layout, and state the wrapper cost beside its physical-width choice. Predicate reuse migration requested by review is tracked in #1509. Validation: all six scripts/check.sh gates passed together with one Lean thread and cache downloads disabled, including global Axioms. Changed-file lints against fd69be7bc, hygiene --check-simp, PR names and diff checks pass. No curated abstract impact. Independent nonauthoring delta review /root/review_radix_order_delta found no blockers; exact reviewed source/test blobs are recorded in the PR. Outside coordinator review remains required. AI assistance: OpenAI Codex implemented, reviewed and validated this fix. --- Complexity/Algorithms/WordRAM/RadixSort.lean | 6 ++++-- .../Algorithms/WordRAM/RadixSort/Deduplicate.lean | 3 ++- ComplexityTests/WordRAMRadixSort.lean | 12 ++++++++++-- 3 files changed, 16 insertions(+), 5 deletions(-) diff --git a/Complexity/Algorithms/WordRAM/RadixSort.lean b/Complexity/Algorithms/WordRAM/RadixSort.lean index d3a1a6584c..a7409b44cd 100644 --- a/Complexity/Algorithms/WordRAM/RadixSort.lean +++ b/Complexity/Algorithms/WordRAM/RadixSort.lean @@ -16,6 +16,7 @@ public import Complexity.Computability.WordRAM.Time The input is the ordinary native list of unsigned words. The code reads its length and the physical word width, and uses that width as an explicit key-bit bound. The lower-level `deduplicateCode` accepts a smaller caller-supplied bit bound when known. +The standard-input wrapper's charged bound is `O(n * w)` for positive `n` and `w`. The key array starts at `1+2n`, the count array at `1+3n`; the return value is their common logical length. The two sorting buffers occupy `[1,1+2n)`. The implementation accepts repeated @@ -78,6 +79,7 @@ theorem input_exec {w : ℕ} (xs : List ℕ) (hinput : InputFits w xs) Block.HoldsList (1 + 2 * xs.length) (rs.map Prod.fst) y.mem ∧ Block.HoldsList (1 + 3 * xs.length) (rs.map Prod.snd) y.mem ∧ (rs.map Prod.fst).Nodup ∧ + (rs.map Prod.fst).Pairwise (· < ·) ∧ (∀ k, k ∈ rs.map Prod.fst ↔ k ∈ xs) ∧ (∀ p ∈ rs, 0 < p.2 ∧ xs.count p.1 = p.2) ∧ Modifies (Block.mk 1 (4 * xs.length)).footprint (Data.init xs : Data 16).mem y.mem := by @@ -97,10 +99,10 @@ theorem input_exec {w : ℕ} (xs : List ℕ) (hinput : InputFits w xs) change (Data.init xs : Data 16).mem (1 + i.val) = xs[i.val] rw [Nat.add_comm, Data.init_mem_succ] simp only [List.getD, List.getElem?_eq_getElem i.isLt, Option.getD_some] - obtain ⟨y, rs, he, hyout, hyrk, hyrc, hnodup, hset, hcount, hframe⟩ := + obtain ⟨y, rs, he, hyout, hyrk, hyrc, hnodup, hordered, hset, hcount, hframe⟩ := deduplicate_exec xs u hn hs hd hb hm hinput.2 hinput.2 hfit (fits_self w) h2 obtain ⟨steps, hsteps, hexec⟩ := hi.seq he - refine ⟨y, rs, steps, hexec, ?_, by omega, hyrk, hyrc, hnodup, hset, hcount, hframe⟩ + refine ⟨y, rs, steps, hexec, ?_, by omega, hyrk, hyrc, hnodup, hordered, hset, hcount, hframe⟩ simpa only [program, inputCode, hyout, State.init, x] using hexec.executes Deduplicate.out end Cslib.Algorithms.WordRAM.RadixSort diff --git a/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean b/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean index eed159b431..49d3ce5bb2 100644 --- a/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean +++ b/Complexity/Algorithms/WordRAM/RadixSort/Deduplicate.lean @@ -66,6 +66,7 @@ theorem deduplicate_exec {w base keyBits : ℕ} (xs : List ℕ) (x : Data 16) Block.HoldsList (base + 2 * xs.length) (rs.map Prod.fst) y.mem ∧ Block.HoldsList (base + 3 * xs.length) (rs.map Prod.snd) y.mem ∧ (rs.map Prod.fst).Nodup ∧ + (rs.map Prod.fst).Pairwise (· < ·) ∧ (∀ k, k ∈ rs.map Prod.fst ↔ k ∈ xs) ∧ (∀ p ∈ rs, 0 < p.2 ∧ xs.count p.1 = p.2) ∧ Modifies (Block.mk base (4 * xs.length)).footprint x.mem y.mem := by @@ -106,7 +107,7 @@ theorem deduplicate_exec {w base keyBits : ℕ} (xs : List ℕ) (x : Data 16) (by rw [hlen]; exact hfit.of_le (by omega)) (h2.of_le (by omega)) have heq : base + 2 * xs.length + xs.length = base + 3 * xs.length := by omega refine ⟨y, rs, ?_, hyout, hyrk, by simpa only [hlen, heq] using hyrc, - hrs.keys_nodup, ?_, ?_, ?_⟩ + hrs.keys_nodup, List.pairwise_map.mpr hrs.ordered, ?_, ?_, ?_⟩ · have hall := hinit.seq (hu.seq (hl.seq hy)) apply hall.mono simp only [deduplicateTime, timeBound, hlen] diff --git a/ComplexityTests/WordRAMRadixSort.lean b/ComplexityTests/WordRAMRadixSort.lean index 55ddb07ef7..06b443a9ad 100644 --- a/ComplexityTests/WordRAMRadixSort.lean +++ b/ComplexityTests/WordRAMRadixSort.lean @@ -46,10 +46,10 @@ private def observe (w fuel : ℕ) (xs : List ℕ) : Option (List ℕ × List #guard_msgs in #eval decide (observe 2 4 [2] = none) --- A raw unsigned input beyond the declared physical word bound faults as well. +-- Key overflow faults even when the singleton's four-buffer layout fits. /-- info: true -/ #guard_msgs in -#eval decide (observe 2 4 [4] = none) +#eval decide (observe 4 5 [16] = none) -- The assembled program is executed using the ordinary native initial state. /-- info: true -/ @@ -65,3 +65,11 @@ example : RadixSort.passes 0 [0, 0] = [0, 0] ∧ -- One finite instruction list is used at every input length and word width. example : RadixSort.program.length = 81 := by decide +kernel + +-- Clients recover strict unsigned key order directly from the complete input contract. +example (w : ℕ) (xs : List ℕ) (hinput : InputFits w xs) + (hfit : Fits w (1 + 4 * xs.length)) (h2 : Fits w 2) : + ∃ keys : List ℕ, keys.Pairwise (· < ·) ∧ ∀ k, k ∈ keys ↔ k ∈ xs := by + obtain ⟨_, rs, _, _, _, _, _, _, _, hordered, hset, _, _⟩ := + RadixSort.input_exec xs hinput hfit h2 + exact ⟨rs.map Prod.fst, hordered, hset⟩ From ec8aa561feae957effa539ac6f84d64c5db903c8 Mon Sep 17 00:00:00 2001 From: Erik Demaine Date: Wed, 28 Oct 2026 15:12:47 -0400 Subject: [PATCH 09/11] test(wordram): drop disconnected radix contract example Address the nonblocking #1465 review finding by deleting the generic sorted-set existential, which did not retain an execution or memory witness. Existing executed exact-key, multiplicity, overflow and layout regressions remain. Native code and public library contracts are unchanged. Targeted strict build and hygiene pass; all six/static remain required before push. AI assistance: OpenAI Codex. --- ComplexityTests/WordRAMRadixSort.lean | 8 -------- 1 file changed, 8 deletions(-) diff --git a/ComplexityTests/WordRAMRadixSort.lean b/ComplexityTests/WordRAMRadixSort.lean index 06b443a9ad..28c8c96cd4 100644 --- a/ComplexityTests/WordRAMRadixSort.lean +++ b/ComplexityTests/WordRAMRadixSort.lean @@ -65,11 +65,3 @@ example : RadixSort.passes 0 [0, 0] = [0, 0] ∧ -- One finite instruction list is used at every input length and word width. example : RadixSort.program.length = 81 := by decide +kernel - --- Clients recover strict unsigned key order directly from the complete input contract. -example (w : ℕ) (xs : List ℕ) (hinput : InputFits w xs) - (hfit : Fits w (1 + 4 * xs.length)) (h2 : Fits w 2) : - ∃ keys : List ℕ, keys.Pairwise (· < ·) ∧ ∀ k, k ∈ keys ↔ k ∈ xs := by - obtain ⟨_, rs, _, _, _, _, _, _, _, hordered, hset, _, _⟩ := - RadixSort.input_exec xs hinput hfit h2 - exact ⟨rs.map Prod.fst, hordered, hset⟩ From c0a41fc9ceb11dc67ae3e4a4f0e97de96641f946 Mon Sep 17 00:00:00 2001 From: Samuel Giuseppe Schlesinger Date: Wed, 28 Oct 2026 22:22:09 -0400 Subject: [PATCH 10/11] fix(imports): deduplicate aggregate after hashing integration --- Complexity.lean | 1 - 1 file changed, 1 deletion(-) diff --git a/Complexity.lean b/Complexity.lean index d05051bb3a..a6efd4ca55 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -4,7 +4,6 @@ public import Complexity.Algorithms.ExactTriangle.Defs public import Complexity.Algorithms.ExactTriangle.Hashing public import Complexity.Algorithms.ExactTriangle.Lopsided public import Complexity.Algorithms.ExactTriangle.Weights -public import Complexity.Algorithms.ExactTriangle.Lopsided public import Complexity.Algorithms.Lean.FusionTree.Capacity public import Complexity.Algorithms.Lean.FusionTree.Compression public import Complexity.Algorithms.Lean.FusionTree.Correction From a1295d519c4494f72a9d8a5ca6af28e734fc219d Mon Sep 17 00:00:00 2001 From: Samuel Giuseppe Schlesinger Date: Wed, 28 Oct 2026 22:52:34 -0400 Subject: [PATCH 11/11] fix(imports): regenerate aggregate after radix-sort integration --- Complexity.lean | 5 ----- 1 file changed, 5 deletions(-) diff --git a/Complexity.lean b/Complexity.lean index e29e03a954..bf3ca017d2 100644 --- a/Complexity.lean +++ b/Complexity.lean @@ -53,11 +53,6 @@ public import Complexity.Algorithms.ThreeSum.ThreeSets public import Complexity.Algorithms.WordRAM.Deduplicate.Count public import Complexity.Algorithms.WordRAM.Deduplicate.Model public import Complexity.Algorithms.WordRAM.Deduplicate.Program -public import Complexity.Algorithms.WordRAM.Dijkstra.BinaryHeap -public import Complexity.Algorithms.WordRAM.Dijkstra.ExactNotation -public import Complexity.Algorithms.WordRAM.Dijkstra.HeapProgram -public import Complexity.Algorithms.WordRAM.Dijkstra.HeapRelax -public import Complexity.Algorithms.WordRAM.Dijkstra.HeapSift public import Complexity.Algorithms.WordRAM.Dijkstra.Fibonacci public import Complexity.Algorithms.WordRAM.Dijkstra.FibonacciStatement public import Complexity.Algorithms.WordRAM.Dijkstra.GenericPQ