feat(ExactTriangle): implement native cyclic-Strassen modular counting

Complexity.lean

@@ -1,5 +1,7 @@
 module  -- shake: keep-all --deprecated_module: ignore
 
+public import Complexity.Algorithms.ExactTriangle.CountPrime
+public import Complexity.Algorithms.ExactTriangle.CountPrime.Extraction
 public import Complexity.Algorithms.ExactTriangle.Defs
 public import Complexity.Algorithms.ExactTriangle.FailureAccounting
 public import Complexity.Algorithms.ExactTriangle.Hashing
@@ -92,6 +94,10 @@ public import Complexity.Algorithms.MatrixMultiplication.Schonhage.Sparse
 public import Complexity.Algorithms.MatrixMultiplication.Schonhage.SplitDigits
 public import Complexity.Algorithms.MatrixMultiplication.Schonhage.ThinProduct
 public import Complexity.Algorithms.MatrixMultiplication.Schonhage.Tiling
+public import Complexity.Algorithms.MatrixMultiplication.Strassen.Coefficients
+public import Complexity.Algorithms.MatrixMultiplication.Strassen.Correctness
+public import Complexity.Algorithms.MatrixMultiplication.Strassen.Defs
+public import Complexity.Algorithms.MatrixMultiplication.Strassen.Phases
 public import Complexity.Algorithms.Reachability.BreadthFirst
 public import Complexity.Algorithms.Reachability.InductiveCounting
 public import Complexity.Algorithms.Reachability.InductiveCounting.Bounds
@@ -151,6 +157,34 @@ public import Complexity.Algorithms.WordRAM.ExactTriangle.Catalogue
 public import Complexity.Algorithms.WordRAM.ExactTriangle.ChunkCatalogue
 public import Complexity.Algorithms.WordRAM.ExactTriangle.Chunks
 public import Complexity.Algorithms.WordRAM.ExactTriangle.Classes
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Access
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Bounds
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.CoefficientMatrices
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Exec
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.ExecReady
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Flat
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Input
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Layout
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Logarithm
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Matrix
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.MatrixCell
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.MatrixReadback
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.MatrixResidue
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Monomial
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Product
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Program
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.ResidueArray
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.ResidueArrays
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Scan
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.ScanCell
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.ScanReadback
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Setup
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Solve
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.SourceBounds
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Spread
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.TimeBound
+public import Complexity.Algorithms.WordRAM.ExactTriangle.CountPrime.Width
 public import Complexity.Algorithms.WordRAM.ExactTriangle.FixedPrime
 public import Complexity.Algorithms.WordRAM.ExactTriangle.FixedPrime.Accounting
 public import Complexity.Algorithms.WordRAM.ExactTriangle.FixedPrime.Controller
@@ -399,6 +433,35 @@ public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Schonhage.Wante
 public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Schonhage.Wanted.ScatterWrite
 public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Schonhage.Wanted.Setup
 public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Schonhage.Wanted.Width
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Accumulate
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.ArithmeticForm
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Bounds
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Child
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Control
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.CopyForm
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Dispatch
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Entry
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Exec
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Finish
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Form
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Frame
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Initialize
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Leaf
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.LeafStep
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Next
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Pending
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Prepare
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Program
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Progress
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Push
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Quadrant
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.QuarterUpdate
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Receive
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Restore
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Resume
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.State
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Traversal
+public import Complexity.Algorithms.WordRAM.MatrixMultiplication.Strassen.Zero
 public import Complexity.Algorithms.WordRAM.NatClog
 public import Complexity.Algorithms.WordRAM.NatGcd
 public import Complexity.Algorithms.WordRAM.NatLog
@@ -2431,6 +2494,9 @@ public import Complexity.Computability.WordRAM.Randomized.UniformLoopExit
 public import Complexity.Computability.WordRAM.RegisterFrame
 public import Complexity.Computability.WordRAM.Sequence
 public import Complexity.Computability.WordRAM.Signed.Arithmetic
+public import Complexity.Computability.WordRAM.Signed.ArrayArithmetic.Body
+public import Complexity.Computability.WordRAM.Signed.ArrayArithmetic.Exec
+public import Complexity.Computability.WordRAM.Signed.ArrayArithmetic.Program
 public import Complexity.Computability.WordRAM.Signed.Comparison
 public import Complexity.Computability.WordRAM.Signed.Compile
 public import Complexity.Computability.WordRAM.Signed.Encoding
@@ -3025,6 +3091,7 @@ public import Complexity.Foundations.Data.Nat.Factorial
 public import Complexity.Foundations.Data.Nat.Fields
 public import Complexity.Foundations.Data.Nat.IntermediateValue
 public import Complexity.Foundations.Data.Nat.Log
