mirror of
https://github.com/allaunthefox/SilverSight.git
synced 2026-07-31 01:25:21 +00:00
allaun
359effdfc6
STEP 1 — integer spiral-index packing (EXACT): pack_spiral_index/ unpack_spiral_index with corpus-wide base=2·max|c|+1 and offset-encoded digits, recorded in the codebook JSON header; round-trip verified over all 250 matrices. Fixes the sign bug in integration_sprint.py's integer path (unoffset signed packing collapsed (-1,1) with (1,0) in base 2). STEP 2 — f(n) layout (DECORATIVE): per-entry radius_sq + angle_frac, computed in decimal arithmetic (indices reach ~1e44 >> 2^52); verified against VERIFICATION_LOG V007 f(20121); min angular separation ≈4.64e-6 of a turn reported. STEP 3 — Cartan Δ-floor (--cartan-floor): Δ=17/1792 exact (CartanConnection.lean:70) as Fisher-distance resolution floor on Δ₇; 13 sub-Δ pairs merge 196 fingerprints into 187 Δ-resolution codewords vs 9 gap-rule clusters; emitted alongside, not replacing, the 3×-median rule. Caveat asserted in tests: the |c|/Σ|c| simplex projection is lossy (one corpus pair at Fisher distance 0) — similarity, not identity. STEP 4 — torus-winding saturation fix: TorusWinding gains a_exact/b_exact (plain ℤ, no Q16.16 clamp); torus_to_spiral_index prefers them, making the round trip lossless (regression-tested at n=1e9). BraidEigensolid.lean needs the same widening — noted, Lean untouched. Also corrects charpoly_codebook.py's stale docstring counts (192/180 → 196/182, cross-checked against the exact Faddeev–LeVerrier fingerprints: identical on all 250 matrices). Tests: 43 passed (19 new). Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
893 lines
39 KiB
Python
893 lines
39 KiB
Python
#!/usr/bin/env python3
|
||
# /// script
|
||
# requires-python = ">=3.10"
|
||
# dependencies = []
|
||
# ///
|
||
"""
|
||
spectral_codebook.py — Gap-aware spectral codebook over the PIST matrix corpus.
|
||
|
||
Implements the "Next Steps" of docs/SPECTRAL_CODEBOOK_ANALYSIS.md, with
|
||
corrections found while building it (see docs/SPECTRAL_CODEBOOK_GENERATOR.md):
|
||
|
||
1. Parse the 250 8×8 integer braid adjacency matrices from
|
||
formal/SilverSight/PIST/Matrices250.lean (regex on the auto-generated
|
||
`def <name> : Array (Array Int)` blocks, plus the `findMatrix` cases).
|
||
2. Compute a full spectral profile per matrix:
|
||
- EXACT integer characteristic-polynomial coefficients via the
|
||
Faddeev–LeVerrier algorithm in fractions.Fraction arithmetic
|
||
(the primary, float-free fingerprint),
|
||
- spectral radius λ as the max root modulus of that exact char poly
|
||
(NumPy eigvals fast path; pure-stdlib Durand–Kerner fallback),
|
||
- power-iteration λ retained only for traceability: the analysis
|
||
doc's 9 "mid band" values in (0.5, 1) are power-iteration
|
||
NON-CONVERGENCE artifacts — all 9 matrices have exact ρ = 1.0
|
||
(peripheral spectrum with several eigenvalues of modulus 1).
|
||
3. Gap-aware quantization with a sample-size guard: dedupe byte-identical
|
||
matrices first (the corpus has 238 distinct matrices over 250 ids),
|
||
find gaps > 3× median gap between consecutive distinct λ, then merge
|
||
any segment with fewer than MIN_SUPPORT distinct matrices into its
|
||
neighbor across the smaller gap. Points isolated above a suppressed
|
||
gap with sub-threshold occupancy are flagged as sparse-tail outliers
|
||
rather than treated as clusters.
|
||
4. Round-trip: Codebook.encode(matrix) → (codeword, index) and
|
||
Codebook.decode(codeword, index) → equation_id are mutually inverse
|
||
over the corpus. Neither λ nor the char poly is injective (nilpotent
|
||
x^8 class alone has 31 distinct members), so the decode key is the
|
||
within-cluster rank index; all collision classes are explicit in the
|
||
emitted JSON.
|
||
|
||
Pure stdlib. NumPy is used only as an optional fast path when importable;
|
||
the pure-Python path is authoritative.
|
||
|
||
Usage:
|
||
python3 python/spectral_codebook.py # build + report
|
||
python3 python/spectral_codebook.py --check-manifold # incl. 278-row corpus
|
||
python3 python/spectral_codebook.py --out data/spectral_codebook.json
|
||
"""
|
||
from __future__ import annotations
|
||
|
||
import argparse
|
||
import json
|
||
import math
|
||
import re
|
||
import sys
|
||
from datetime import date
|
||
from fractions import Fraction
|
||
from pathlib import Path
|
||
from typing import Dict, List, Optional, Sequence, Tuple
|
||
|
||
try: # optional fast path only; both paths agree to float rounding
|
||
import numpy as _np
|
||
except ImportError: # pragma: no cover
|
||
_np = None
|
||
|
||
REPO_ROOT = Path(__file__).resolve().parent.parent
|
||
MATRICES_LEAN = REPO_ROOT / "formal" / "SilverSight" / "PIST" / "Matrices250.lean"
|
||
MANIFOLD_LEAN = REPO_ROOT / "formal" / "SilverSight" / "RRC" / "Q16_16Manifold.lean"
|
||
DEFAULT_OUT = REPO_ROOT / "data" / "spectral_codebook.json"
|
||
|
||
POWER_ITERS = 300 # matches docs/SPECTRAL_CODEBOOK_ANALYSIS.md
|
||
GAP_FACTOR = 3.0 # boundary threshold = GAP_FACTOR × median gap
|
||
MIN_SUPPORT = 10 # min distinct matrices per cluster (sample-size guard)
|
||
LAMBDA_DP = 4 # rounding used for "unique λ" counting (as in doc)
|
||
|
||
Matrix = Tuple[Tuple[int, ...], ...]
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# Parsing Matrices250.lean
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
|
||
_DEF_RE = re.compile(
|
||
r"^def\s+(\w+)\s*:\s*Array\s*\(Array\s+Int\)\s*:=\s*\n\s*#\[(.*?)\n\s*\]",
|
||
re.MULTILINE | re.DOTALL,
|
||
)
|
||
_ROW_RE = re.compile(r"#\[([^\]]*)\]")
|
||
_FIND_CASE_RE = re.compile(r'\|\s*"([^"]+)"\s*=>\s*some\s+(\w+)')
|
||
|
||
|
||
def parse_matrices_lean(path: Path = MATRICES_LEAN) -> Dict[str, Matrix]:
|
||
"""Parse `def <name> : Array (Array Int) := #[#[...], ...]` blocks.
|
||
|
||
Returns an equation_id → matrix mapping (row-major tuples of ints),
|
||
resolving ids through the `findMatrix` match cases when present so the
|
||
mapping survives Lean-identifier sanitization of ids.
|
||
"""
|
||
text = path.read_text()
|
||
defs: Dict[str, Matrix] = {}
|
||
for name, body in _DEF_RE.findall(text):
|
||
rows = tuple(
|
||
tuple(int(x) for x in row.split(",") if x.strip())
|
||
for row in _ROW_RE.findall(body)
|
||
)
|
||
if rows:
|
||
defs[name] = rows
|
||
|
||
cases = _FIND_CASE_RE.findall(text)
|
||
if cases:
|
||
mats = {eid: defs[name] for eid, name in cases if name in defs}
|
||
else: # fall back: def name == equation id
|
||
mats = dict(defs)
|
||
|
||
for eid, mat in mats.items():
|
||
n = len(mat)
|
||
if any(len(r) != n for r in mat):
|
||
raise ValueError(f"non-square matrix for {eid}")
|
||
return mats
|
||
|
||
|
||
def parse_manifold_ids(path: Path = MANIFOLD_LEAN) -> List[str]:
|
||
"""All equationId occurrences (with multiplicity) in the fixture corpus."""
|
||
text = path.read_text()
|
||
return re.findall(r'equationId\s*:=\s*"([^"]+)"', text)
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# Spectral profile
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
|
||
def charpoly_coeffs(mat: Sequence[Sequence[int]]) -> Tuple[int, ...]:
|
||
"""Exact characteristic-polynomial coefficients via Faddeev–LeVerrier.
|
||
|
||
For an n×n integer matrix M returns (c1, ..., cn) with
|
||
|
||
det(λI − M) = λ^n + c1·λ^(n−1) + ... + cn.
|
||
|
||
Runs in exact fractions.Fraction arithmetic; every ck is provably an
|
||
integer (char poly of an integer matrix), which is asserted.
|
||
"""
|
||
n = len(mat)
|
||
m = [[Fraction(x) for x in row] for row in mat]
|
||
aux = [[Fraction(0)] * n for _ in range(n)] # M_0 = 0
|
||
coeffs: List[int] = []
|
||
c = Fraction(1)
|
||
for k in range(1, n + 1):
|
||
# M_k = M · (M_{k−1} + c_{k−1}·I)
|
||
shifted = [
|
||
[aux[i][j] + (c if i == j else 0) for j in range(n)]
|
||
for i in range(n)
|
||
]
|
||
aux = [
|
||
[sum(m[i][t] * shifted[t][j] for t in range(n)) for j in range(n)]
|
||
for i in range(n)
|
||
]
|
||
c = -sum(aux[i][i] for i in range(n)) / k
|
||
assert c.denominator == 1, "char-poly coefficient must be integral"
|
||
coeffs.append(int(c))
|
||
return tuple(coeffs)
|
||
|
||
|
||
def _durand_kerner_max_modulus(coeffs: Sequence[int]) -> float:
|
||
"""Max root modulus of the monic polynomial x^n + c1 x^(n-1) + ... + cn.
|
||
|
||
Pure-stdlib root finder (Durand–Kerner). Zero roots are stripped
|
||
exactly first, so nilpotent classes return 0.0 with no iteration.
|
||
"""
|
||
poly = [1] + [int(c) for c in coeffs]
|
||
# strip exact zero roots (trailing zero coefficients)
|
||
while len(poly) > 1 and poly[-1] == 0:
|
||
poly.pop()
|
||
n = len(poly) - 1
|
||
if n == 0:
|
||
return 0.0
|
||
# Cauchy bound and standard (0.4 + 0.9i)^k seeds
|
||
bound = 1.0 + max(abs(c) for c in poly[1:])
|
||
roots = [complex(0.4, 0.9) ** k * (bound / 1.1) for k in range(n)]
|
||
|
||
def ev(z: complex) -> complex:
|
||
r = complex(0)
|
||
for c in poly:
|
||
r = r * z + c
|
||
return r
|
||
|
||
for _ in range(1000):
|
||
delta = 0.0
|
||
new_roots = []
|
||
for i, z in enumerate(roots):
|
||
den = complex(1)
|
||
for j, w in enumerate(roots):
|
||
if j != i:
|
||
den *= (z - w)
|
||
dz = ev(z) / den if den != 0 else complex(0)
|
||
new_roots.append(z - dz)
|
||
delta = max(delta, abs(dz))
|
||
roots = new_roots
|
||
if delta < 1e-13:
|
||
break
|
||
return max(abs(z) for z in roots)
|
||
|
||
|
||
def spectral_radius(mat: Sequence[Sequence[int]],
|
||
coeffs: Optional[Sequence[int]] = None) -> float:
|
||
"""Exact-spectrum spectral radius ρ(M) = max |eigenvalue|.
|
||
|
||
NumPy eigvals fast path when available; otherwise Durand–Kerner on the
|
||
exact integer char poly. Do NOT use power iteration here: it fails to
|
||
converge on peripheral spectra (several eigenvalues of equal modulus),
|
||
which is exactly what the ρ = 1 matrices in this corpus have.
|
||
"""
|
||
if _np is not None:
|
||
m = _np.asarray(mat, dtype=float)
|
||
if m.size == 0:
|
||
return 0.0
|
||
return float(max(abs(_np.linalg.eigvals(m))))
|
||
if coeffs is None:
|
||
coeffs = charpoly_coeffs(mat)
|
||
return _durand_kerner_max_modulus(coeffs)
|
||
|
||
|
||
def power_iteration(mat: Sequence[Sequence[int]], iters: int = POWER_ITERS) -> float:
|
||
"""Power-iteration λ estimate — UNTRUSTED, kept only for traceability
|
||
against data/spectral_codebook_raw.json. Does not converge for
|
||
peripheral spectra (see spectral_radius)."""
|
||
n = len(mat)
|
||
if n == 0:
|
||
return 0.0
|
||
v = [1.0 / math.sqrt(n)] * n
|
||
lam = 0.0
|
||
for _ in range(iters):
|
||
w = [sum(mat[i][j] * v[j] for j in range(n)) for i in range(n)]
|
||
norm = math.sqrt(sum(x * x for x in w))
|
||
if norm < 1e-12:
|
||
return 0.0
|
||
lam = norm
|
||
v = [x / norm for x in w]
|
||
return lam
|
||
|
||
|
||
def density(mat: Sequence[Sequence[int]]) -> float:
|
||
"""Total entry weight / n² (matches SpectralN.lean and the raw JSON)."""
|
||
n = len(mat)
|
||
return sum(sum(row) for row in mat) / (n * n) if n else 0.0
|
||
|
||
|
||
def frobenius(mat: Sequence[Sequence[int]]) -> float:
|
||
return math.sqrt(sum(x * x for row in mat for x in row))
|
||
|
||
|
||
def spectral_profile(mat: Sequence[Sequence[int]]) -> dict:
|
||
coeffs = charpoly_coeffs(mat)
|
||
return {
|
||
"spectral_radius": round(spectral_radius(mat, coeffs), 6),
|
||
"lambda_power_iteration": round(power_iteration(mat), 6),
|
||
"charpoly": list(coeffs),
|
||
"density": round(density(mat), 4),
|
||
"frobenius": round(frobenius(mat), 4),
|
||
"trace": sum(mat[i][i] for i in range(len(mat))),
|
||
}
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# Φ-corkscrew: integer spiral-index packing (exact) + f(n) layout
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
#
|
||
# The packing layer is EXACT and provable: with a corpus-wide base
|
||
# B = 2·max|c| + 1 and offset o = max|c|,
|
||
# every offset digit c_i + o lies in [0, B), so
|
||
# n = Σ (c_i + o)·B^i
|
||
# is the base-B positional numeral of the digit tuple — injective on
|
||
# char-poly fingerprints by uniqueness of positional representation, and
|
||
# unpack_spiral_index inverts it exactly. It inherits charpoly's
|
||
# NON-injectivity on matrices (cospectral classes share one index).
|
||
#
|
||
# The layout layer f(n) = (√n·cos(nψ), √n·sin(nψ)), ψ = 2π/φ², is
|
||
# DECORATIVE placement only: injective because 1/φ² is irrational
|
||
# (distinct n give distinct angles), information-neutral (everything is
|
||
# recoverable from n alone), low-discrepancy by the three-distance
|
||
# theorem. It carries no information beyond spiral_index.
|
||
|
||
def packing_params(fingerprints: Sequence[Sequence[int]]) -> Tuple[int, int]:
|
||
"""Corpus-wide (base, offset) for spiral-index packing.
|
||
|
||
offset = max|c| over every coefficient of every fingerprint;
|
||
base = 2·offset + 1, so offset-encoded digits fill [0, base).
|
||
"""
|
||
offset = max((abs(c) for fp in fingerprints for c in fp), default=0)
|
||
return 2 * offset + 1, offset
|
||
|
||
|
||
def pack_spiral_index(coeffs: Sequence[int], base: int, offset: int) -> int:
|
||
"""Little-endian base-B packing of offset-encoded coefficients.
|
||
|
||
Exact integer arithmetic; injective for digits in range (raises
|
||
otherwise). Replaces the float phinary packing of
|
||
integration_sprint.phi_corkscrew_index, which was not injective.
|
||
"""
|
||
n = 0
|
||
for i, c in enumerate(coeffs):
|
||
d = c + offset
|
||
if not 0 <= d < base:
|
||
raise ValueError(
|
||
f"coefficient {c} out of packing range for "
|
||
f"base={base}, offset={offset}")
|
||
n += d * base ** i
|
||
return n
|
||
|
||
|
||
def unpack_spiral_index(n: int, base: int, offset: int,
|
||
length: int = 8) -> Tuple[int, ...]:
|
||
"""Exact inverse of pack_spiral_index."""
|
||
if n < 0:
|
||
raise ValueError("spiral index must be non-negative")
|
||
out: List[int] = []
|
||
for _ in range(length):
|
||
n, d = divmod(n, base)
|
||
out.append(d - offset)
|
||
if n:
|
||
raise ValueError("spiral index exceeds length digits in this base")
|
||
return tuple(out)
|
||
|
||
|
||
# 1/φ² = (3 − √5)/2: the fractional rotation per step of the corkscrew.
|
||
def angle_frac(n: int) -> float:
|
||
"""frac(n/φ²) — the corkscrew angle of index n as a fraction of 2π.
|
||
|
||
Computed in decimal arithmetic at a precision that covers the integer
|
||
part of n·(3−√5)/2, because float64 loses the fractional part
|
||
entirely once n exceeds 2^52 (spiral indices here reach base^8).
|
||
"""
|
||
from decimal import Decimal, getcontext
|
||
getcontext().prec = len(str(abs(n))) + 25
|
||
alpha = (Decimal(3) - Decimal(5).sqrt()) / 2
|
||
return float((Decimal(n) * alpha) % 1)
|
||
|
||
|
||
def corkscrew_xy(n: int) -> Tuple[float, float]:
|
||
"""f(n) = (√n·cos(nψ), √n·sin(nψ)), ψ = 2π/φ². Placement only."""
|
||
theta = 2.0 * math.pi * angle_frac(n)
|
||
r = math.sqrt(n)
|
||
return r * math.cos(theta), r * math.sin(theta)
|
||
|
||
|
||
def min_angular_separation(indices: Sequence[int]) -> Optional[float]:
|
||
"""Minimum circular distance (as a fraction of a full turn) between
|
||
the corkscrew angles of distinct spiral indices."""
|
||
uniq = sorted({angle_frac(n) for n in set(indices)})
|
||
if len(uniq) < 2:
|
||
return None
|
||
gaps = [uniq[i + 1] - uniq[i] for i in range(len(uniq) - 1)]
|
||
gaps.append(1.0 - uniq[-1] + uniq[0]) # wrap-around
|
||
return min(gaps)
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# Cartan ∆-floor quantization (Fisher metric on Δ₇)
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
#
|
||
# ∆ = 17/1792 is the minimum eigenvalue of the Cartan crossing blocks,
|
||
# proven exact in formal/.../CartanConnection.lean:70 (see
|
||
# docs/cartan_fingerprint.md). It is adopted here as the resolution
|
||
# floor on the Fisher metric: two codewords closer than ∆ are treated
|
||
# as indistinguishable at operator resolution. This gives a PRINCIPLED
|
||
# alternative to the 3×-median-gap heuristic; both are emitted so the
|
||
# rules can be compared, neither replaces the other.
|
||
|
||
CARTAN_DELTA = Fraction(17, 1792) # exact; CartanConnection.lean:70
|
||
|
||
|
||
def simplex_project(coeffs: Sequence[int]) -> Tuple[float, ...]:
|
||
"""Project a fingerprint onto Δ₇: p_i = |c_i| / Σ|c_j|.
|
||
|
||
Convention: the all-zero fingerprint (nilpotent x⁸ class) maps to
|
||
the uniform distribution — flagged in the report, since that choice
|
||
is a convention, not a theorem.
|
||
"""
|
||
total = sum(abs(c) for c in coeffs)
|
||
k = len(coeffs)
|
||
if total == 0:
|
||
return tuple(1.0 / k for _ in range(k))
|
||
return tuple(abs(c) / total for c in coeffs)
|
||
|
||
|
||
def fisher_distance(p: Sequence[float], q: Sequence[float]) -> float:
|
||
"""d_F(p, q) = 2·arccos(Σ√(p_i·q_i)) — Fisher geodesic distance on
|
||
the probability simplex (Bhattacharyya angle × 2)."""
|
||
bc = sum(math.sqrt(a * b) for a, b in zip(p, q))
|
||
return 2.0 * math.acos(min(1.0, bc))
|
||
|
||
|
||
def cartan_floor_report(profiles: Dict[str, dict]) -> dict:
|
||
"""All-pairs Fisher distances between distinct fingerprints; pairs
|
||
below ∆ = 17/1792 are merged into ∆-resolution components (union-
|
||
find), giving the codeword count the Cartan floor supports."""
|
||
fp_rep: Dict[Tuple[int, ...], str] = {}
|
||
for eid in sorted(profiles):
|
||
fp_rep.setdefault(tuple(profiles[eid]["charpoly"]), eid)
|
||
fps = list(fp_rep)
|
||
pts = [simplex_project(fp) for fp in fps]
|
||
delta = float(CARTAN_DELTA)
|
||
|
||
parent = list(range(len(fps)))
|
||
|
||
def find(i: int) -> int:
|
||
while parent[i] != i:
|
||
parent[i] = parent[parent[i]]
|
||
i = parent[i]
|
||
return i
|
||
|
||
sub_delta: List[dict] = []
|
||
min_d, min_pair = None, None
|
||
for i in range(len(fps)):
|
||
for j in range(i + 1, len(fps)):
|
||
d = fisher_distance(pts[i], pts[j])
|
||
if min_d is None or d < min_d:
|
||
min_d, min_pair = d, (fp_rep[fps[i]], fp_rep[fps[j]])
|
||
if d < delta:
|
||
sub_delta.append({
|
||
"representatives": [fp_rep[fps[i]], fp_rep[fps[j]]],
|
||
"fisher_distance": round(d, 8),
|
||
})
|
||
ri, rj = find(i), find(j)
|
||
if ri != rj:
|
||
parent[rj] = ri
|
||
components = len({find(i) for i in range(len(fps))})
|
||
return {
|
||
"delta_exact": "17/1792",
|
||
"delta_source": "CartanConnection.lean:70 (docs/cartan_fingerprint.md)",
|
||
"delta_float": delta,
|
||
"metric": "Fisher on Δ₇: d_F = 2·arccos(Σ√(p_i·q_i)); "
|
||
"p_i = |c_i|/Σ|c_j| (all-zero fingerprint → uniform, "
|
||
"by convention)",
|
||
"distinct_fingerprints": len(fps),
|
||
"pairs_below_delta": len(sub_delta),
|
||
"sub_delta_pairs": sub_delta[:50],
|
||
"delta_resolution_codewords": components,
|
||
"min_fisher_distance": round(min_d, 8) if min_d is not None else None,
|
||
"min_fisher_pair": list(min_pair) if min_pair else None,
|
||
}
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# Gap-aware quantization (deduped, sample-size guarded)
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
|
||
def gap_quantize(lambda_counts: Dict[float, int],
|
||
factor: float = GAP_FACTOR,
|
||
min_support: int = MIN_SUPPORT) -> dict:
|
||
"""Gap-aware boundaries over distinct-λ values with occupancy counts.
|
||
|
||
`lambda_counts` maps λ → number of DISTINCT matrices at that λ (dedupe
|
||
byte-identical matrices before calling; 12 of the 250 corpus matrices
|
||
are exact duplicates of another).
|
||
|
||
Candidate boundaries sit at midpoints of gaps > factor × median gap.
|
||
Sample-size guard: segments holding < min_support distinct matrices are
|
||
merged into the neighbor across the smaller bounding gap, so no cluster
|
||
is asserted on a handful of points. The suppressed high-end boundaries
|
||
define a sparse-tail cutoff: everything above the lowest suppressed
|
||
boundary with < min_support mass above it is an outlier region, not a
|
||
cluster.
|
||
"""
|
||
uniq = sorted(lambda_counts)
|
||
if len(uniq) < 2:
|
||
return {"boundaries": [], "median_gap": 0.0, "threshold": 0.0,
|
||
"candidates": [], "suppressed": [], "sparse_tail_cutoff": None}
|
||
gaps = [uniq[i + 1] - uniq[i] for i in range(len(uniq) - 1)]
|
||
sizes = sorted(gaps)
|
||
mid = len(sizes) // 2
|
||
median = sizes[mid] if len(sizes) % 2 else 0.5 * (sizes[mid - 1] + sizes[mid])
|
||
threshold = factor * median
|
||
# candidate boundary index i ↔ gap between uniq[i] and uniq[i+1]
|
||
candidates = [i for i, g in enumerate(gaps) if g > threshold]
|
||
|
||
def seg_counts(bnds: List[int]) -> List[int]:
|
||
counts, prev = [], 0
|
||
for b in bnds + [len(uniq) - 1]:
|
||
counts.append(sum(lambda_counts[uniq[k]] for k in range(prev, b + 1)))
|
||
prev = b + 1
|
||
return counts
|
||
|
||
kept = list(candidates)
|
||
suppressed: List[int] = []
|
||
while kept:
|
||
counts = seg_counts(kept)
|
||
weak = [s for s, c in enumerate(counts) if c < min_support]
|
||
if not weak:
|
||
break
|
||
s = min(weak, key=lambda i: counts[i])
|
||
# bounding boundaries of segment s in `kept`
|
||
options = []
|
||
if s > 0:
|
||
options.append(kept[s - 1])
|
||
if s < len(kept):
|
||
options.append(kept[s])
|
||
drop = min(options, key=lambda b: gaps[b]) # merge across smaller gap
|
||
kept.remove(drop)
|
||
suppressed.append(drop)
|
||
|
||
def midpoint(i: int) -> float:
|
||
return 0.5 * (uniq[i] + uniq[i + 1])
|
||
|
||
total = sum(lambda_counts.values())
|
||
sparse_cutoff: Optional[float] = None
|
||
for b in sorted(suppressed):
|
||
above = sum(lambda_counts[v] for v in uniq if v > midpoint(b))
|
||
if above < min_support:
|
||
sparse_cutoff = midpoint(b)
|
||
break
|
||
|
||
return {
|
||
"boundaries": [round(midpoint(b), 6) for b in sorted(kept)],
|
||
"median_gap": median,
|
||
"threshold": threshold,
|
||
"candidates": [round(midpoint(b), 6) for b in sorted(candidates)],
|
||
"suppressed": [round(midpoint(b), 6) for b in sorted(suppressed)],
|
||
"sparse_tail_cutoff": sparse_cutoff,
|
||
"total_support": total,
|
||
}
|
||
|
||
|
||
class Codebook:
|
||
"""Gap-aware spectral codebook: matrix → (codeword, index) → equation_id."""
|
||
|
||
def __init__(self, matrices: Dict[str, Matrix],
|
||
gap_factor: float = GAP_FACTOR,
|
||
min_support: int = MIN_SUPPORT) -> None:
|
||
self.matrices = matrices
|
||
self.gap_factor = gap_factor
|
||
self.min_support = min_support
|
||
self.profiles: Dict[str, dict] = {
|
||
eid: spectral_profile(mat) for eid, mat in matrices.items()
|
||
}
|
||
|
||
# Φ-corkscrew integer packing (corpus-wide base/offset, exact)
|
||
self.pack_base, self.pack_offset = packing_params(
|
||
[p["charpoly"] for p in self.profiles.values()])
|
||
self.spiral_indices: Dict[str, int] = {
|
||
eid: pack_spiral_index(p["charpoly"], self.pack_base, self.pack_offset)
|
||
for eid, p in self.profiles.items()
|
||
}
|
||
|
||
# dedupe byte-identical matrices BEFORE gap analysis
|
||
self.matrix_groups: Dict[Matrix, List[str]] = {}
|
||
for eid in sorted(self.profiles):
|
||
self.matrix_groups.setdefault(matrices[eid], []).append(eid)
|
||
lambda_counts: Dict[float, int] = {}
|
||
for mat, ids in self.matrix_groups.items():
|
||
lam = self.profiles[ids[0]]["spectral_radius"]
|
||
lambda_counts[lam] = lambda_counts.get(lam, 0) + 1
|
||
|
||
self.quant = gap_quantize(lambda_counts, gap_factor, min_support)
|
||
self.boundaries: List[float] = self.quant["boundaries"]
|
||
|
||
# cluster assignment + deterministic within-cluster rank
|
||
by_cluster: Dict[int, List[str]] = {}
|
||
for eid, p in self.profiles.items():
|
||
by_cluster.setdefault(self.cluster_of(p["spectral_radius"]), []).append(eid)
|
||
self._assign: Dict[str, Tuple[str, int]] = {}
|
||
self._decode: Dict[Tuple[str, int], str] = {}
|
||
self._by_matrix: Dict[Matrix, Tuple[str, int]] = {}
|
||
for ci in sorted(by_cluster):
|
||
members = sorted(
|
||
by_cluster[ci],
|
||
key=lambda e: (
|
||
self.profiles[e]["spectral_radius"],
|
||
self.profiles[e]["charpoly"],
|
||
self.matrices[e],
|
||
e,
|
||
),
|
||
)
|
||
for idx, eid in enumerate(members):
|
||
cw = f"C{ci}"
|
||
self._assign[eid] = (cw, idx)
|
||
self._decode[(cw, idx)] = eid
|
||
# first (canonical) entry wins for byte-identical duplicates
|
||
self._by_matrix.setdefault(self.matrices[eid], (cw, idx))
|
||
|
||
def cluster_of(self, lam: float) -> int:
|
||
return sum(1 for b in self.boundaries if lam > b)
|
||
|
||
def is_sparse_tail(self, lam: float) -> bool:
|
||
cut = self.quant["sparse_tail_cutoff"]
|
||
return cut is not None and lam > cut
|
||
|
||
def codeword_of(self, equation_id: str) -> Tuple[str, int]:
|
||
return self._assign[equation_id]
|
||
|
||
def encode(self, mat: Sequence[Sequence[int]]) -> Tuple[str, Optional[int]]:
|
||
"""(codeword, index) for a matrix. Known matrices get their exact
|
||
index (canonical entry for byte-identical duplicates); novel
|
||
matrices get a cluster codeword with index None."""
|
||
key = tuple(tuple(int(x) for x in row) for row in mat)
|
||
if key in self._by_matrix:
|
||
return self._by_matrix[key]
|
||
return f"C{self.cluster_of(spectral_radius(key))}", None
|
||
|
||
def decode(self, codeword: str, index: int) -> str:
|
||
return self._decode[(codeword, index)]
|
||
|
||
# ── reporting ────────────────────────────────────────────────────
|
||
|
||
def cluster_table(self) -> List[dict]:
|
||
clusters: Dict[str, List[str]] = {}
|
||
for eid, (cw, _) in self._assign.items():
|
||
clusters.setdefault(cw, []).append(eid)
|
||
rows = []
|
||
for cw in sorted(clusters, key=lambda c: int(c[1:])):
|
||
ids = clusters[cw]
|
||
lams = [self.profiles[e]["spectral_radius"] for e in ids]
|
||
distinct = len({self.matrices[e] for e in ids})
|
||
rows.append({
|
||
"codeword": cw,
|
||
"count": len(ids),
|
||
"distinct_matrices": distinct,
|
||
"lambda_min": round(min(lams), 6),
|
||
"lambda_max": round(max(lams), 6),
|
||
"sparse_tail_members": sum(1 for l in lams if self.is_sparse_tail(l)),
|
||
})
|
||
return rows
|
||
|
||
def collision_report(self) -> dict:
|
||
lam_groups: Dict[float, List[str]] = {}
|
||
cp_groups: Dict[Tuple[int, ...], List[str]] = {}
|
||
for eid, p in self.profiles.items():
|
||
lam_groups.setdefault(round(p["spectral_radius"], LAMBDA_DP), []).append(eid)
|
||
cp_groups.setdefault(tuple(p["charpoly"]), []).append(eid)
|
||
|
||
def collide(groups): # groups with >1 member
|
||
return {k: sorted(v) for k, v in groups.items() if len(v) > 1}
|
||
|
||
lam_c, cp_c = collide(lam_groups), collide(cp_groups)
|
||
mat_c = {m: ids for m, ids in self.matrix_groups.items() if len(ids) > 1}
|
||
# cospectral = same char poly but NOT byte-identical matrices
|
||
dup_sets = [frozenset(v) for v in mat_c.values()]
|
||
cospectral = {
|
||
",".join(map(str, k)): v
|
||
for k, v in cp_c.items() if frozenset(v) not in dup_sets
|
||
}
|
||
pi_mismatch = [
|
||
eid for eid, p in self.profiles.items()
|
||
if abs(p["spectral_radius"] - p["lambda_power_iteration"]) > 1e-3
|
||
]
|
||
n = len(self.profiles)
|
||
return {
|
||
"corpus_size": n,
|
||
"distinct_matrices": len(self.matrix_groups),
|
||
"duplicate_matrix_groups": len(mat_c),
|
||
"duplicate_matrix_members": sum(len(v) for v in mat_c.values()),
|
||
"duplicate_matrix_classes": [sorted(v) for v in mat_c.values()],
|
||
"unique_lambda_4dp": len(lam_groups),
|
||
"lambda_collision_groups": len(lam_c),
|
||
"lambda_collision_members": sum(len(v) for v in lam_c.values()),
|
||
"unique_charpoly": len(cp_groups),
|
||
"charpoly_collision_groups": len(cp_c),
|
||
"charpoly_collision_members": sum(len(v) for v in cp_c.values()),
|
||
"charpoly_collision_classes": {
|
||
",".join(map(str, k)): v for k, v in cp_c.items()
|
||
},
|
||
"cospectral_nonidentical_groups": len(cospectral),
|
||
"cospectral_nonidentical": cospectral,
|
||
"identified_by_lambda_4dp": sum(1 for v in lam_groups.values() if len(v) == 1),
|
||
"identified_by_charpoly": sum(1 for v in cp_groups.values() if len(v) == 1),
|
||
"power_iteration_artifacts": sorted(pi_mismatch),
|
||
"power_iteration_artifact_count": len(pi_mismatch),
|
||
}
|
||
|
||
def to_json(self) -> dict:
|
||
entries = []
|
||
for eid in sorted(self.profiles):
|
||
cw, idx = self._assign[eid]
|
||
p = self.profiles[eid]
|
||
n = self.spiral_indices[eid]
|
||
entry = {
|
||
"equation_id": eid,
|
||
"codeword": cw,
|
||
"index": idx,
|
||
"spectral_radius": p["spectral_radius"],
|
||
"lambda_power_iteration": p["lambda_power_iteration"],
|
||
"charpoly": p["charpoly"],
|
||
"density": p["density"],
|
||
"frobenius": p["frobenius"],
|
||
"trace": p["trace"],
|
||
"spiral_index": n,
|
||
"layout": {
|
||
"radius_sq": n,
|
||
"angle_frac": round(angle_frac(n), 12),
|
||
},
|
||
}
|
||
if self.is_sparse_tail(p["spectral_radius"]):
|
||
entry["sparse_tail"] = True
|
||
twins = [t for t in self.matrix_groups[self.matrices[eid]] if t != eid]
|
||
if twins:
|
||
entry["identical_matrix_ids"] = twins
|
||
entries.append(entry)
|
||
return {
|
||
"schema": "spectral_codebook_v2",
|
||
"generated": date.today().isoformat(),
|
||
"source": "formal/SilverSight/PIST/Matrices250.lean",
|
||
"matrix_count": len(entries),
|
||
"lambda_semantics": (
|
||
"spectral_radius is the exact max root modulus of the integer "
|
||
"characteristic polynomial (charpoly); lambda_power_iteration is "
|
||
"the legacy estimate and is unreliable on peripheral spectra"
|
||
),
|
||
"phi_corkscrew": {
|
||
"packing": {
|
||
"base": self.pack_base,
|
||
"offset": self.pack_offset,
|
||
"digits": 8,
|
||
"order": "little-endian: charpoly c1 is digit 0",
|
||
"semantics": (
|
||
"spiral_index = Σ (c_i + offset)·base^i — exact "
|
||
"integer, injective on charpoly fingerprints "
|
||
"(unpack_spiral_index inverts it exactly); inherits "
|
||
"charpoly's non-injectivity on matrices, so identity "
|
||
"still needs the within-cluster index"
|
||
),
|
||
},
|
||
"layout": {
|
||
"psi": "2π/φ² (golden angle)",
|
||
"angle_frac": "frac(n·(3−√5)/2), decimal-exact precision",
|
||
"min_angular_separation_frac": (
|
||
min_angular_separation(list(self.spiral_indices.values()))
|
||
),
|
||
"semantics": (
|
||
"decorative placement only: f(n) is injective because "
|
||
"1/φ² is irrational, and carries no information beyond "
|
||
"spiral_index"
|
||
),
|
||
},
|
||
},
|
||
"quantization": {
|
||
"method": "gap_aware_deduped_min_support",
|
||
"gap_factor": self.gap_factor,
|
||
"min_support": self.min_support,
|
||
"median_gap": round(self.quant["median_gap"], 6),
|
||
"threshold": round(self.quant["threshold"], 6),
|
||
"boundaries": self.boundaries,
|
||
"candidate_boundaries": self.quant["candidates"],
|
||
"suppressed_boundaries": self.quant["suppressed"],
|
||
"sparse_tail_cutoff": self.quant["sparse_tail_cutoff"],
|
||
},
|
||
"clusters": self.cluster_table(),
|
||
"collisions": self.collision_report(),
|
||
"entries": entries,
|
||
}
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# 278-row manifold cross-check
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
|
||
def manifold_check(codebook: Codebook,
|
||
manifold_path: Path = MANIFOLD_LEAN) -> dict:
|
||
"""Run the pipeline over the full 278-row fixture corpus.
|
||
|
||
The manifold repeats equation ids (278 rows over 250 unique ids), so the
|
||
28 extra rows reuse matrices already in the codebook. After dedupe the
|
||
λ multiset is unchanged, so boundaries must match the 250-matrix
|
||
codebook exactly; this is reported rather than assumed.
|
||
"""
|
||
ids = parse_manifold_ids(manifold_path)
|
||
missing = sorted({e for e in ids if e not in codebook.profiles})
|
||
seen: Dict[Matrix, str] = {}
|
||
lambda_counts: Dict[float, int] = {}
|
||
for eid in ids:
|
||
if eid in codebook.profiles:
|
||
mat = codebook.matrices[eid]
|
||
if mat not in seen: # dedupe, matching the codebook build
|
||
seen[mat] = eid
|
||
lam = codebook.profiles[eid]["spectral_radius"]
|
||
lambda_counts[lam] = lambda_counts.get(lam, 0) + 1
|
||
quant = gap_quantize(lambda_counts, codebook.gap_factor, codebook.min_support)
|
||
return {
|
||
"rows": len(ids),
|
||
"unique_ids": len(set(ids)),
|
||
"extra_rows": len(ids) - len(set(ids)),
|
||
"ids_without_matrix": missing,
|
||
"distinct_matrices": len(seen),
|
||
"boundaries": quant["boundaries"],
|
||
"median_gap": round(quant["median_gap"], 6),
|
||
"boundaries_match_250": quant["boundaries"] == codebook.boundaries,
|
||
"clusters": len(quant["boundaries"]) + 1,
|
||
}
|
||
|
||
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
# CLI
|
||
# ──────────────────────────────────────────────────────────────────────
|
||
|
||
def build_codebook(matrices_path: Path = MATRICES_LEAN,
|
||
gap_factor: float = GAP_FACTOR,
|
||
min_support: int = MIN_SUPPORT) -> Codebook:
|
||
return Codebook(parse_matrices_lean(matrices_path), gap_factor, min_support)
|
||
|
||
|
||
def main(argv: Optional[List[str]] = None) -> int:
|
||
ap = argparse.ArgumentParser(description="Gap-aware spectral codebook generator")
|
||
ap.add_argument("--matrices", type=Path, default=MATRICES_LEAN,
|
||
help="Matrices250.lean path")
|
||
ap.add_argument("--out", type=Path, default=DEFAULT_OUT,
|
||
help="Output codebook JSON (default data/spectral_codebook.json)")
|
||
ap.add_argument("--gap-factor", type=float, default=GAP_FACTOR,
|
||
help="Boundary threshold = factor × median gap (default 3)")
|
||
ap.add_argument("--min-support", type=int, default=MIN_SUPPORT,
|
||
help="Minimum distinct matrices per cluster (default 10)")
|
||
ap.add_argument("--check-manifold", action="store_true",
|
||
help="Also run over the 278-row RRC/Q16_16Manifold corpus")
|
||
ap.add_argument("--cartan-floor", action="store_true",
|
||
help="Also emit the Cartan ∆=17/1792 Fisher-floor report "
|
||
"(alongside, not replacing, the gap rule)")
|
||
ap.add_argument("--no-write", action="store_true",
|
||
help="Report only; do not write the JSON")
|
||
ap.add_argument("--sync-db", action="store_true",
|
||
help="DRY-RUN sync preview to ene.rrc_predictions on the "
|
||
"neon-64gb Postgres ($NEON_PG). Always a dry run from "
|
||
"this entry point; use python/spectral_codebook_db.py "
|
||
"--apply to actually insert.")
|
||
args = ap.parse_args(argv)
|
||
|
||
cb = build_codebook(args.matrices, args.gap_factor, args.min_support)
|
||
doc = cb.to_json()
|
||
|
||
col = doc["collisions"]
|
||
q = doc["quantization"]
|
||
print(f"Spectral codebook over {doc['matrix_count']} matrices "
|
||
f"({col['distinct_matrices']} distinct)")
|
||
print(f" median gap {q['median_gap']:.6f} × {args.gap_factor:g} → "
|
||
f"threshold {q['threshold']:.6f}; min support {args.min_support} "
|
||
f"(candidates {len(q['candidate_boundaries'])}, "
|
||
f"kept {len(q['boundaries'])}, "
|
||
f"suppressed {len(q['suppressed_boundaries'])})")
|
||
print(f" boundaries: {q['boundaries']}")
|
||
print(f" sparse-tail cutoff: {q['sparse_tail_cutoff']}")
|
||
print("\n codeword ids distinct lambda range sparse-tail")
|
||
for row in doc["clusters"]:
|
||
print(f" {row['codeword']:<8} {row['count']:>4} {row['distinct_matrices']:>8}"
|
||
f" [{row['lambda_min']:.4f}, {row['lambda_max']:.4f}]"
|
||
f" {row['sparse_tail_members'] or ''}")
|
||
print(f"\n fingerprint uniqueness (of {col['corpus_size']} ids / "
|
||
f"{col['distinct_matrices']} distinct matrices):")
|
||
print(f" λ @4dp : {col['unique_lambda_4dp']} unique values, "
|
||
f"{col['identified_by_lambda_4dp']} ids uniquely identified, "
|
||
f"{col['lambda_collision_groups']} collision groups "
|
||
f"({col['lambda_collision_members']} ids)")
|
||
print(f" char-poly ℤ : {col['unique_charpoly']} unique fingerprints, "
|
||
f"{col['identified_by_charpoly']} ids uniquely identified, "
|
||
f"{col['charpoly_collision_groups']} collision groups "
|
||
f"({col['charpoly_collision_members']} ids)")
|
||
print(f" exact duplicate matrices: {col['duplicate_matrix_groups']} groups "
|
||
f"({col['duplicate_matrix_members']} ids)")
|
||
print(f" cospectral non-identical groups: "
|
||
f"{col['cospectral_nonidentical_groups']}")
|
||
print(f" power-iteration artifacts (|λ_pi − ρ| > 1e-3): "
|
||
f"{col['power_iteration_artifact_count']}")
|
||
|
||
if args.check_manifold:
|
||
mc = manifold_check(cb)
|
||
doc["manifold_278"] = mc
|
||
print(f"\n 278-row manifold: {mc['rows']} rows / {mc['unique_ids']} unique ids "
|
||
f"({mc['extra_rows']} duplicate rows, "
|
||
f"{mc['distinct_matrices']} distinct matrices)")
|
||
print(f" boundaries match 250-matrix codebook: {mc['boundaries_match_250']} "
|
||
f"({mc['clusters']} clusters)")
|
||
|
||
if args.cartan_floor:
|
||
cf = cartan_floor_report(cb.profiles)
|
||
doc["cartan_floor"] = cf
|
||
gap_clusters = len(q["boundaries"]) + 1
|
||
print(f"\n Cartan ∆-floor (∆ = {cf['delta_exact']} ≈ "
|
||
f"{cf['delta_float']:.6f}, {cf['delta_source']}):")
|
||
print(f" distinct fingerprints on Δ₇ : {cf['distinct_fingerprints']}")
|
||
print(f" pairs below ∆ (Fisher) : {cf['pairs_below_delta']}")
|
||
print(f" min Fisher distance : {cf['min_fisher_distance']} "
|
||
f"{cf['min_fisher_pair']}")
|
||
print(f" rule comparison — codewords : "
|
||
f"∆-floor {cf['delta_resolution_codewords']} vs "
|
||
f"3×-median-gap {gap_clusters} (λ clusters); both emitted, "
|
||
f"neither replaces the other")
|
||
|
||
if args.sync_db:
|
||
import spectral_codebook_db
|
||
print("\n DB sync (dry run):")
|
||
spectral_codebook_db.sync_db(codebook=cb) # never applies from here
|
||
|
||
if not args.no_write:
|
||
args.out.parent.mkdir(parents=True, exist_ok=True)
|
||
args.out.write_text(json.dumps(doc, indent=2) + "\n")
|
||
print(f"\nWrote {args.out}")
|
||
return 0
|
||
|
||
|
||
if __name__ == "__main__":
|
||
sys.exit(main())
|