diff --git a/research/a1-isotypic-128-charpoly/PROOF.md b/research/a1-isotypic-128-charpoly/PROOF.md new file mode 100644 index 0000000000..e598bd522d --- /dev/null +++ b/research/a1-isotypic-128-charpoly/PROOF.md @@ -0,0 +1,166 @@ +# Characteristic polynomial on the 128-dimensional isotypic component + +This is an exact computation in linear algebra over `Q`. It determines one +monic characteristic polynomial on one invariant subspace of the historical +`c=1` centre block. The selected conditional kappa stays +`781451/250000000`. The degree-881 minimal polynomial of the full one-slot +operator is not replayed. + +## Operator and group + +The input is the 8-by-112 query in `INTERFACE.json`, SHA256 +`5c6c51fbffc07a683e922f17532fd90262ae9e23cc3560a696002073dd234ef3`. +Its columns lie in `{-1,0,1}^8`, each has weight 4, the 112 port labels are +distinct odd-weight vectors in `F_2^9`, and `R R^T = 56 I_8`. + +On `Q^8 ⊗ Q[F_2^9]`, of dimension 4096, + +```text +A_1 = (1/4) sum_q r_q r_q^T ⊗ P_q(1), +B = 4 A_1 = sum_q r_q r_q^T ⊗ P_q(1). +``` + +`P_q(1)` sends `e_a` to `e_{a+q}` when the least bit of `q` is 1, and to +`(-1)^{q·a} e_a` otherwise. The basis index `i*512+a` is query coordinate +`i` and address `a`. + +`E` is the joint `(+1,+1)` space of the order-16 sign stabilizer of sector 1 +and the involutions `g` and `t`. The trace formula gives dimension 216. The +recomputed signed-orbit basis equals the basis pinned by +`research/a1-joint-216-charpoly/evidence/basis.json`. On that basis, the +port formula for `B` reproduces `evidence/matrix.json`. + +The dihedral image `G` has order 12. Its generators are the signed +permutations + +| name | `inv` | output signs `-1` | +|---|---|---| +| `u` | `(0, 1, 2, 3, 7, 6, 5, 4)` | none | +| `r` | `(0, 1, 6, 7, 3, 2, 5, 4)` | coordinates 0 and 7 | + +`u^2 = 1`, `r^6 = 1`, `r` has order 6, and `u r u = r^{-1}`. Both preserve +sector 1 and `E`. Each power of `r` acts on the pinned basis of `E` by a +signed permutation. `B` commutes with `r` and with `u` on the ambient +4096-dimensional space. + +## The isotypic component + +Let `ψ` be the rational character with values `2, 1, -1, -2, -1, 1` on +`r^k` for `k = 0..5` and value 0 on the reflections. These are the traces of +the rational representation whose rotation satisfies `x^2 - x + 1 = 0`. The +central idempotent is + +```text +P = (1/6) sum_{k=0}^{5} ψ(r^k) ρ(r^k), +6P = 2 I + ρ(r) - ρ(r^2) - 2 ρ(r^3) - ρ(r^4) + ρ(r^5). +``` + +On the basis of `E`, `6P` has trace 768, rank 128, and `(6P)^2 = 6(6P)`. +The rational image `W` therefore has dimension 128, and the multiplicity of +`ψ` is `128/2 = 64`. The other two-dimensional rational character, +`2, -1, -1, 2, -1, -1` on `r^k`, has trace 0 on `E`. + +`B` commutes with `r`, so `B` preserves the column lattice `L = (6P) Z^216`. +A unimodular Hermite process produces a `Z`-basis of `L`. Its 128 diagonals +lie in `{2, 6}`, so `L` is a finite-index sublattice of the saturated +lattice `{v in Z^216 : (6P)v = 6v}`. Every Hermite basis vector satisfies +`(6P)v = 6v`, and there are 128 of them, so they are a basis of `W` over +`Q`. The matrix `X` of `B` on this basis is an integer matrix: each image +`B v` reduces through the Hermite pivots with exact divisions. + +## Polynomial + +The monic characteristic polynomial is `chi_B(t) = det(t I - X)`. It lies in +`Z[t]`, every coefficient of `t^k` is divisible by `4^{128-k}`, and + +```text +chi_{A_1, W}(x) = 4^{-128} chi_B(4x) +``` + +is monic in `Z[x]`. Both polynomials are even. Evaluations of `chi_B` at +`0, 1, 2, 3, 5` agree with `det(t I - X)` modulo the first three CRT primes. + +Write + +```text +f(x) = x^11 + 2x^10 - 108x^9 - 192x^8 + 3208x^7 + 2128x^6 - 33888x^5 + 19264x^4 + 115776x^3 - 196352x^2 + 102144x - 12288 +g(x) = x^11 - 2x^10 - 108x^9 + 192x^8 + 3208x^7 - 2128x^6 - 33888x^5 - 19264x^4 + 115776x^3 + 196352x^2 + 102144x + 12288 +h(x) = x^14 - 2x^13 - 112x^12 + 136x^11 + 4376x^10 - 3632x^9 - 76384x^8 + 44928x^7 + 621184x^6 - 260736x^5 - 2179840x^4 + 704000x^3 + 3049472x^2 - 702464x - 1167360 +k(x) = x^14 + 2x^13 - 112x^12 - 136x^11 + 4376x^10 + 3632x^9 - 76384x^8 - 44928x^7 + 621184x^6 + 260736x^5 - 2179840x^4 - 704000x^3 + 3049472x^2 + 702464x - 1167360 +``` + +Then `g(x) = -f(-x)`, `k(x) = h(-x)`, and each of `f, g, h, k` is irreducible +over `Z` by the Rabin test modulo one prime. The factorization on `W` is + +```text +chi_{A_1, W}(x) = q(x)^2, +q(x) = x^6 (x-2)^4 (x+2)^4 f(x) g(x) h(x) k(x). +``` + +Equivalently, `q(x) = x^6 r(x^2)` where `r` is the monic degree-29 polynomial +with ascending coefficients + +```text +-52675902681576269414400 +2303306134076870941999104 +-18827961748477891024257024 +73124429304073607471693824 +-171456969792808134744473600 +271456076002980399446753280 +-309537536542443194852311040 +264951255789256454094979072 +-175106550821218960407527424 +91149737035716567039475712 +-37905952006491036612296704 +12722870692816032484032512 +-3471386398306750712250368 +773676955090062483128320 +-141267596592764868886528 +21162534505722224312320 +-2601171087886668267520 +262024005207787307008 +-21579308813538295808 +1447508150155345920 +-78661706621779968 +3437815337254912 +-119644063870976 +3271915571200 +-69036760064 +1095226368 +-12571072 +98112 +-464 +1 +``` + +The product of these factors is the expanded polynomial. Ascending-coefficient +SHA256 values, as canonical JSON, are + +```text +chi_A bf10fad08aa427ff461bc3a7a7377f5b90ab7a810a1004d81e550a2a328bbb65 +chi_B 33f15f1bcc5d9a57fec18f40653e2362792a77b2c82bdbcf4bbd469694dd03a2 +q d895a734501995193ce1a820dfe8d7cf9121fdf12a1abd1620a344e3315f5762 +``` + +Each factor above occurs in the published degree-216 factorization from +`research/a1-joint-216-charpoly/`, and the multiplicity on `W` does not exceed +the multiplicity on that 216-dimensional space. The factors which occur there +with multiplicity one do not divide `chi_{A_1, W}`. That comparison reads the +published coefficient lists. It does not rebuild the degree-216 matrix +characteristic polynomial. + +The square is proved by exhibiting `q` in `Z[x]`. It is the polynomial +predicted by Schur's lemma for this character: the representation is +absolutely irreducible, so an operator commuting with `G` acts as `I_2 ⊗ T` +on a 64-dimensional multiplicity space, and `q` is the characteristic +polynomial of `T`. This package does not construct a basis in which that +tensor form is visible. + +## Evidence class + +Finite exact linear algebra over `Q`. The verifier rebuilds the query, the +group, the basis of `E`, the operator `6P`, a `Z`-basis of its column +lattice, and the characteristic polynomial, then checks the factorization, +the square, the determinant samples, and the rejection controls. A one-entry +change of `X`, a dropped character coefficient, a dropped splitting +involution, and several nearby permutations fail those controls. diff --git a/research/a1-isotypic-128-charpoly/README.md b/research/a1-isotypic-128-charpoly/README.md new file mode 100644 index 0000000000..bbb045bfd4 --- /dev/null +++ b/research/a1-isotypic-128-charpoly/README.md @@ -0,0 +1,99 @@ +# Characteristic polynomial of one 128-dimensional isotypic component + +Prepared by unpaidmoon with substantial Cursor Grok assistance. Apache-2.0. + +For the historical 8-by-112 centre query, let `E` be the 216-dimensional +joint `(+1,+1)` space of the `c=1` sign stabilizer and the involutions `g` +and `t`. Let `G` be the order-12 dihedral image generated by the signed +permutations `r` and `u` inside that space. The rational character with +values `2, 1, -1, -2, -1, 1` on `r^k` and value 0 on the reflections cuts +out a 128-dimensional `A_1`-invariant subspace `W`. Its monic characteristic +polynomial is the square recorded in [PROOF.md](PROOF.md). + +The selected conditional kappa stays `781451/250000000`. + +## Result + +`B = 4 A_1` preserves `W`. On a Hermite `Z`-basis of the column lattice of +the integer projector `6P`, + +```text +chi_{A_1, W}(x) = q(x)^2, +q(x) = x^6 (x-2)^4 (x+2)^4 f(x) g(x) h(x) k(x), +``` + +where `f` and `g` are the published degree-11 pair and `h` and `k` are the +published degree-14 pair that already appears squared in the degree-216 +factorization. `g(x) = -f(-x)` and `k(x) = h(-x)`. The four higher factors +are irreducible over `Z`. Ascending-coefficient SHA256 of `chi_{A_1, W}` is + +```text +bf10fad08aa427ff461bc3a7a7377f5b90ab7a810a1004d81e550a2a328bbb65 +``` + +Eigenvalue 0 has algebraic multiplicity 12 on `W`. The multiplicity-one +factors of the published degree-216 polynomial do not divide this +characteristic polynomial. + +## Scope + +Catalogue card `formal-interface-lemma` asks for a formal caller binding. +The computation answers the second open question of PR #7023: the +characteristic polynomial of `A_1` on this 128-dimensional component. The +degree-216 polynomial and basis are the comparison data in +`research/a1-joint-216-charpoly/`. This package recomputes `W` from the +pinned query. + +Selected-baseline orientation returned `PASS_RESEARCH_BASELINE_BINDINGS_ONLY` +for kappa `781451/250000000`. + +## Reproduce + +From the repository root, with CPython 3.12 and the standard library: + +```sh +python3 -B -W error research/a1-isotypic-128-charpoly/verify.py +python3 -B -O research/a1-isotypic-128-charpoly/verify.py +python3 -B scripts/research_package.py --package research/a1-isotypic-128-charpoly --check-only +``` + +The marker is `PASS_A1_ISOTYPIC_128_CHARPOLY`. Acceptance uses explicit +checks, so `-O` does not skip a theorem or a rejection control. No Lean, +Docker, degree-881 replay, or selected cold replay was run. + +## Overlap + +`scripts/find_overlap.py --offline` was run on 2026-11-10 against the +snapshot fetched at 10:16:36Z for open PRs and 10:17:43Z for discussions. +Coverage: local notes 10180 paths against HEAD +`17469fff78c22941280f04d39dc5a1faee09c840` (10178 committed paths, complete); +339 open PRs, pagination complete; 7013 discussion records, pagination +complete. The phrases `128-dimensional isotypic`, `two-dimensional rational`, +`isotypic component` and `multiplicity 64` matched PR #7023 and discussion +comment +[19514162](https://github.com/CrocSwap/integer-mult-bounds/discussions/1166#discussioncomment-19514162). +Both state this polynomial as remaining work. A direct open-PR search for +`128-dimensional isotypic` returned only PR #7023. Closed work was searched +through local READMEs, TASK metadata and the ledger. + +## Open next questions + +The structured list is in [TASK.json](TASK.json). + +- Eigenvalue 0 has algebraic multiplicity 12 on `W`. The geometric + multiplicity, and the ranks of `B ∓ 8I` for eigenvalues `±2`, are not + computed. +- `q` is the characteristic polynomial of the multiplicity-space operator + `T`. An explicit integer matrix of `T` is not constructed. + +## Credit + +JoyboyBrian recorded the sign-character dimensions. PR #6944 supplied the +216-dimensional blocks. PR #6964 computed the degree-216 characteristic +polynomial, the joint basis and the matrix of `B`. PR #7008 identified the +signed permutations `r` and `u` and the multiplicity 64. PR #7023 proved the +8-dimensional summand is reducible and asked for this polynomial. Ali +Romman's PR #4958 and PR #4720 retain the query and the Bell action. Nick +Petalas's PR #5120 is the degree-410 lower bound. Krzysztof Wojt's +`research/exact-one-slot-gram/` is the degree-881 annihilator. The doubled +representation comes from PR #1719. New work is Apache-2.0. diff --git a/research/a1-isotypic-128-charpoly/SOURCES.json b/research/a1-isotypic-128-charpoly/SOURCES.json new file mode 100644 index 0000000000..f6d600f7e9 --- /dev/null +++ b/research/a1-isotypic-128-charpoly/SOURCES.json @@ -0,0 +1,36 @@ +{ + "schema_version": 1, + "base_commit": "17469fff78c22941280f04d39dc5a1faee09c840", + "files": [ + { + "path": "certificates/selected-result.json", + "usage": "provenance", + "git_blob": "d59b5035e7b4cbaa9d688a1d32dd0d6eb79c466a", + "sha256": "9b8d0eea95b3bf0dfcb7938a732a6439b87a2685dec2a5ab9e525b3b34debbfa" + }, + { + "path": "docs/research/agent-task-catalog.json", + "usage": "provenance", + "git_blob": "3404fbdd5e7529a0674eef26114b5071e4586cec", + "sha256": "ead6cdc9d0af6d58e8a49744d8ec8c7c15e13d2ba93c276c434bcf9b51b93a81" + }, + { + "path": "research/a1-joint-216-charpoly/evidence/basis.json", + "usage": "input", + "git_blob": "5e524a41ed1c65017b7d902b5445aedd2eb606ff", + "sha256": "e198cd71f4100d2629ecc25fb40b135a8c18c7756be8bc80bf24d7b349af418a" + }, + { + "path": "research/a1-joint-216-charpoly/evidence/matrix.json", + "usage": "input", + "git_blob": "23bccaef6f24e6e644bdba81f58e9217787ab16c", + "sha256": "44f40ad5a12daf7c8493b42e1e9ddc4ef338b317ee5fdfc6de1b413642ca97aa" + }, + { + "path": "research/centre-polynomial-cache-limits/frozen/work/centre-encoding-campaign-20261014/lanes/circuit-interface/INTERFACE.json", + "usage": "input", + "git_blob": "4201e9a04c12273b919e2cae60a6b6482c8150f3", + "sha256": "5c6c51fbffc07a683e922f17532fd90262ae9e23cc3560a696002073dd234ef3" + } + ] +} diff --git a/research/a1-isotypic-128-charpoly/TASK.json b/research/a1-isotypic-128-charpoly/TASK.json new file mode 100644 index 0000000000..15220c2bb2 --- /dev/null +++ b/research/a1-isotypic-128-charpoly/TASK.json @@ -0,0 +1,118 @@ +{ + "schema_version": 1, + "task_id": "formal-interface-lemma", + "priority_revision": 8, + "base_commit": "17469fff78c22941280f04d39dc5a1faee09c840", + "original_card": { + "id": "formal-interface-lemma", + "question": "Can one uncovered interface implication be formalized with exact caller and certificate bindings?", + "starting_points": [ + "docs/research/research-family-status.md", + "certificates/selected-result.json", + "research/gate-operand-parity-lean/README.md", + "research/staged544-certified-endpoint/README.md", + "research/unconditional-batched-crt/README.md" + ], + "methods": [ + "Formalize a missing boundary lemma identified by the selected review.", + "Audit a formal theorem-to-certificate binding independently." + ], + "acceptance_gate": "State the exact formalized proposition and remaining assumptions; avoid duplicate arithmetic formalizations already present.", + "first_step": { + "command": [ + "python3", + "-B", + "scripts/research_selected_baseline.py" + ], + "expected": "PASS_RESEARCH_BASELINE_BINDINGS_ONLY", + "inspect": "Resolve the current selected review, separate complex/ordinary BIT/admitted BIT/outer budgets, and compare the card against live handoffs before starting. Historical ceilings require exact source/model equivalence.", + "scope": "Source-binding orientation only; not a scientific replay or all-size proof. Materialize selected inputs if this checkout is sparse." + }, + "first_experiment": "Write one missing statement with quantified hypotheses, actual caller and existing coverage. A compact-switch lemma is a candidate only after its mathematical contract is reviewed. Establish a bounded proof plan before a separately budgeted proof-assistant build; do not rerun large Lean jobs as orientation.", + "bottleneck_link": "A missing interface proof can reduce the conditional gap even if it does not change kappa.", + "allowed_variations": [ + "Narrow to one formal lemma or one theorem-to-certificate binding." + ], + "scope_change_triggers": [ + "Formalizing already-covered arithmetic and claiming the full multiplication theorem." + ], + "overlap_terms": [ + "hyperplane end", + "formal interface", + "Lean binding", + "3807", + "3857", + "4100", + "gate-operand-parity-lean", + "staged544-certified-endpoint", + "unconditional-batched-crt" + ], + "existing_checks": [ + { + "scope": "PR #3807: scoped gate-operand parity formalization.", + "evidence": "research/gate-operand-parity-lean/README.md", + "remaining": "Exact theorem-to-caller binding and genuinely uncovered interface lemmas.", + "status": "author-checked", + "source_commit": "4c3530f42a6abcc48e7a6bcffd6c901911dc8269", + "selected_record_sha256": "bb2eaa58647b353f88d0aec2eafb9d30df1be53589231b23e12ceef3e4bad129" + }, + { + "scope": "PR #3857: certified real-power moment/endpoint arithmetic.", + "evidence": "research/staged544-certified-endpoint/README.md", + "remaining": "Uniform algorithmic implications, not another endpoint digit.", + "status": "author-checked", + "source_commit": "ef4315b8caa605e47a2515004845eb2d2ed71a86", + "selected_record_sha256": "bb2eaa58647b353f88d0aec2eafb9d30df1be53589231b23e12ceef3e4bad129" + }, + { + "scope": "PR #4100: compact-switch formalization proposed as future validation.", + "evidence": "research/unconditional-batched-crt/README.md", + "remaining": "Audit the statement and assumptions before a formal proof campaign.", + "status": "author-checked", + "source_commit": "337d04939571a637aeab54a0f58f54fde811663a", + "selected_record_sha256": "bb2eaa58647b353f88d0aec2eafb9d30df1be53589231b23e12ceef3e4bad129" + } + ] + }, + "scope_status": "changed", + "actual_scope": "Characteristic polynomial, over Q, of A_1 on the 128-dimensional isotypic component of the two-dimensional rational representation of the dihedral image inside the c=1 joint (+1,+1) space.", + "scope_change_reason": "The catalogue card asks for a formal caller or certificate binding. PR #7023 asks for the characteristic polynomial of A_1 on the 128-dimensional isotypic component of the same dihedral image, and that polynomial is not on main. This package computes it. It does not add a Lean binding.", + "bottleneck_link": "Scoped side result on the historical centre Gram. The selected conditional kappa remains 781451/250000000. The one-slot minimal polynomial remains degree 881 in research/exact-one-slot-gram/; this polynomial locates one invariant subspace of the c=1 block and does not replace that certificate.", + "budget_minutes": 60, + "stop_rule": "Stop once a basis of the 128-dimensional component and the monic characteristic polynomial of A_1 on it are recomputed from the pinned query, with the square and the comparison to the published degree-216 factors checked. Do not rerun the degree-881 certificate.", + "outcome": "completed", + "evidence": "verify.py rebuilds the joint basis, the integer projector 6P, a Hermite basis of its column lattice, and the characteristic polynomial. Marker PASS_A1_ISOTYPIC_128_CHARPOLY. chi_A SHA256 bf10fad08aa427ff461bc3a7a7377f5b90ab7a810a1004d81e550a2a328bbb65. Selected baseline returned PASS_RESEARCH_BASELINE_BINDINGS_ONLY.", + "author": "unpaidmoon", + "ai_assistance": "Prepared by unpaidmoon with substantial Cursor Grok assistance. The assistant derived the projector basis and checked the polynomial. This is not an independent human or proof-kernel review.", + "open_questions": [ + { + "question": "What is the geometric multiplicity of eigenvalue 0 for A_1 on this 128-dimensional component, and what are the rational ranks of B-8I and B+8I?", + "why_it_matters": "The characteristic polynomial gives algebraic multiplicities 12, 8 and 8 for eigenvalues 0, 2 and -2. Those ranks decide whether the repeated factors are semisimple on this component.", + "first_step": "Compute the exact rank over Q of the integer matrix X of B and of X-8I and X+8I.", + "blocked_by": "A modular rank is only a lower bound on the rational rank, and this package does not run that exact elimination." + }, + { + "question": "What is an explicit integer matrix of the multiplicity-space operator T for which B acts as I_2 tensor T?", + "why_it_matters": "The degree-64 polynomial q is the characteristic polynomial of T. A matrix would determine the rational canonical form on a multiplicity basis, which the characteristic polynomial alone does not list as linear maps.", + "first_step": "Choose a G-equivariant splitting of the 128-dimensional image into 64 copies of the rational representation with minimal polynomial x^2-x+1, then read the action of B on the multiplicity coordinate.", + "blocked_by": "The Hermite basis diagonalizes the projector, not the dihedral action. This package does not construct that splitting." + } + ], + "environment": { + "python": "3.12.3", + "platform": "linux" + }, + "overlap": { + "checked_at": "2026-11-10T10:17:43Z", + "canonical_main": "17469fff78c22941280f04d39dc5a1faee09c840", + "open_pr_records": 339, + "discussion_records": 7013, + "remote_pagination_complete": true, + "local_records": 10180, + "head_records": 10178, + "local_coverage_complete": true, + "nearest_open_pr": "https://github.com/CrocSwap/integer-mult-bounds/pull/7023", + "nearest_discussion": "https://github.com/CrocSwap/integer-mult-bounds/discussions/1166#discussioncomment-19514162", + "limits": "Closed PRs were searched through local READMEs, TASK metadata and the ledger. The cached open-PR and discussion snapshot was fetched at 2026-11-10T10:16:36Z and 2026-11-10T10:17:43Z. A later open-PR search for the phrase 128-dimensional isotypic returned only PR #7023." + } +} diff --git a/research/a1-isotypic-128-charpoly/VERIFY.json b/research/a1-isotypic-128-charpoly/VERIFY.json new file mode 100644 index 0000000000..516e2756b6 --- /dev/null +++ b/research/a1-isotypic-128-charpoly/VERIFY.json @@ -0,0 +1,34 @@ +{ + "schema_version": 1, + "status": "ready", + "python": "3.12", + "requirements": [], + "inputs": [ + "research/a1-joint-216-charpoly/evidence/basis.json", + "research/a1-joint-216-charpoly/evidence/matrix.json", + "research/centre-polynomial-cache-limits/frozen/work/centre-encoding-campaign-20261014/lanes/circuit-interface/INTERFACE.json" + ], + "commands": [ + { + "argv": [ + "{python}", + "-B", + "-W", + "error", + "{package}/verify.py" + ], + "timeout_seconds": 120, + "checks": "Rebuild the pinned query, the 216-dimensional joint basis, the dihedral action, the rank-128 projector and the characteristic polynomial of A_1; check the square, the factorisation, the determinant samples and the rejection controls." + }, + { + "argv": [ + "{python}", + "-B", + "-O", + "{package}/verify.py" + ], + "timeout_seconds": 120, + "checks": "Repeat the same certificate with assertions removed. Acceptance uses explicit checks, so optimization does not skip a theorem or a rejection control." + } + ] +} diff --git a/research/a1-isotypic-128-charpoly/verify.py b/research/a1-isotypic-128-charpoly/verify.py new file mode 100644 index 0000000000..314388d998 --- /dev/null +++ b/research/a1-isotypic-128-charpoly/verify.py @@ -0,0 +1,1390 @@ +#!/usr/bin/env python3 +"""Characteristic polynomial of A_1 on one 128-dimensional dihedral isotypic component. + +The pinned query is the historical 8-by-112 centre matrix. E is the joint +(+1,+1) eigenspace of the c=1 sign stabilizer and the involutions g and t. +Inside E, G is the order-12 dihedral image generated by the signed +permutations r and u. The character on the rotations r^k, k = 0..5, is +2, 1, -1, -2, -1, 1, and the character vanishes on the reflections. Its +isotypic component W has dimension 128. This verifier computes the monic +characteristic polynomial of A_1 on W. + +The computation is exact linear algebra over Q: one characteristic +polynomial, recomputed from the pinned query. It does not change the +selected conditional kappa. Acceptance uses explicit checks, so -O does not +skip a theorem or a rejection control. +""" + +import hashlib +import json +import pathlib + +ROOT = pathlib.Path(__file__).resolve().parents[2] +INTERFACE = ( + ROOT + / "research/centre-polynomial-cache-limits/frozen/work" + / "centre-encoding-campaign-20261014/lanes/circuit-interface/INTERFACE.json" +) +BASIS_PATH = ROOT / "research/a1-joint-216-charpoly/evidence/basis.json" +MATRIX_PATH = ROOT / "research/a1-joint-216-charpoly/evidence/matrix.json" +SELECTED = ROOT / "certificates/selected-result.json" +INTERFACE_SHA256 = "5c6c51fbffc07a683e922f17532fd90262ae9e23cc3560a696002073dd234ef3" +BASIS_SHA256 = "e198cd71f4100d2629ecc25fb40b135a8c18c7756be8bc80bf24d7b349af418a" +MATRIX_SHA256 = "44f40ad5a12daf7c8493b42e1e9ddc4ef338b317ee5fdfc6de1b413642ca97aa" +SELECTED_SHA256 = "9b8d0eea95b3bf0dfcb7938a732a6439b87a2685dec2a5ab9e525b3b34debbfa" +SELECTED_KAPPA = "781451/250000000" +C1 = 1 +G_INV = (0, 2, 1, 3, 4, 6, 5, 7) +G_EPS = (1, 1, 1, 1, -1, 1, 1, -1) +T_INV = (0, 5, 6, 3, 4, 1, 2, 7) +T_EPS = (1, 1, 1, -1, 1, 1, 1, -1) +U_INV = (0, 1, 2, 3, 7, 6, 5, 4) +U_EPS = (1, 1, 1, 1, 1, 1, 1, 1) +R_INV = (0, 1, 6, 7, 3, 2, 5, 4) +R_EPS = (-1, 1, 1, 1, 1, 1, 1, -1) +# chi(r^k) for the two-dimensional rational representation. Reflections have +# character 0, so the integer combination 6P uses only the rotations. +CHI = (2, 1, -1, -2, -1, 1) +OTHER_CHI = (2, -1, -1, 2, -1, -1) +# Monic factors of chi_{A_1} on W, ascending coefficients, with multiplicity. +# These are the factors already identified on the 216-dimensional joint space, +# kept here with the multiplicities that lie on W. f and g are the degree-11 +# pair, and h and k are the degree-14 pair that appears squared in that +# factorization. g(x) = -f(-x) and k(x) = h(-x). +F11 = (-12288, 102144, -196352, 115776, 19264, -33888, 2128, 3208, -192, -108, 2, 1) +G11 = (12288, 102144, 196352, 115776, -19264, -33888, -2128, 3208, 192, -108, -2, 1) +H14 = (-1167360, -702464, 3049472, 704000, -2179840, -260736, 621184, 44928, -76384, -3632, 4376, 136, -112, -2, 1) +K14 = (-1167360, 702464, 3049472, -704000, -2179840, 260736, 621184, -44928, -76384, 3632, 4376, -136, -112, 2, 1) +COMPONENT_FACTORS = ( + (12, (0, 1)), + (8, (-2, 1)), + (8, (2, 1)), + (2, F11), + (2, G11), + (2, H14), + (2, K14), +) +SQUARE_FACTORS = ( + (6, (0, 1)), + (4, (-2, 1)), + (4, (2, 1)), + (1, F11), + (1, G11), + (1, H14), + (1, K14), +) +# Published multiplicities on the whole 216-dimensional joint space. The +# comparison uses these coefficient lists; it does not recompute that polynomial. +PUBLISHED_216 = ( + (24, (0, 1)), + (12, (-2, 1)), + (12, (2, 1)), + (2, (-6, 1)), + (2, (6, 1)), + (1, (-4, -2, 1)), + (1, (-4, 2, 1)), + (2, F11), + (2, G11), + (1, (-498401280, -739049472, 132964352, 453709824, 8403968, -86923264, -97024, 6293376, -173568, -202144, 10112, 2824, -180, -14, 1)), + (1, (-498401280, 739049472, 132964352, -453709824, 8403968, 86923264, -97024, -6293376, -173568, 202144, 10112, -2824, -180, 14, 1)), + (2, H14), + (2, K14), + (1, (-1702232064, -112656384, 3560341504, 1656389632, -1259880448, -773404672, 108563456, 105269248, -321280, -6136576, -277632, 173824, 11488, -2344, -180, 12, 1)), + (1, (-1702232064, 112656384, 3560341504, -1656389632, -1259880448, 773404672, 108563456, -105269248, -321280, 6136576, -277632, -173824, 11488, 2344, -180, -12, 1)), +) +# q(x) = x^6 r(x^2), ascending coefficients of the degree-29 polynomial r. +R_COEFFICIENTS = ( + -52675902681576269414400, + 2303306134076870941999104, + -18827961748477891024257024, + 73124429304073607471693824, + -171456969792808134744473600, + 271456076002980399446753280, + -309537536542443194852311040, + 264951255789256454094979072, + -175106550821218960407527424, + 91149737035716567039475712, + -37905952006491036612296704, + 12722870692816032484032512, + -3471386398306750712250368, + 773676955090062483128320, + -141267596592764868886528, + 21162534505722224312320, + -2601171087886668267520, + 262024005207787307008, + -21579308813538295808, + 1447508150155345920, + -78661706621779968, + 3437815337254912, + -119644063870976, + 3271915571200, + -69036760064, + 1095226368, + -12571072, + 98112, + -464, + 1, +) +PIVOT_ROWS = ( + 2, 3, 4, 5, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, + 24, 25, 28, 29, 30, 31, 40, 41, 44, 45, 46, 47, 48, 49, 50, 51, 56, 57, 60, 61, + 78, 79, 80, 81, 84, 85, 86, 87, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, + 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, + 126, 127, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, + 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, + 178, 179, 180, 181, 182, 183, 184, 185, 186, 187, 188, 189, 190, 191, +) +DIAGONALS = ( + 2, 2, 2, 2, 6, 6, 2, 2, 2, 2, 2, 2, 6, 6, 2, 2, + 6, 6, 6, 6, 2, 2, 6, 6, 6, 6, 2, 2, 2, 2, 2, 2, + 6, 6, 6, 6, 6, 6, 6, 6, 2, 2, 6, 6, 2, 2, 6, 6, + 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 6, 6, 6, 6, + 6, 6, 6, 6, 2, 2, 2, 2, 6, 6, 6, 6, 6, 6, 6, 6, + 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, + 2, 2, 2, 2, 2, 2, 2, 2, 6, 6, 6, 6, 6, 6, 6, 6, + 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, +) +CHI_A_SHA256 = "bf10fad08aa427ff461bc3a7a7377f5b90ab7a810a1004d81e550a2a328bbb65" +CHI_B_SHA256 = "33f15f1bcc5d9a57fec18f40653e2362792a77b2c82bdbcf4bbd469694dd03a2" +CHI_SQRT_SHA256 = "d895a734501995193ce1a820dfe8d7cf9121fdf12a1abd1620a344e3315f5762" + +# The c=1 symmetry, lattice, port action and modular characteristic-polynomial +# routines below are the ones in research/a1-joint-216-charpoly/verify.py. +# This package adds the dihedral projector and the restriction to W. + +def fail(message): + raise SystemExit("FAIL: " + message) + + +def expect(condition, message): + if not condition: + fail(message) + + +def sha256_bytes(data): + return hashlib.sha256(data).hexdigest() + + +def sha256_file(path): + return sha256_bytes(path.read_bytes()) + + +def canonical(obj): + return json.dumps(obj, separators=(",", ":"), ensure_ascii=True) + + +def load_query(): + raw = INTERFACE.read_bytes() + expect(sha256_bytes(raw) == INTERFACE_SHA256, "INTERFACE.json SHA256 does not match the pinned query") + data = json.loads(raw) + rows = data["R"] + ports = [int(q) for q in data["ports"]] + expect(len(rows) == 8 and len(ports) == 112 and len(set(ports)) == 112, "query shape") + columns = [] + for j in range(112): + columns.append(tuple(int(rows[i].get(str(j), "0")) for i in range(8))) + return columns, ports + + +def norm(vector): + for entry in vector: + if entry: + return tuple(entry * item for item in vector) + return vector + + +def rank_solve(vectors, target=None): + width = len(vectors) + matrix = [[(vector >> row) & 1 for vector in vectors] for row in range(9)] + image = None if target is None else [(target >> row) & 1 for row in range(9)] + pivot_at = [-1] * 9 + rank = 0 + for col in range(width): + pivot = None + for row in range(rank, 9): + if matrix[row][col]: + pivot = row + break + if pivot is None: + continue + matrix[rank], matrix[pivot] = matrix[pivot], matrix[rank] + if image is not None: + image[rank], image[pivot] = image[pivot], image[rank] + for row in range(9): + if row != rank and matrix[row][col]: + for j in range(width): + matrix[row][j] ^= matrix[rank][j] + if image is not None: + image[row] ^= image[rank] + pivot_at[rank] = col + rank += 1 + if rank == 9: + break + if target is None: + return rank + if any(image[row] for row in range(rank, 9)): + return None + coef = [0] * width + for row in range(rank): + coef[pivot_at[row]] = image[row] + return coef + + +def prepare(columns, ports): + port_of = {columns[j]: ports[j] for j in range(112)} + lookup = {norm(column): port_of[column] for column in columns} + basis = [] + for port in ports: + if rank_solve(basis + [port]) == len(basis) + 1: + basis.append(port) + if len(basis) == 9: + break + expect(len(basis) == 9, "the 112 ports do not span F2^9") + coords = {address: tuple(rank_solve(basis, address)) for address in range(512)} + gram = [[0] * 8 for _ in range(8)] + for column in columns: + for i in range(8): + for k in range(8): + gram[i][k] += column[i] * column[k] + expect( + all(gram[i][k] == (56 if i == k else 0) for i in range(8) for k in range(8)), + "RR^T is not 56 I", + ) + expect(all(bin(port).count("1") % 2 == 1 for port in ports), "a port label has even parity") + expect(all(entry in (-1, 0, 1) for column in columns for entry in column), "a root entry is outside {-1,0,1}") + expect(all(sum(entry != 0 for entry in column) == 4 for column in columns), "a root does not have weight 4") + return port_of, lookup, basis, coords + + +def apply_linear(matrix, address): + image = 0 + for bit in range(9): + if (address >> bit) & 1: + image ^= matrix[bit] + return image + + +def fixed_addresses(matrix): + rows = [ + [((matrix[col] >> row) & 1) ^ (1 if row == col else 0) for col in range(9)] + for row in range(9) + ] + rank = 0 + for col in range(9): + pivot = None + for row in range(rank, 9): + if rows[row][col]: + pivot = row + break + if pivot is None: + continue + rows[rank], rows[pivot] = rows[pivot], rows[rank] + for row in range(9): + if row != rank and rows[row][col]: + rows[row] = [left ^ right for left, right in zip(rows[row], rows[rank])] + rank += 1 + return 1 << (9 - rank) + + +def trace_s(inv, eps): + return sum(eps[index] for index in range(8) if inv[index] == index) + + +def operator(inv, eps, port_of, lookup, basis, coords): + mapping = {} + for column, port in port_of.items(): + image = tuple(eps[index] * column[inv[index]] for index in range(8)) + matched = lookup.get(norm(image)) + if matched is None: + return None + mapping[port] = matched + matrix = [] + for bit in range(9): + image = 0 + for index, flag in enumerate(coords[1 << bit]): + if flag: + image ^= mapping[basis[index]] + matrix.append(image) + for port, image in mapping.items(): + if apply_linear(matrix, port) != image: + return None + images = [mapping[vector] for vector in basis] + for left in range(9): + for right in range(left, 9): + if (bin(basis[left] & basis[right]).count("1") & 1) != ( + bin(images[left] & images[right]).count("1") & 1 + ): + return None + return {"inv": tuple(inv), "eps": tuple(eps), "T": tuple(matrix)} + + +def compose(left, right): + inv = tuple(right["inv"][left["inv"][index]] for index in range(8)) + eps = tuple(left["eps"][index] * right["eps"][left["inv"][index]] for index in range(8)) + matrix = tuple(apply_linear(left["T"], column) for column in right["T"]) + return {"inv": inv, "eps": eps, "T": matrix} + + +def same_operator(left, right): + return left["inv"] == right["inv"] and left["eps"] == right["eps"] and left["T"] == right["T"] + + +def sign_group(port_of, lookup, basis, coords, sector): + identity = tuple(range(8)) + group = [] + for mask in range(128): + eps = tuple([1] + [-1 if (mask >> bit) & 1 else 1 for bit in range(7)]) + built = operator(identity, eps, port_of, lookup, basis, coords) + if built is None: + continue + if apply_linear(built["T"], sector) == sector: + group.append(built) + return group + + +def action_tables(element): + inv, eps, transform = element["inv"], element["eps"], element["T"] + sigma = [0] * 8 + for out_index, in_index in enumerate(inv): + sigma[in_index] = out_index + dest = [0] * 4096 + sign = [0] * 4096 + for coordinate in range(8): + image_coordinate = sigma[coordinate] + image_sign = eps[image_coordinate] + for address in range(512): + src = coordinate * 512 + address + dest[src] = image_coordinate * 512 + apply_linear(transform, address) + sign[src] = image_sign + return dest, sign + + +def action_trace(dest, sign): + return sum(sign[index] for index in range(4096) if dest[index] == index) + + +def check_involution(dest, sign, name): + for src in range(4096): + mid = dest[src] + if dest[mid] != src or sign[src] * sign[mid] != 1: + fail(name + " is not a signed involution") + + +class Lattice: + """Solution space of v[dst] = sign * v[src] equations on Q^4096.""" + + def __init__(self): + self.parent = list(range(4096)) + self.rel = [1] * 4096 + self.dead = [False] * 4096 + + def find(self, index): + sign = 1 + seen = [] + while self.parent[index] != index: + seen.append(index) + sign *= self.rel[index] + index = self.parent[index] + accumulated = 1 + for step in reversed(seen): + accumulated = self.rel[step] * accumulated + self.rel[step] = accumulated + self.parent[step] = index + return index, sign + + def constrain(self, dst, src, sgn): + rd, sd = self.find(dst) + rs, ss = self.find(src) + if rd == rs: + if sd != sgn * ss: + self.dead[rd] = True + return + self.parent[rd] = rs + self.rel[rd] = sgn * ss * sd + if self.dead[rd] or self.dead[rs]: + self.dead[rs] = True + + def add_operator(self, dest, sign, eigenvalue): + for src in range(4096): + self.constrain(dest[src], src, sign[src] * eigenvalue) + + def basis(self): + groups = {} + for index in range(4096): + root, sign = self.find(index) + if self.dead[root]: + continue + groups.setdefault(root, []).append((index, sign)) + vectors = [] + for items in groups.values(): + items.sort() + pivot = items[0][1] + vectors.append((items[0][0], [(index, coeff * pivot) for index, coeff in items])) + vectors.sort() + return vectors + + +def fixed_by(dest, sign, items, eigenvalue): + values = {index: coeff for index, coeff in items} + for index, coeff in items: + image = dest[index] + if values.get(image, 0) != sign[index] * eigenvalue * coeff: + return False + for index, coeff in items: + # Every preimage that the operator sends onto this support must agree. + # The support check above is complete because the operator is a bijection + # and the vector vanishes off the support: if rho v differed at an + # outside index, some support source would map there. + if coeff == 0: + return False + return True + + +def column_table(columns, ports): + raw = [{} for _ in range(4096)] + for port, root in zip(ports, columns): + support = [k for k in range(8) if root[k]] + if port & 1: + for coordinate in support: + weight = root[coordinate] + for address in range(512): + src = coordinate * 512 + address + target = address ^ port + bucket = raw[src] + for out_coordinate in support: + index = out_coordinate * 512 + target + bucket[index] = bucket.get(index, 0) + root[out_coordinate] * weight + else: + for coordinate in support: + weight = root[coordinate] + for address in range(512): + src = coordinate * 512 + address + parity = -1 if (port & address).bit_count() & 1 else 1 + term = weight * parity + bucket = raw[src] + for out_coordinate in support: + index = out_coordinate * 512 + address + bucket[index] = bucket.get(index, 0) + root[out_coordinate] * term + return [{index: value for index, value in column.items() if value} for column in raw] + + +def literal_image(items, columns, ports): + image = {} + for port, root in zip(ports, columns): + support = [k for k in range(8) if root[k]] + buckets = {} + for index, coeff in items: + coordinate = index // 512 + address = index - 512 * coordinate + if root[coordinate]: + buckets[address] = buckets.get(address, 0) + root[coordinate] * coeff + if port & 1: + for address, acc in buckets.items(): + target = address ^ port + for coordinate in support: + index = coordinate * 512 + target + image[index] = image.get(index, 0) + root[coordinate] * acc + else: + for address, acc in buckets.items(): + if (port & address).bit_count() & 1: + acc = -acc + for coordinate in support: + index = coordinate * 512 + address + image[index] = image.get(index, 0) + root[coordinate] * acc + return {index: value for index, value in image.items() if value} + + +def column_image(columns_b, items): + image = {} + for index, coeff in items: + for target, value in columns_b[index].items(): + image[target] = image.get(target, 0) + coeff * value + return {index: value for index, value in image.items() if value} + + +def operators_commute(columns_b, dest, sign): + for src in range(4096): + left_sign = sign[src] + left = columns_b[dest[src]] + right = {} + for index, value in columns_b[src].items(): + image = dest[index] + right[image] = right.get(image, 0) + sign[index] * value + keys = set(left) | set(right) + for key in keys: + if left.get(key, 0) * left_sign != right.get(key, 0): + return False + return True + + +def restriction_matrix(columns_b, columns, ports, vectors): + roots = [root for root, _items in vectors] + position = {root: index for index, root in enumerate(roots)} + members = [] + for _root, items in vectors: + members.append(items) + owners = {} + grouped = [[] for _ in vectors] + for component, items in enumerate(members): + for index, coeff in items: + owners[index] = (component, coeff) + grouped[component].append((index, coeff)) + matrix = [[0] * len(vectors) for _ in vectors] + sizes = [len(items) for _root, items in vectors] + for column, (_root, items) in enumerate(vectors): + by_columns = column_image(columns_b, items) + by_literal = literal_image(items, columns, ports) + expect(by_columns == by_literal, "the two actions of B disagree on a basis vector") + alphas = {} + for index, value in by_columns.items(): + info = owners.get(index) + expect(info is not None, "B sends a basis vector outside the joint space") + component, coeff = info + alpha = value * coeff + previous = alphas.get(component) + expect(previous is None or previous == alpha, "B image is not constant on a component") + alphas[component] = alpha + for component, alpha in alphas.items(): + for index, coeff in grouped[component]: + expect(by_columns.get(index, 0) == coeff * alpha, "B image misses a component coordinate") + matrix[component][column] = alpha + for index, value in by_columns.items(): + component, _coeff = owners[index] + expect(component in alphas, "an image coordinate has no component value") + for left in range(len(vectors)): + for right in range(left + 1, len(vectors)): + expect( + sizes[left] * matrix[left][right] == sizes[right] * matrix[right][left], + "the Gram identity size_r M_rs = size_s M_sr failed", + ) + return matrix, sizes + + +def character_joint_dimension(group, actors): + """Trace-formula dimension of the trivial character's (+1,+1) space.""" + identity = {"inv": tuple(range(8)), "eps": (1, 1, 1, 1, 1, 1, 1, 1), "T": tuple(1 << bit for bit in range(9))} + + def character_trace(element): + total = 0 + for symmetry in group: + composed = compose(element, symmetry) + dest, sign = action_tables(composed) + closed = trace_s(composed["inv"], composed["eps"]) * fixed_addresses(composed["T"]) + expect(action_trace(dest, sign) == closed, "permutation trace disagrees with the closed trace formula") + total += closed + expect(total % len(group) == 0, "trivial character trace is not integral") + return total // len(group) + + traces = [character_trace(element) for element in (identity, actors[0], actors[1], compose(actors[0], actors[1]))] + numerator = traces[0] + traces[1] + traces[2] + traces[3] + expect(numerator % 4 == 0, "joint trace numerator is not divisible by 4") + return traces[0], numerator // 4 + + +def coeff_bound(size, row_sum): + best = 1 + binomial = 1 + power = 1 + for degree in range(size + 1): + term = binomial * power + if term > best: + best = term + if degree < size: + binomial = binomial * (size - degree) // (degree + 1) + power *= row_sum + return best + + +def is_prime(number): + if number < 2: + return False + if number % 2 == 0: + return number == 2 + divisor = 3 + while divisor * divisor <= number: + if number % divisor == 0: + return False + divisor += 2 + return True + + +def primes_past(limit): + modulus = 1 + primes = [] + candidate = 1_000_003 + while modulus <= limit: + if is_prime(candidate): + primes.append(candidate) + modulus *= candidate + candidate += 2 + return primes, modulus + + +def charpoly_mod(matrix, prime): + size = len(matrix) + hess = [[entry % prime for entry in row] for row in matrix] + for column in range(size - 2): + pivot = None + for row in range(column + 1, size): + if hess[row][column]: + pivot = row + break + if pivot is None: + continue + if pivot != column + 1: + hess[column + 1], hess[pivot] = hess[pivot], hess[column + 1] + for row in hess: + row[column + 1], row[pivot] = row[pivot], row[column + 1] + inverse = pow(hess[column + 1][column], -1, prime) + for row in range(column + 2, size): + factor = hess[row][column] + if not factor: + continue + factor = factor * inverse % prime + target = hess[row] + source = hess[column + 1] + for index in range(column, size): + target[index] = (target[index] - factor * source[index]) % prime + for host in hess: + host[column + 1] = (host[column + 1] + factor * host[row]) % prime + minors = [[1], [(prime - hess[0][0]) % prime, 1]] + for end in range(1, size): + previous = minors[end] + diagonal = hess[end][end] + output = [0] * (end + 2) + for index, coefficient in enumerate(previous): + output[index + 1] = coefficient + output[index] = (output[index] - diagonal * coefficient) % prime + product = 1 + for start in range(end - 1, -1, -1): + product = product * hess[start + 1][start] % prime + if product == 0: + break + scale = hess[start][end] * product % prime + if scale == 0: + continue + for index, coefficient in enumerate(minors[start]): + output[index] = (output[index] - scale * coefficient) % prime + expect(output[-1] % prime == 1, "Hessenberg characteristic polynomial is not monic") + minors.append(output) + return minors[size] + + +def det_mod(matrix, prime): + size = len(matrix) + work = [[entry % prime for entry in row] for row in matrix] + determinant = 1 + for column in range(size): + pivot = None + for row in range(column, size): + if work[row][column]: + pivot = row + break + if pivot is None: + return 0 + if pivot != column: + work[column], work[pivot] = work[pivot], work[column] + determinant = (-determinant) % prime + determinant = determinant * work[column][column] % prime + inverse = pow(work[column][column], -1, prime) + for row in range(column + 1, size): + factor = work[row][column] + if not factor: + continue + factor = factor * inverse % prime + target = work[row] + source = work[column] + for index in range(column, size): + target[index] = (target[index] - factor * source[index]) % prime + return determinant + + +def crt_signed(residues, primes): + value = 0 + modulus = 1 + for residue, prime in zip(residues, primes): + increment = ((residue - value) % prime) * pow(modulus % prime, -1, prime) % prime + value += increment * modulus + modulus *= prime + if value > modulus // 2: + value -= modulus + return value + + +def poly_mul(left, right): + output = [0] * (len(left) + len(right) - 1) + for i, a in enumerate(left): + for j, b in enumerate(right): + output[i + j] += a * b + return output + + +def poly_pow(base, exponent): + result = [1] + while exponent: + if exponent & 1: + result = poly_mul(result, base) + exponent >>= 1 + if exponent: + base = poly_mul(base, base) + return result + + +def product_factors(factors): + value = [1] + for multiplicity, coefficients in factors: + value = poly_mul(value, poly_pow(list(coefficients), multiplicity)) + return value + + +def scale_to_b(coefficients): + degree = len(coefficients) - 1 + return [coefficient * 4 ** (degree - index) for index, coefficient in enumerate(coefficients)] + + +def monic_divides(divisor, target): + if len(divisor) == 1: + return divisor == [1] and target == [0] or divisor == [1] + work = list(target) + degree = len(divisor) - 1 + while len(work) - 1 >= degree and not (len(work) == 1 and work[0] == 0): + shift = len(work) - 1 - degree + lead = work[-1] + for index, coefficient in enumerate(divisor): + work[shift + index] -= lead * coefficient + while len(work) > 1 and work[-1] == 0: + work.pop() + return work == [0] or work == [] + + +def trim(coefficients, prime): + output = [coefficient % prime for coefficient in coefficients] + while len(output) > 1 and output[-1] == 0: + output.pop() + return output + + +def monic_mod(coefficients, prime): + inverse = pow(coefficients[-1], -1, prime) + return [(coefficient * inverse) % prime for coefficient in coefficients] + + +def remainder_mod(numerator, divisor, prime): + work = trim(numerator, prime) + modulus = monic_mod(trim(divisor, prime), prime) + degree = len(modulus) - 1 + while len(work) - 1 >= degree and not (len(work) == 1 and work[0] == 0): + shift = len(work) - 1 - degree + lead = work[-1] + for index, coefficient in enumerate(modulus): + work[shift + index] = (work[shift + index] - lead * coefficient) % prime + while len(work) > 1 and work[-1] == 0: + work.pop() + return work + + +def gcd_mod(left, right, prime): + first = trim(left, prime) + second = trim(right, prime) + while not (len(second) == 1 and second[0] == 0): + first, second = second, remainder_mod(first, second, prime) + if len(first) == 1 and first[0] == 0: + return [0] + return monic_mod(first, prime) + + +def pow_mod_poly(exponent, modulus, prime): + result = [1] + base = [0, 1] + degree = len(modulus) - 1 + while exponent: + if exponent & 1: + result = remainder_mod(poly_mul(result, base), modulus, prime) + exponent >>= 1 + if exponent: + base = remainder_mod(poly_mul(base, base), modulus, prime) + expect(len(result) <= degree + 1, "polynomial remainder grew past the modulus") + return result + + +def is_square(number): + if number < 0: + return False + root = int(number ** 0.5) + while (root + 1) * (root + 1) <= number: + root += 1 + while root * root > number: + root -= 1 + return root * root == number + + +def x_minus_power(exponent, modulus, prime): + power = pow_mod_poly(exponent, modulus, prime) + if len(power) < 2: + power = power + [0] * (2 - len(power)) + power = power[:] + power[1] = (power[1] - 1) % prime + return trim(power, prime) + + +def irreducible_mod(coefficients, prime): + reduced = trim(coefficients, prime) + if len(reduced) <= 2: + return len(reduced) == 2 + if reduced[-1] % prime == 0: + return False + modulus = monic_mod(reduced, prime) + degree = len(modulus) - 1 + derivative = [((index * coefficient) % prime) for index, coefficient in enumerate(modulus)][1:] + if gcd_mod(modulus, derivative or [0], prime) != [1]: + return False + primes = [] + value = degree + divisor = 2 + while divisor * divisor <= value: + if value % divisor == 0: + primes.append(divisor) + while value % divisor == 0: + value //= divisor + divisor += 1 + if value > 1: + primes.append(value) + if x_minus_power(prime ** degree, modulus, prime) != [0]: + return False + for factor in primes: + if gcd_mod(x_minus_power(prime ** (degree // factor), modulus, prime), modulus, prime) != [1]: + return False + return True + + +def prove_irreducible(coefficients): + degree = len(coefficients) - 1 + expect(coefficients[-1] == 1, "a factor is not monic") + if degree == 1: + return + if degree == 2: + discriminant = coefficients[1] * coefficients[1] - 4 * coefficients[0] + expect(not is_square(discriminant), "a quadratic factor is reducible over Z") + return + candidate = 1009 + tried = 0 + while tried < 80: + if is_prime(candidate) and irreducible_mod(coefficients, candidate): + return + if is_prime(candidate): + tried += 1 + candidate += 2 + fail("no prime in the first 80 primes at or above 1009 proved a factor irreducible") + + +def descending(coefficients): + terms = [] + for power, coefficient in enumerate(coefficients): + if coefficient: + terms.append((power, coefficient)) + parts = [] + for power, coefficient in reversed(terms): + if power == 0: + body = str(abs(coefficient)) + if not parts: + parts.append(str(coefficient)) + else: + parts.append(("+ " if coefficient > 0 else "- ") + body) + continue + magnitude = abs(coefficient) + letter = "x" if power == 1 else "x^" + str(power) + if magnitude == 1: + body = letter + else: + body = str(magnitude) + letter + if not parts: + parts.append(body if coefficient > 0 else "-" + body) + else: + parts.append(("+ " if coefficient > 0 else "- ") + body) + return " ".join(parts) + + +def row_sum(matrix): + return max(sum(abs(entry) for entry in row) for row in matrix) + + +def evaluate(coefficients, point, prime): + value = 0 + power = 1 + for coefficient in coefficients: + value = (value + coefficient * power) % prime + power = power * point % prime + return value + + +def shifted_det(matrix, point, prime): + size = len(matrix) + work = [[(-entry) % prime for entry in row] for row in matrix] + for index in range(size): + work[index][index] = (work[index][index] + point) % prime + return det_mod(work, prime) + + +def power(element, count, identity): + value = identity + for _ in range(count): + value = compose(value, element) + return value + + +def egcd(left, right): + if right == 0: + return abs(left), 1 if left >= 0 else -1, 0 + gcd_value, bezout_x, bezout_y = egcd(right, left % right) + return gcd_value, bezout_y, bezout_x - (left // right) * bezout_y + + +def coordinates_216(sparse, owners, vectors): + alphas = {} + for index, value in sparse.items(): + if not value: + continue + info = owners.get(index) + expect(info is not None, "a vector leaves the 216-dimensional joint space") + component, base = info + alpha = value * base + previous = alphas.get(component) + expect(previous is None or previous == alpha, "a joint component receives two coefficients") + alphas[component] = alpha + for component, alpha in alphas.items(): + for index, coeff in vectors[component][1]: + expect(sparse.get(index, 0) == coeff * alpha, "a joint component is only partly present") + return alphas + + +def permute_component(dest, sign, items, owners, vectors): + image = {} + for index, coeff in items: + target = dest[index] + image[target] = image.get(target, 0) + sign[index] * coeff + image = {index: value for index, value in image.items() if value} + return coordinates_216(image, owners, vectors) + + +def monomial(element, vectors, owners): + dest, sign = action_tables(element) + images = [] + for _root, items in vectors: + coords = permute_component(dest, sign, items, owners, vectors) + expect(len(coords) == 1, "a rotation does not permute the joint basis by a signed coordinate") + (target, coef), = coords.items() + expect(coef in (-1, 1), "a joint-basis coefficient is not a unit") + images.append((target, coef)) + expect(len({target for target, _coef in images}) == 216, "a rotation is not bijective on the joint basis") + return images + + +def integer_projector(monos, coefficients): + columns = [] + trace = 0 + for src in range(216): + acc = {} + for power_index, coef in enumerate(coefficients): + if coef == 0: + continue + target, sign = monos[power_index][src] + acc[target] = acc.get(target, 0) + coef * sign + acc = {key: value for key, value in acc.items() if value} + columns.append(acc) + trace += acc.get(src, 0) + return columns, trace + + +def apply_columns(columns, sparse): + acc = {} + for index, coef in sparse.items(): + for target, value in columns[index].items(): + acc[target] = acc.get(target, 0) + coef * value + return {key: value for key, value in acc.items() if value} + + +def scale_sparse(sparse, factor): + return {key: factor * value for key, value in sparse.items() if factor * value} + + +def gcd_many(values): + value = 0 + for item in values: + item = abs(item) + while item: + value, item = item, value % item + return value + + +def lattice_rank(columns): + pivots = [] + for raw in columns: + vector = {key: value for key, value in raw.items() if value} + for pivot_index, pivot in pivots: + coef = vector.get(pivot_index, 0) + if not coef: + continue + scale = pivot[pivot_index] + updated = {} + for key in set(vector) | set(pivot): + value = scale * vector.get(key, 0) - coef * pivot.get(key, 0) + if value: + updated[key] = value + vector = updated + if vector: + common = gcd_many(vector.values()) + if common > 1: + vector = {key: value // common for key, value in vector.items()} + if not vector: + continue + pivot_index = min(vector) + if vector[pivot_index] < 0: + vector = {key: -value for key, value in vector.items()} + pivots.append((pivot_index, vector)) + return len(pivots) + + +def hermite_basis(columns): + """Z-basis of the column lattice of 6P. + + Column operations are unimodular. A finished pivot column is never mixed + by a later gcd step, so it stays zero above its pivot row. Entries below + that row may be nonzero, and the pivot submatrix is lower triangular. + """ + width = len(columns) + height = 216 + expect(width == 216, "the projector is not presented on the 216-dimensional basis") + dense = [[columns[col].get(row, 0) for col in range(width)] for row in range(height)] + used = [False] * width + chosen = [] + pivot_rows = [] + for row in range(height): + pivot = None + for col in range(width): + if not used[col] and dense[row][col]: + pivot = col + break + if pivot is None: + continue + for col in range(width): + if used[col] or col == pivot or not dense[row][col]: + continue + gcd_value, s_coef, t_coef = egcd(dense[row][pivot], dense[row][col]) + diagonal = dense[row][pivot] + other = dense[row][col] + expect(s_coef * diagonal + t_coef * other == gcd_value, "extended gcd failed") + expect( + s_coef * (diagonal // gcd_value) - t_coef * (-(other // gcd_value)) == 1, + "a Hermite column operation is not unimodular", + ) + for host in dense: + left = host[pivot] + right = host[col] + host[pivot] = s_coef * left + t_coef * right + host[col] = (diagonal // gcd_value) * right - (other // gcd_value) * left + if dense[row][pivot] < 0: + for host in dense: + host[pivot] = -host[pivot] + diagonal = dense[row][pivot] + expect(diagonal > 0, "a Hermite pivot is not positive") + for col in range(width): + if col == pivot: + continue + quotient = dense[row][col] // diagonal + if quotient: + for host in dense: + host[col] -= quotient * host[pivot] + expect( + all(dense[row][col] in range(diagonal) for col in range(width) if col != pivot), + "a Hermite row was not reduced", + ) + expect( + all(dense[earlier][pivot] == 0 for earlier in range(row)), + "a pivot column gained an entry above its pivot", + ) + used[pivot] = True + chosen.append(pivot) + pivot_rows.append(row) + for col in range(width): + if used[col]: + continue + expect(all(dense[row][col] == 0 for row in range(height)), "a nonpivot Hermite column survived") + basis = [] + for col in chosen: + basis.append({index: dense[index][col] for index in range(height) if dense[index][col]}) + return basis, pivot_rows + + +def express_in_hermite(basis_columns, pivot_rows, images): + size = len(basis_columns) + matrix = [[0] * len(images) for _ in range(size)] + for col, image in enumerate(images): + residual = {index: value for index, value in enumerate(image) if value} + for index in range(size): + pivot = pivot_rows[index] + diagonal = basis_columns[index][pivot] + value = residual.get(pivot, 0) + expect(value % diagonal == 0, "B does not preserve the column lattice of 6P") + coef = value // diagonal + matrix[index][col] = coef + if coef: + for key, entry in basis_columns[index].items(): + residual[key] = residual.get(key, 0) - coef * entry + if residual[key] == 0: + del residual[key] + expect(not residual, "an image vector leaves the Hermite span") + return matrix + + +def identity_element(): + return { + "inv": tuple(range(8)), + "eps": (1, 1, 1, 1, 1, 1, 1, 1), + "T": tuple(1 << bit for bit in range(9)), + } + + +def poly_sha(polynomial): + return sha256_bytes(canonical(polynomial).encode("ascii")) + + +def integer_square_root(polynomial): + degree = len(polynomial) - 1 + expect(degree % 2 == 0 and polynomial[-1] == 1, "chi is not monic of even degree") + half = degree // 2 + leading = [0] * (half + 1) + leading[0] = 1 + for step in range(1, half + 1): + target = polynomial[degree - step] + mixed = sum(leading[index] * leading[step - index] for index in range(1, step)) + expect((target - mixed) % 2 == 0, "chi is not a square in Z[x]") + leading[step] = (target - mixed) // 2 + root = [leading[half - index] for index in range(half + 1)] + expect(poly_mul(root, root) == polynomial, "the square root does not square back to chi") + return root + + +def charpoly_integer(matrix): + size = len(matrix) + bound = coeff_bound(size, row_sum(matrix)) + primes, _modulus = primes_past(2 * bound) + residues = [charpoly_mod(matrix, prime) for prime in primes] + polynomial = [ + crt_signed([residue[degree] for residue in residues], primes) for degree in range(size + 1) + ] + expect(polynomial[-1] == 1, "the characteristic polynomial is not monic") + for point in (0, 1, 2, 3, 5): + for prime in primes[:3]: + expect( + evaluate(polynomial, point, prime) == shifted_det(matrix, point, prime), + "a Hessenberg characteristic polynomial disagrees with the determinant", + ) + return polynomial + + +def negated_argument(coefficients): + return tuple(coef if index % 2 == 0 else -coef for index, coef in enumerate(coefficients)) + + +def main(): + columns, ports = load_query() + selected_raw = SELECTED.read_bytes() + expect(sha256_bytes(selected_raw) == SELECTED_SHA256, "selected-result.json bytes changed") + selected = json.loads(selected_raw) + expect(selected["kappa"] == SELECTED_KAPPA, "selected kappa is not 781451/250000000") + basis_raw = BASIS_PATH.read_bytes() + expect(sha256_bytes(basis_raw) == BASIS_SHA256, "the published 216-basis bytes changed") + matrix_raw = MATRIX_PATH.read_bytes() + expect(sha256_bytes(matrix_raw) == MATRIX_SHA256, "the published 216-matrix bytes changed") + published_basis = json.loads(basis_raw) + published_matrix = json.loads(matrix_raw)["rows"] + + port_of, lookup, port_basis, coords = prepare(columns, ports) + group = sign_group(port_of, lookup, port_basis, coords, C1) + expect(len(group) == 16, "the c=1 sign stabilizer does not have order 16") + first = operator(G_INV, G_EPS, port_of, lookup, port_basis, coords) + second = operator(T_INV, T_EPS, port_of, lookup, port_basis, coords) + reflection = operator(U_INV, U_EPS, port_of, lookup, port_basis, coords) + rotation = operator(R_INV, R_EPS, port_of, lookup, port_basis, coords) + expect(first is not None and second is not None, "g or t is not a symmetry") + expect(reflection is not None and rotation is not None, "u or r is not a symmetry") + identity = identity_element() + expect(same_operator(compose(first, first), identity), "g is not an involution") + expect(same_operator(compose(second, second), identity), "t is not an involution") + expect(same_operator(compose(first, second), compose(second, first)), "g and t do not commute") + expect(apply_linear(first["T"], C1) == C1 and apply_linear(second["T"], C1) == C1, "a splitter moves sector 1") + expect(same_operator(compose(reflection, reflection), identity), "u is not an involution") + expect(same_operator(power(rotation, 6, identity), identity), "r does not have order dividing 6") + for count in (1, 2, 3): + expect(not same_operator(power(rotation, count, identity), identity), "r has order less than 6") + expect( + same_operator(compose(reflection, compose(rotation, reflection)), power(rotation, 5, identity)), + "u does not invert r", + ) + words = [] + for index in range(6): + words.append(power(rotation, index, identity)) + words.append(compose(power(rotation, index, identity), reflection)) + expect(len({(word["inv"], word["eps"], word["T"]) for word in words}) == 12, "the twelve dihedral words are not distinct") + for word in words: + expect(apply_linear(word["T"], C1) == C1, "a dihedral word moves sector 1") + trivial_dimension, joint_dimension = character_joint_dimension(group, (first, second)) + expect(trivial_dimension == 864 and joint_dimension == 216, "the trace formula does not give a 216-dimensional ++ space") + + actions = [action_tables(element) for element in group] + g_dest, g_sign = action_tables(first) + t_dest, t_sign = action_tables(second) + lattice = Lattice() + for dest, sign in actions: + lattice.add_operator(dest, sign, 1) + lattice.add_operator(g_dest, g_sign, 1) + lattice.add_operator(t_dest, t_sign, 1) + vectors = lattice.basis() + expect(len(vectors) == 216, "the explicit ++ basis does not have dimension 216") + recomputed = { + "index": "coordinate*512+address", + "eigenvalues": {"H": 1, "g": 1, "t": 1}, + "roots": [root for root, _items in vectors], + "vectors": [[[index, coeff] for index, coeff in items] for _root, items in vectors], + } + expect(recomputed == published_basis, "the recomputed basis is not the published 216-basis") + owners = {} + for component, (_root, items) in enumerate(vectors): + expect(len(items) in (1, 4, 8, 16), "a joint basis vector has an unexpected support size") + for index, coeff in items: + expect(coeff in (-1, 1), "a joint basis coefficient is not a unit") + expect(index not in owners, "joint basis supports overlap") + owners[index] = (component, coeff) + + rotations = [identity] + current = identity + for _ in range(5): + current = compose(current, rotation) + rotations.append(current) + monos = [monomial(element, vectors, owners) for element in rotations] + expect(monos[0] == [(index, 1) for index in range(216)], "the identity rotation does not act as the identity") + columns_p, trace = integer_projector(monos, CHI) + expect(trace == 768, "6P does not have trace 768") + for src, column in enumerate(columns_p): + expect(apply_columns(columns_p, column) == scale_sparse(column, 6), "6P is not idempotent up to the factor 6") + expect(lattice_rank(columns_p) == 128, "the image of 6P does not have rank 128") + _other_columns, other_trace = integer_projector(monos, OTHER_CHI) + expect(other_trace == 0, "the other two-dimensional rational character does not have trace 0") + + basis_columns, pivot_rows = hermite_basis(columns_p) + expect(len(basis_columns) == 128, "the Hermite basis does not have rank 128") + expect(tuple(pivot_rows) == PIVOT_ROWS, "the Hermite pivot rows changed") + diagonals = tuple(basis_columns[index][pivot_rows[index]] for index in range(128)) + expect(diagonals == DIAGONALS, "the Hermite diagonals changed") + generator_rows = [] + for column in columns_p: + row = [0] * 216 + for key, value in column.items(): + row[key] = value + generator_rows.append(row) + express_in_hermite(basis_columns, pivot_rows, generator_rows) + for column in basis_columns: + expect(apply_columns(columns_p, column) == scale_sparse(column, 6), "a Hermite basis vector is not fixed by 6P") + + columns_b = column_table(columns, ports) + rotation_dest, rotation_sign = action_tables(rotation) + reflection_dest, reflection_sign = action_tables(reflection) + expect(operators_commute(columns_b, rotation_dest, rotation_sign), "B does not commute with r") + expect(operators_commute(columns_b, reflection_dest, reflection_sign), "B does not commute with u") + matrix, _sizes = restriction_matrix(columns_b, columns, ports, vectors) + expect(matrix == published_matrix, "the recomputed matrix of B is not the published 216-by-216 matrix") + + images = [] + for column in basis_columns: + image = [0] * 216 + for index, coef in column.items(): + for row in range(216): + image[row] += coef * matrix[row][index] + images.append(image) + x_matrix = express_in_hermite(basis_columns, pivot_rows, images) + for row in range(216): + for col in range(128): + got = sum(basis_columns[k].get(row, 0) * x_matrix[k][col] for k in range(128)) + expect(got == images[col][row], "the integer matrix does not represent B on the Hermite basis") + + polynomial_b = charpoly_integer(x_matrix) + expect(len(polynomial_b) == 129, "chi_B does not have degree 128") + expect( + all(coef % (4 ** (128 - index)) == 0 for index, coef in enumerate(polynomial_b)), + "chi_A is not an integer polynomial", + ) + polynomial_a = [coef // (4 ** (128 - index)) for index, coef in enumerate(polynomial_b)] + expect(scale_to_b(polynomial_a) == polynomial_b, "the relation chi_B(t) = 4^128 chi_A(t/4) failed") + expect(product_factors(COMPONENT_FACTORS) == polynomial_a, "the factored form is not chi_A") + published = {coefficients: multiplicity for multiplicity, coefficients in PUBLISHED_216} + for multiplicity, coefficients in COMPONENT_FACTORS: + expect(coefficients in published, "a factor of chi_A is not one of the published degree-216 factors") + expect(multiplicity <= published[coefficients], "a factor occurs more often on W than on the 216-dimensional space") + expect(monic_divides(list(coefficients), polynomial_a), "a listed factor does not divide chi_A") + for multiplicity, coefficients in PUBLISHED_216: + if multiplicity == 1: + expect(not monic_divides(list(coefficients), polynomial_a), "a multiplicity-one factor of the degree-216 polynomial divides chi_A") + for coefficients in (F11, G11, H14, K14): + prove_irreducible(list(coefficients)) + expect(tuple(-coef for coef in negated_argument(F11)) == G11, "g(x) is not -f(-x)") + expect(negated_argument(H14) == K14, "k(x) is not h(-x)") + root = integer_square_root(polynomial_a) + expect(product_factors(SQUARE_FACTORS) == root, "the degree-64 square root is not the displayed product") + expect(all(root[index] == 0 for index in range(6)), "the square root is not divisible by x^6") + expect(all(root[index] == 0 for index in range(65) if index % 2), "the square root has an odd-degree term") + expect(tuple(root[6 + 2 * index] for index in range(30)) == R_COEFFICIENTS, "the degree-29 compression changed") + expect(poly_sha(polynomial_a) == CHI_A_SHA256, "chi_A SHA256 changed") + expect(poly_sha(polynomial_b) == CHI_B_SHA256, "chi_B SHA256 changed") + expect(poly_sha(root) == CHI_SQRT_SHA256, "the square-root SHA256 changed") + + # Rejection controls. Each must fail closed. + swapped = list(range(8)) + swapped[0], swapped[1] = 1, 0 + expect( + operator(tuple(swapped), (1,) * 8, port_of, lookup, port_basis, coords) is None, + "a coordinate transposition was accepted as a symmetry", + ) + partial = list(range(8)) + partial[4], partial[7] = 7, 4 + expect( + operator(tuple(partial), (1,) * 8, port_of, lookup, port_basis, coords) is None, + "swapping coordinates 4 and 7 alone still looked like a symmetry", + ) + moved = dict(port_of) + moved[columns[0]], moved[columns[1]] = port_of[columns[1]], port_of[columns[0]] + moved_lookup = {norm(column): moved[column] for column in columns} + expect( + operator(G_INV, G_EPS, moved, moved_lookup, port_basis, coords) is None, + "a swapped port label still looked like a symmetry", + ) + nearby = list(R_INV) + nearby[2], nearby[3] = nearby[3], nearby[2] + expect( + operator(tuple(nearby), R_EPS, port_of, lookup, port_basis, coords) is None, + "swapping two values of r still looked like a symmetry", + ) + _dropped, dropped_trace = integer_projector(monos, (2, 1, -1, -2, -1, 0)) + expect(dropped_trace != 768, "dropping the r^5 coefficient still left trace 768") + short = Lattice() + for dest, sign in actions: + short.add_operator(dest, sign, 1) + expect(len(short.basis()) != 216, "dropping g and t still produced the 216-dimensional basis") + corrupted = [row[:] for row in x_matrix] + corrupted[0][0] += 1 + prime = 1_000_003 + genuine = charpoly_mod(x_matrix, prime) + changed = charpoly_mod(corrupted, prime) + expect(genuine[127] != changed[127], "a one-entry change of the integer matrix left chi unchanged") + + print("PASS_A1_ISOTYPIC_128_CHARPOLY") + print("dihedral image order 12") + print("isotypic rank 128") + print("trace of 6P 768") + print("chi_A = q(x)^2") + print("q(x) = x^6 (x-2)^4 (x+2)^4 f(x) g(x) h(x) k(x)") + print("chi_A sha256 " + CHI_A_SHA256) + print("selected kappa 781451/250000000") + + +if __name__ == "__main__": + main()