+public import Complexity.Foundations.Data.Nat.Morton
 public import Complexity.Foundations.Data.Nat.OddMulMod
 public import Complexity.Foundations.Data.Nat.OfDigits
 public import Complexity.Foundations.Data.Nat.PolynomialBound.Closure

Complexity/Algorithms/ExactTriangle/CountPrime.lean

@@ -0,0 +1,147 @@
+/-
+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.Hashing
+public import Complexity.Algorithms.MatrixMultiplication.Strassen.Correctness
+
+/-!
+# Modular triangle counts from the actual seven-product coefficient recursion
+
+The coefficient arrays contain zero padding and one monomial for each real edge residue.
+`CountPrime.count` applies the executable finite Strassen recursion and selects one product
+coefficient for each original AB pair. It equals the existing `countZeroMod`, including
+modular false positives. This is a mathematical specification and correctness bridge for
+the native counter, not a native execution or running-time theorem.
+
+The construction selectively adapts Spec/Sec3/Theorem17/Strassen.lean, the `matPList`,
+`matQList`, `countBy` and `countOf_eq` sections, at revision
+16a5ce41304fcdaea4e0be7c132e1c361dc26da5. Existing `matP`, `matQ`, cyclic coefficients and
+modular triangle counts are reused; no parallel triangle problem or polynomial ring is added.
+-/
+
+public section
+
+namespace Cslib.Algorithms.ExactTriangle.CountPrime
+
+open Cslib.CyclicIntRing Cslib.Algorithms.MatrixMultiplication.Strassen
+
+/-- Native P coefficient specification, including true zero padding outside the original side. -/
+@[expose] def pCoefficients {n : ℕ} (p : ℕ) (wAC : Fin n → Fin n → ℤ) (i : ℕ) : ℤ :=
+  if h : Nat.Morton.row (i / p) < n ∧ Nat.Morton.column (i / p) < n then
+    if i % p = (wAC ⟨_, h.1⟩ ⟨_, h.2⟩ % (p : ℤ)).toNat then 1 else 0
+  else 0
+
+/-- Native Q coefficient specification; the original BC coordinates are transposed. -/
+@[expose] def qCoefficients {n : ℕ} (p : ℕ) (wBC : Fin n → Fin n → ℤ) (i : ℕ) : ℤ :=
+  if h : Nat.Morton.row (i / p) < n ∧ Nat.Morton.column (i / p) < n then
+    if i % p = (wBC ⟨_, h.2⟩ ⟨_, h.1⟩ % (p : ℤ)).toNat then 1 else 0
+  else 0
+
+/-- Every actual monomial coefficient, including padding, has magnitude at most one. -/
+theorem pCoefficients_natAbs_le_one {n : ℕ} (p : ℕ) (wAC : Fin n → Fin n → ℤ) (i : ℕ) :
+    (pCoefficients p wAC i).natAbs ≤ 1 := by
+  unfold pCoefficients
+  split_ifs <;> decide
+
+/-- The transposed monomial matrix has the same unit coefficient bound. -/
+theorem qCoefficients_natAbs_le_one {n : ℕ} (p : ℕ) (wBC : Fin n → Fin n → ℤ) (i : ℕ) :
+    (qCoefficients p wBC i).natAbs ≤ 1 := by
+  unfold qCoefficients
+  split_ifs <;> decide
+
+private noncomputable def pad {n p : ℕ} (M : Matrix (Fin n) (Fin n) (CyclicIntRing p))
+    (a b : ℕ) : CyclicIntRing p := if h : a < n ∧ b < n then M ⟨a, h.1⟩ ⟨b, h.2⟩ else 0
+
+private theorem pCoefficients_eq {n p : ℕ} (hp : p ≠ 0) (wAC : Fin n → Fin n → ℤ) :
+    pCoefficients p wAC = coefficients hp (pad (matP wAC p)) := by
+  funext i
+  by_cases h : Nat.Morton.row (i / p) < n ∧ Nat.Morton.column (i / p) < n
+  · simp only [pCoefficients, coefficients, pad, dite_eq_left h, matP, Matrix.of_apply,
+      coeff_x_pow]
+    rw [Nat.mod_eq_of_lt (show (wAC ⟨_, h.1⟩ ⟨_, h.2⟩ % (p : ℤ)).toNat < p from
+      Int.natMod_lt hp)]
+  · simp only [pCoefficients, coefficients, pad, dite_eq_right h, map_zero]
+
+private theorem qCoefficients_eq {n p : ℕ} (hp : p ≠ 0) (wBC : Fin n → Fin n → ℤ) :
+    qCoefficients p wBC = coefficients hp (pad (matQ wBC p)) := by
+  funext i
+  by_cases h : Nat.Morton.row (i / p) < n ∧ Nat.Morton.column (i / p) < n
+  · simp only [qCoefficients, coefficients, pad, dite_eq_left h, matQ, Matrix.of_apply,
+      coeff_x_pow]
+    rw [Nat.mod_eq_of_lt (show (wBC ⟨_, h.2⟩ ⟨_, h.1⟩ % (p : ℤ)).toNat < p from
+      Int.natMod_lt hp)]
+  · simp only [qCoefficients, coefficients, pad, dite_eq_right h, map_zero]
+
+private theorem sum_pad {n p N : ℕ} (hN : n ≤ N)
+    (M K : Matrix (Fin n) (Fin n) (CyclicIntRing p)) (a b : Fin n) :
+    ∑ c ∈ Finset.range N, pad M a c * pad K c b = (M * K) a b := by
+  obtain ⟨k, rfl⟩ := Nat.exists_eq_add_of_le hN
+  have ht : ∑ c ∈ Finset.range k, pad M a (n + c) * pad K (n + c) b = 0 :=
+    Finset.sum_eq_zero fun c _ => by simp [pad]
+  rw [Finset.sum_range_add, ht, add_zero, Matrix.mul_apply, Finset.sum_range]
+  exact Finset.sum_congr rfl fun c _ => by simp [pad]
+
+/-- Every original pair and residue in the executable finite Strassen coefficient array
+has exactly the existing cyclic matrix-product meaning, including power-of-two padding. -/
+theorem multiply_entry_eq_coeff {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wBC wAC : Fin n → Fin n → ℤ) (a b : Fin n) (r : ℕ) (hr : r < p) :
+    multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+      (Nat.Morton.index a b * p + r) =
+      CyclicIntRing.coeff hp r ((matP wAC p * matQ wBC p) a b) := by
+  have hi := Nat.Morton.index_lt (a.isLt.trans_le hn) (b.isLt.trans_le hn)
+  have hmul := Nat.mul_le_mul_right p (Nat.succ_le_of_lt hi)
+  rw [Nat.succ_mul] at hmul
+  have hidx : Nat.Morton.index a b * p + r < 4 ^ k * p := by omega
+  rw [multiply_eq_coefficients hp k (pad (matP wAC p)) (pad (matQ wBC p))
+    (pCoefficients p wAC) (qCoefficients p wBC)
+    (fun i _ => congrFun (pCoefficients_eq hp wAC) i)
+    (fun i _ => congrFun (qCoefficients_eq hp wBC) i) hidx,
+    coefficients_at hp _ a b hr, sum_pad hn]
+
+/-- The selected real product coefficient is a nonnegative count of at most n middle
+vertices. Native signed-to-unsigned extraction uses this bound before every sum. -/
+theorem multiply_entry_bounds {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wBC wAC : Fin n → Fin n → ℤ) (a b : Fin n) (r : ℕ) (hr : r < p) :
+    0 ≤ multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+        (Nat.Morton.index a b * p + r) ∧
+      multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+        (Nat.Morton.index a b * p + r) ≤ (n : ℤ) := by
+  classical
+  rw [multiply_entry_eq_coeff hp hn wBC wAC a b r hr,
+    coeff_matP_mul_matQ wBC wAC hp a b r hr]
+  have hc := Finset.card_le_card (Finset.filter_subset (s := Finset.univ)
+    (p := fun c : Fin n => wAC a c + wBC b c ≡ (r : ℤ) [ZMOD (p : ℤ)]))
+  simp only [Finset.card_univ, Fintype.card_fin] at hc
+  constructor
+  · positivity
+  · exact_mod_cast hc
+
+/-- Select the (-AB)-residue coefficient for each original pair after the real seven-product
+finite recursion. The depth controls zero padding; it is not a supplied multiplication. -/
+@[expose] def count {n : ℕ} (p k : ℕ) (wAB wBC wAC : Fin n → Fin n → ℤ) : ℤ :=
+  ∑ a : Fin n, ∑ b : Fin n,
+    multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+      (Nat.Morton.index a b * p + ((-wAB a b) % (p : ℤ)).toNat)
+
+/-- The finite coefficient algorithm counts exactly the existing modular-zero predicate,
+for every nonzero modulus and every power-of-two padding covering the input. -/
+theorem count_eq {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) :
+    count p k wAB wBC wAC = (countZeroMod wAB wBC wAC p : ℤ) := by
+  have hc := falsePositiveCount_add_zeroCount_eq_sum_coeff wAB wBC wAC hp
+  rw [← countZeroMod_eq] at hc
+  rw [count, hc]
+  apply Finset.sum_congr rfl
+  intro a _
+  apply Finset.sum_congr rfl
+  intro b _
+  have hr : ((-wAB a b) % (p : ℤ)).toNat < p := Int.natMod_lt hp
+  exact multiply_entry_eq_coeff hp hn wBC wAC a b _ hr
+
+end Cslib.Algorithms.ExactTriangle.CountPrime

Complexity/Algorithms/ExactTriangle/CountPrime/Extraction.lean

@@ -0,0 +1,91 @@
+/-
+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.ExactTriangle.CountPrime
+
+/-!
+# Nonnegative selected coefficients for modular triangle counting
+
+These bridges justify the real counter's unsigned extraction and summation from the existing
+cyclic product specification. They do not assume a matrix multiplication or native execution.
+-/
+
+public section
+
+namespace Cslib.Algorithms.ExactTriangle.CountPrime
+
+open Cslib.Algorithms.MatrixMultiplication.Strassen
+
+/-- The unsigned value extracted from one actual original-coordinate product coefficient. -/
+@[expose] def pairCoefficient {n : ℕ} (p k : ℕ)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) (a b : Fin n) : ℕ :=
+  (multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+    (Nat.Morton.index a b * p + ((-wAB a b) % (p : ℤ)).toNat)).natAbs
+
+/-- Nonnegativity makes signed decoding by division by two an exact natural count. -/
+theorem pairCoefficient_cast {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) (a b : Fin n) :
+    (pairCoefficient p k wAB wBC wAC a b : ℤ) =
+      multiply p k (pCoefficients p wAC) (qCoefficients p wBC)
+        (Nat.Morton.index a b * p + ((-wAB a b) % (p : ℤ)).toNat) := by
+  exact Int.natAbs_of_nonneg
+    (multiply_entry_bounds hp hn wBC wAC a b _ (Int.natMod_lt hp)).1
+
+/-- Each selected coefficient counts at most the original side many middle vertices. -/
+theorem pairCoefficient_le {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) (a b : Fin n) :
+    pairCoefficient p k wAB wBC wAC a b ≤ n := by
+  have he := pairCoefficient_cast hp hn wAB wBC wAC a b
+  have hb := (multiply_entry_bounds hp hn wBC wAC a b
+    (((-wAB a b) % (p : ℤ)).toNat) (Int.natMod_lt hp)).2
+  omega
+
+/-- All selected original-coordinate coefficients sum to the existing modular triangle count. -/
+theorem sum_pairCoefficient {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) :
+    ∑ a : Fin n, ∑ b : Fin n, pairCoefficient p k wAB wBC wAC a b =
+      countZeroMod wAB wBC wAC p := by
+  apply Int.natCast_inj.mp
+  push_cast
+  simp_rw [pairCoefficient_cast hp hn]
+  exact count_eq hp hn wAB wBC wAC
+
+/-- The selected counts in the exact original row-major order traversed by native code.
+Out-of-range entries are zero, so empty input requires no arbitrary coordinate witness. -/
+@[expose] def rowMajorCoefficient {n : ℕ} (p k : ℕ)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) (i : ℕ) : ℕ :=
+  if h : i < n * n then
+    let q := finProdFinEquiv.symm (⟨i, h⟩ : Fin (n * n))
+    pairCoefficient p k wAB wBC wAC q.1 q.2
+  else 0
+
+/-- The real row-major coefficient bound supplies the native prefix-accumulator admission. -/
+theorem rowMajorCoefficient_le {n p k i : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) :
+    rowMajorCoefficient p k wAB wBC wAC i ≤ n := by
+  unfold rowMajorCoefficient
+  split_ifs
+  · exact pairCoefficient_le hp hn wAB wBC wAC _ _
+  · omega
+
+/-- Summation in the actual native scan order preserves the existing modular count. -/
+theorem sum_rowMajorCoefficient {n p k : ℕ} (hp : p ≠ 0) (hn : n ≤ 2 ^ k)
+    (wAB wBC wAC : Fin n → Fin n → ℤ) :
+    ∑ i ∈ Finset.range (n * n), rowMajorCoefficient p k wAB wBC wAC i =
+      countZeroMod wAB wBC wAC p := by
+  rw [Finset.sum_range]
+  simp only [rowMajorCoefficient, Fin.isLt, dite_true]
+  rw [← Equiv.sum_comp (finProdFinEquiv : Fin n × Fin n ≃ Fin (n * n))]
+  rw [Fintype.sum_prod_type]
+  change (∑ a : Fin n, ∑ b : Fin n, pairCoefficient p k wAB wBC wAC
+    (finProdFinEquiv.symm (finProdFinEquiv (a, b))).1
+    (finProdFinEquiv.symm (finProdFinEquiv (a, b))).2) = _
+  simpa only [Equiv.symm_apply_apply] using sum_pairCoefficient hp hn wAB wBC wAC
+
+end Cslib.Algorithms.ExactTriangle.CountPrime