From 54f4bc0e86d0322b685b47e9dae3575a660f4d49 Mon Sep 17 00:00:00 2001 From: Claude Date: Sun, 12 Jul 2026 06:09:24 +0000 Subject: [PATCH 1/4] =?UTF-8?q?motivic:=20A=5FC,=20the=20C-motivic=20Steen?= =?UTF-8?q?rod=20algebra=20engine=20over=20F=E2=82=82[=CF=84]?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Foundation layer for computing the C-motivic Adams E₂ by deformation. Adds the coefficient ring and the product engine, with no wiring into the resolution engine yet — the mod-τ reduction that the engine resolves comes in a follow-up. - `tau`: F₂[τ] as a small homogeneous scalar (`Tau`). Every structure constant in the motivic world is a single power of τ, so coefficients are a one-integer valuation rather than a polynomial; this type carries the whole τ-tower. - `milnor`: `MotivicMilnorAlgebra` = A_C as a free F₂[τ]-module on the Milnor basis Q(E)P(R). The product is computed two ways — a duality oracle (dualize the coproduct ψ) and the closed-form Kong–Lin Theorem 5.1 (arXiv:2411.12890, ρ = 0) — and the fast path is validated exhaustively against the oracle. It is intentionally not an `Algebra` (that trait is over F_p); it is the F₂[τ] engine the deformation lift builds on. `MOTIVIC_PLAN.md` documents the overall deformation approach; `MOTIVIC_GPU_HANDOFF.md` documents the product engine's batching boundary. Tests: rewrite_tau identities, product associativity, weight-homogeneity, and closed-form-vs-duality agreement over a range. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013ePtYD7Bt4iPeCtmqtqvZE --- MOTIVIC_PLAN.md | 275 +++ ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md | 329 +++ ext/crates/algebra/src/algebra/mod.rs | 3 + .../algebra/src/algebra/motivic/milnor.rs | 1772 +++++++++++++++++ ext/crates/algebra/src/algebra/motivic/mod.rs | 28 + ext/crates/algebra/src/algebra/motivic/tau.rs | 143 ++ 6 files changed, 2550 insertions(+) create mode 100644 MOTIVIC_PLAN.md create mode 100644 ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md create mode 100644 ext/crates/algebra/src/algebra/motivic/milnor.rs create mode 100644 ext/crates/algebra/src/algebra/motivic/mod.rs create mode 100644 ext/crates/algebra/src/algebra/motivic/tau.rs diff --git a/MOTIVIC_PLAN.md b/MOTIVIC_PLAN.md new file mode 100644 index 0000000000..cccc1ec22d --- /dev/null +++ b/MOTIVIC_PLAN.md @@ -0,0 +1,275 @@ +# C-motivic Ext: end-to-end implementation plan + +Computing the **C-motivic Adams E₂ page** — `Ext` over the C-motivic mod-2 +Steenrod algebra, `p = 2`, over `ℂ` — as a first-class computation that reuses +the existing `ext` engine without disturbing the classical path. + +> Base: this branch starts fresh from `master`. The earlier exploratory work +> (a generic-ground-ring engine + `FpTau` `Solvable`) is on +> `73de90ee4671068b3a018640e993ad1a61f0d272` and is recoverable for porting the +> validated *algebra* (`motivic_milnor.rs`), which is the only piece we reuse. + +> **Status (implemented).** Phases 0–4 are done. Phase 0: `algebra::motivic` +> (`Tau`, `MotivicMilnorAlgebra`, `CTauAlgebra`). Phase 1: +> `examples/resolve_motivic_ctau.rs` + `tests/motivic_ctau.rs`. Phases 2–3: +> `ext::motivic` (weight-homogeneous lift, `H(δ)`). Phase 4: +> `examples/resolve_motivic.rs` + the `examples/benchmarks/motivic-S_2` golden +> fixture. **All three anchors pass**: invert `τ` reproduces classical `Ext_A` +> rank-for-rank; `τ = 0` gives the algebraic Novikov ranks; keeping `τ` exhibits +> the `h₁`-tower `τ`-torsion (`h₁⁴`, `h₁⁵`). The classical path is untouched and +> the 69 resolution benchmarks stay bit-identical. One deviation from §2.5: the +> lift and `H(δ)` live in `ext::motivic` (they consume the resolution's +> quasi-inverses, populated by the `ext` engine) rather than the `algebra` crate; +> the pure `Tau` linear algebra could still be factored down later. + +--- + +## 0. Success criteria (the three anchors) + +The motivic E₂ is pinned by three views of one object, so a correct +implementation must reproduce all three: + +1. **invert τ → classical Adams E₂.** Regress against + `ext/examples/benchmarks/num_gens-S_2`. +2. **kill τ (τ = 0) → algebraic Novikov E₂.** The generator counts of the + resolution; equivalently the motivic Adams E₂ of `Cτ` (Gheorghe–Wang–Xu). +3. **keep τ → motivic Adams E₂ with τ-torsion.** The full bigraded answer, + including genuine τ-torsion classes (e.g. `h₁⁵` at `(s,t,w)`), regressed + against a published C-motivic chart (Isaksen–Wang–Xu). + +**Guardrail:** the classical path stays **bit-identical** — all 69 resolution +benchmarks unchanged. We achieve this by *not modifying the `ext` engine at all* +(see §2). + +--- + +## 1. The mathematics we are implementing + +Three algebras, not two. Let `A_C` be the C-motivic dual Steenrod algebra +(Voevodsky): + +``` +A_C = F₂[τ][ξ₁, ξ₂, …, τ₀, τ₁, …] / (τᵢ² = τ·ξᵢ₊₁), |τ| = (0, −1) (stem, weight) +``` + +- **`A_C`** — over `F₂[τ]`, bigraded by (stem, weight). Resolving the sphere + `F₂[τ]` over it gives the **motivic Adams E₂**. +- **`A_C/τ`** — set `τ = 0`: `A_C/τ ≅ F₂[ξᵢ] ⊗ Λ(τᵢ)`. This is a *connected, + finite-type `F₂`-algebra* (every generator has positive stem once τ is gone), + and it is **strictly bigger than the classical `A_* = F₂[ξᵢ]`**. Its Ext is + the **algebraic Novikov E₂**. +- **`A`** — the classical Steenrod algebra. A *further* collapse of `A_C` + (invert τ, or run the algebraic Novikov SS off `A_C/τ`). This is what `ext` + resolves today. + +### Why this makes motivic first-class *and* fast + +The minimal resolution over `A_C` is a genuine object over the motivic Steenrod +algebra. We compute it the way you compute over any τ-adically filtered ring — +**associated graded first, then deform** — which here means: + +- **Leading term:** resolve `A_C/τ` over the field `F₂`. This is a *motivic* + algebra (not classical `A`), resolved by the ordinary field engine. Fast. +- **Deform:** lift the differentials to `A_C` over `F₂[τ]`; the lifted `d∘d` is + zero only mod τ, and a small τ-power (Gaussian-elimination) correction + restores `d² = 0`. This is the honest `A_C` minimal resolution. + +Because the mod-τ object is `A_C/τ` (motivic), this is **not** "compute +classical and correct" — both objects are motivic. First-class semantics, and +the heavy combinatorics live over a field. (This is Guozhen Wang's Adams–Novikov +algorithm, module side; see §7.) + +### Representation: one `F₂[τ]`-matrix per stem, valuation entries + +A homogeneous element of `F₂[τ]` in a fixed bidegree is exactly `τᵏ` for a +unique `k` — store it as the integer `k` (`τ`-valuation), with a sentinel for 0. +Homogeneity fixes every differential entry: generator `j` at weight `wⱼ` maps to +generator `i` at weight `wᵢ` via `τ^{wᵢ − wⱼ}`. So the resolution differential at +each stem is **one matrix over `F₂[τ]`, the same size as the classical `F₂` +matrix**, with integer entries — *not* a sequence of matrices per weight. The +τ-tower is stored once. Linear algebra over it is `F₂` elimination with graded-PID +pivoting (pick the lowest-valuation entry in a column) plus integer bookkeeping — +a small constant over `F₂`, never `×(#weights)`, never polynomial arithmetic. + +### The E₂ is a cohomology, not a generator count + +`Ext_{A_C}(F₂[τ], F₂[τ]) = H(Hom_{A_C}(F_•, F₂[τ]), δ)`, where +`Hom_{A_C}(F_s, F₂[τ]) = Vₛ*` is a **free `F₂[τ]`-module** (rank = # generators) +and `δ` is the **identity-operation part** of the resolution differentials — +apply the augmentation `A_C → F₂[τ]`, which kills positive Steenrod operations +and leaves the τ-power on the unit. Every entry of `δ` is `τᵏ` (or 0), and +`δ² = 0`. The three anchors are three operations on this one complex: + +| operation on τ | result | +|---|---| +| invert τ | classical Adams E₂ | +| set τ = 0 | `δ ≡ 0` ⟹ generator counts = algebraic Novikov E₂ | +| keep τ, take `H(δ)` | motivic Adams E₂, with `F₂[τ]/τᵏ` torsion | + +The first (exploratory) implementation printed the middle row — generator counts +— and never applied `δ`. The real work is `H(δ)`, organised efficiently by the +**τ-Bockstein spectral sequence / Curtis table** (Wang §6): leading-term +cycle/tag table, `a` = permanent cycle, `a→b` = differential whose *length* is +the τ-filtration difference (= torsion order). This keeps the heavy pass over +`F₂` while producing honest `F₂[τ]`-cohomology. + +--- + +## 2. Architecture decisions (committed) + +1. **Do not genericize the `ext` engine.** Motivic factors through the *field* + engine via `A_C/τ`, so the engine stays `F₂`-only and classical stays + bit-identical by construction. (This supersedes / closes PR #263 — the + ground-ring genericization is unnecessary under this factoring.) +2. **Bigrading = internal degree `t` (engine) + weight carried alongside.** + `A_C/τ` is single-`t`-graded with weight a bounded per-basis label; `A_C` + carries weight in its `F₂[τ]` coefficients (the τ-tower). No ℤ²-graded engine. +3. **The deformation lives in `algebra`/the motivic layer, never in `ext`.** +4. **Module picture now; comodules deferred.** BP/Adams–Novikov would add a + *parallel* comodule engine later (§8) — additive, not a rewrite. +5. **Crate placement.** New code in the `algebra` crate (algebras, coefficient + rep, lift) plus, on the `ext` side, the cohomology (`H(δ)`) as a + generalization of the existing `ExtAlgebra` bridge (§6), and a thin + driver/example. The `ext` *resolution engine* and all classical modules are + untouched; `ExtAlgebra` is *extended*, with the classical `δ = 0` path + preserved. + +--- + +## 3. Phase 0 — Coefficient ring and the two algebras + +**Goal:** the `F₂[τ]` valuation ring and the algebras `A_C` (over `F₂[τ]`) and +`A_C/τ` (over `F₂`), with no resolution yet. + +- `algebra/src/algebra/motivic/tau.rs` — `Tau`: homogeneous `F₂[τ]` as a + `(coeff ∈ F₂, valuation ∈ u32)` with a zero sentinel. Ops: add (same + valuation), multiply (add valuations), `valuation`, `shift`, `pow`. Vector / + matrix helpers over `Tau` for the motivic layer only (small, self-contained; + **not** an `fp` replacement threaded through `ext`). +- `algebra/src/algebra/motivic/milnor.rs` — port `MotivicMilnorAlgebra` from + `73de90ee46` (validated: bigrading, Milnor product, coproduct, antipode, + cross-checked against classical). Strip its dependence on the old + `Ring`/`Solvable` tower; it needs only `Tau` arithmetic and the (stem, weight) + bigrading. +- `A_C/τ`: expose the mod-τ reduction as an ordinary `F₂` `Algebra` + (implements the existing `algebra::Algebra` trait), basis by internal degree + `t`, each basis element tagged with its weight. This is the object `ext` + resolves in Phase 1. + +**Tests (port + extend):** product vs classical Milnor after inverting τ; +associativity; coproduct is an algebra map; antipode; `A_C/τ` basis dimensions +per (t, weight) against hand computation for small t. + +--- + +## 4. Phase 1 — Resolve `A_C/τ` over `F₂` → algebraic Novikov E₂ + +**Goal:** first end-to-end numbers, using the **unmodified** `ext` engine. + +- Present `A_C/τ` to `ext` as an `Algebra`; resolve the trivial module `F₂` over + it with the existing minimal-resolution engine. Generators per (s, t) = + algebraic Novikov E₂ ranks; carry the weight tag per generator (needed by the + lift in Phase 2). +- `ext/examples/resolve_motivic_ctau.rs` — a `resolve`-style driver printing the + `A_C/τ` chart. + +**Milestone value:** a real, publishable computation (algebraic Novikov E₂ / +motivic Ext of `Cτ`) with essentially zero engine risk, and it validates the +`A_C/τ` algebra structure before the harder τ work. + +**Tests:** generator counts contain the classical ranks (each classical class +survives); spot-check low-degree algebraic Novikov ranks against Wang/known +values; `d² = 0` in the `F₂` resolution (engine already guarantees this). + +--- + +## 5. Phase 2 — Lift to `A_C` over `F₂[τ]` → the motivic minimal resolution + +**Goal:** the honest minimal resolution over `A_C`, in the valuation rep. + +- Take the Phase-1 `A_C/τ` resolution as the **model**: the generator skeleton + (ranks, weights, and the mod-τ differentials) is fixed. +- Lift each differential entry from `A_C/τ` to `A_C` (choose the `F₂[τ]` lift of + the mod-τ operation). Compute `d∘d`; it is `≡ 0 mod τ`, so it equals τ·(error). + Solve the small `F₂[τ]` linear system (graded-PID elimination on valuation + matrices) that adjusts the lift so `d² = 0` exactly — Wang §4–5. Optionally use + the "project to cogenerators" size reduction later as an optimization. +- Store differentials as per-stem `Tau` matrices (§1 representation). + +**Tests:** generator counts identical to Phase 1 (the lift preserves the +skeleton); `d² = 0` over `F₂[τ]` on every generator through a fixed range; +reducing the lifted differentials mod τ reproduces Phase 1 exactly. + +--- + +## 6. Phase 3 — `H(δ)` in `ExtAlgebra` → the motivic Adams E₂ + +`ExtAlgebra` (`ext/src/ext_algebra/`) is already the intended **resolution → Ext +bridge**: it wraps a resolution and presents `Ext(M, k)` as a bigraded algebra. +Today that bridge is *trivial* — `dimension(b) = number_of_gens_in_bidegree(b)` +— because over a field the Hom-complex differential vanishes (`δ = 0`), so +`Ext = generators`, read off. The motivic case is the **same bridge finally doing +the real work**: `Ext = H(δ)`, with the classical field case (`δ = 0`) recovering +today's behavior. Phase 3 is therefore a *generalization of `ExtAlgebra`*, not a +new extractor. + +**Goal:** teach `ExtAlgebra` that "Ext from a resolution" means taking cohomology, +then read the motivic E₂ off it. + +- Give `ExtAlgebra` a δ-aware Ext: extract `δ` (the identity-operation, τ-power + part of the resolution differentials — apply the augmentation `A_C → F₂[τ]`), a + complex of free `F₂[τ]`-modules with `δ² = 0`, and compute `Ext = H(δ)` as + f.g. `F₂[τ]`-modules (free ⊕ `F₂[τ]/τᵏ`). The current `dimension`/`basis` + become the `δ = 0` specialization. +- Organize the cohomology by the τ-Bockstein Curtis table (Wang §6): leading-term + cycle/tag, `diff_length` = τ-filtration difference = torsion order. This lands + next to the existing `secondary` (d₂) submodule — the τ-BSS / algebraic-Novikov + differentials are *more differentials on Ext*, the same family `ExtAlgebra` + already hosts — so products, d₂, and Massey products keep operating on the + motivic Ext through one surface. +- The two derived views fall out: τ = 0 (drop δ) = algebraic Novikov (Phase 1); + invert τ (rank of δ over `F₂(τ)`) = classical. + +**Tests — the three anchors:** +1. invert τ reproduces `num_gens-S_2` exactly, *including killing τ-torsion* + (`h₁⁵` disappears). +2. τ = 0 reproduces the Phase-1 algebraic Novikov ranks. +3. the full E₂ matches a published C-motivic chart over a fixed range: the + `h₁`-tower is τ-torsion-free/periodic as expected, and named τ-torsion classes + (`h₁⁵` etc.) appear with the correct `F₂[τ]/τᵏ` order. + +--- + +## 7. Phase 4 — Presentation, fixtures, docs + +- `ext/examples/resolve_motivic.rs` — print the bigraded motivic E₂ (ASCII, or + emit a SeqSee/JSON chart with weight and τ-torsion annotations). +- Golden fixtures for the motivic chart (small range) under + `ext/examples/benchmarks/` so the anchors regress in CI. +- Doc pass: a `motivic` module doc explaining the three algebras, the three + pages, and where the deformation lives; keep the model identifier out of any + committed text. + +--- + +## 8. Out of scope (future) + +- **BP / Adams–Novikov.** Wang's algorithm is natively about `BP∗BP`-comodules; + motivic is the `A_C/τ ↔ P` special case. Reaching BP means adding a *parallel* + comodule / Hopf-algebroid resolver (cofree resolutions, cobar `δ`) that shares + the coefficient layer (extended to p-adic `Z/pᵏ`) and the Curtis/SS extractor. + Additive; the module engine here is untouched. +- **McTague's Steenrod-algebra data structure.** A faster cycle/tag reduction + for the `F₂` skeleton (Phase 1). Slot it behind the leading-term-reduction + interface if/when Phase 1 becomes the bottleneck; not on the critical path. + +--- + +## 9. Reference + +- G. Wang, *Computations of the Adams–Novikov E₂-term* (2020) — the algorithm + (`MinimalResolution.pdf` in `joeybf/minimalresolution`; §4–5 model-and-lift, + §6 Curtis table). +- Gheorghe–Wang–Xu, *The special fiber of the motivic deformation…* + (arXiv:1809.09290) — algebraic Novikov E₂ = motivic Adams E₂ of `Cτ`. +- Isaksen–Wang–Xu — C-motivic charts for the anchor regressions. diff --git a/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md b/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md new file mode 100644 index 0000000000..3d54b83ba4 --- /dev/null +++ b/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md @@ -0,0 +1,329 @@ +# Handoff: CubeCL kernel for the C-motivic Milnor multiply (ℂ, p = 2) + +A self-contained brief for building a GPU kernel that offloads the product in the +C-motivic Steenrod algebra `A_C` (and its mod-τ reduction `A_C/τ`). It mirrors the +classical Nassau handoff (`GPU_KERNEL_HANDOFF.md`, branch `nassau_gpu`) and reuses +as much of that machinery as the math allows. Written from a **CPU-only container +that cannot run GPU code**, so everything a GPU box needs — the exact CPU +reference to port, the output encoding, the validation oracle — is captured here. + +The motivic layer lives on branch **`claude/awesome-davinci-khkr8x`** of +`JoeyBF/sseq` (`ext/` workspace), in `ext/crates/algebra/src/algebra/motivic/`. +Pull that first. Read the classical `GPU_KERNEL_HANDOFF.md` alongside this — the +device plumbing (CubeCL toolchain, batching boundary, atomic-XOR output, hash-free +indexing, staging, gotchas) is shared and not repeated in full here. + +--- + +## 0. TL;DR — what to build, and the two facts that make it tractable + +Build a CubeCL kernel that computes, in one launch, the batch of `Q(E₁)P(R₁) · s` +products the motivic resolution issues at a differential-assembly step (the same +"fixed operation × general element" unit as classical; §3). Port the CPU +**two-matrix closed form** `multiply_closed` (Kong–Lin Theorem 5.1, ρ = 0) in +`motivic/milnor.rs`. Two facts collapse most of the apparent extra complexity: + +1. **The device output stays F₂ atomic-XOR — the τ never touches it.** Every + surviving `(X, Y)` pair contributes coefficient exactly `1`, and by weight- + homogeneity a fixed output monomial `z` carries a *single, determined* τ-power + `k(z)` — the weight difference between `a·b` and `z` (magnitude + `|w(a)+w(b)−w(z)|`; `bidegree(...).1` uses a negated-weight convention, so pin + the exact expression against `multiply_closed`, which the §7 oracle checks). So + the kernel computes **F₂ presence** of each output index (atomic XOR, identical + to classical); the τ-power is a host-side weight lookup, not part of the device + accumulator. **You do not need `Tau` arithmetic on the GPU.** + +2. **The resolution's products (the bottleneck) are the mod-τ product, which is + odd-primary-Milnor-structured at p = 2 — so the classical admissible-matrix + kernel ports almost directly.** `A_C/τ ≅ 𝔽₂[ξᵢ] ⊗ E(τᵢ)` has the odd-primary + dual's shape (with `2ⁱ` powers, not `pⁱ`), and mod τ the X-matrix constraint + collapses from `S(X) ≤ R₁` to `S(X) = R₁` exactly — **the classical admissible + matrix**. The only additions over the classical kernel are the second matrix + `Y` (the `τᵢ` / Q-part interaction) and `2ⁱ` weights. The *full* `A_C` product + (with τ, needed only by the Phase-2 lift) needs both matrices and the τ-power. + +So: start from the classical `nassau_gpu` kernel, keep its admissible-matrix X +enumeration and atomic-XOR output, add the Y enumeration and the weight table. +Validate bit-for-bit against `multiply_closed`. + +--- + +## 1. Environment reality + +Same as the classical handoff §1 — everything there applies verbatim (CubeCL +`0.10` + `wgpu` compiles headless; `cubecl-cpu` needs a network fetch that was +403-blocked in the sandbox but works on a normal box and gives a GPU-less +correctness backend; pick `wgpu`/`cuda`/`hip` per hardware; ignore the sandbox's +`-C target-cpu=native` SIGILL). **Nothing in this doc was executed on a GPU.** + +--- + +## 2. Where things live (branch `claude/awesome-davinci-khkr8x`) + +Key file: **`ext/crates/algebra/src/algebra/motivic/milnor.rs`** + +| Item | Symbol | Line (approx.) | +|---|---|---| +| **CPU reference to port** | `multiply_closed` (returns `DualElement`; the fuzz-validated product) | ~851 | +| X enumeration (ξ-matrix) | `Closed::enum_x` | ~702 | +| X→Y bridge (S′, RY targets, T(X)) | `Closed::on_x` | ~734 | +| Y enumeration (τ/Q-matrix) | `enum_y` | ~788 | +| Incremental b / T / S state | `struct Acc` | ~ (near `Closed`) | +| τ-rewrite coefficient `c(S, R)` | `c_coeff` | ~82 | +| Column enumerators | `columns_eq`, `columns_le`, `col0_x_options` | ~588, ~609, ~632 | +| Buffer bound (rows/anti-diagonals) | `const NB` (= 64) | ~656 | +| Product used by the engine | `MotivicMilnorAlgebra::product_indexed` | (returns `Vec<(Tau, usize)>`) | +| Basis index (binary search) | `index_of`, `index_of_ref` | ~945 | +| **Weight of a basis element** | `MotivicMilnorAlgebra::bidegree` (`.1` = weight) | ~950 | +| Product cache (replace on GPU) | `products: RwLock>` | struct field | +| CPU profiling counters | `PRODUCT_HITS/MISSES/NANOS` (env `MOT_PROFILE`) | ~52 | + +Mod-τ reduction (**the resolution algebra**): +**`ext/crates/algebra/src/algebra/motivic/ctau.rs`** + +| Item | Symbol | Line | +|---|---|---| +| `A_C/τ` product = τ⁰ part of `A_C` | `CTauAlgebra::multiply_basis_elements` | ~98 | +| Weight label per basis element | `CTauAlgebra::weight` | ~64 | +| Underlying `A_C` engine | `CTauAlgebra::engine` | ~55 | + +Coefficient ring (host-side only; **not** needed on device): `motivic/tau.rs` +(`Tau` = a τ-valuation, `Option`). + +Full-product consumer (Phase-2 lift): **`ext/src/motivic/mod.rs`** — +`MotivicResolution::lift`/`correct`/`composite` call `engine().product_indexed` +and *do* use the τ-powers. Driver + profiling: `ext/examples/resolve_motivic.rs`. + +Validation oracle (already in the crate): `multiply` (duality product) in +`motivic/milnor.rs`, and the test `test_closed_form_matches_duality_oracle`. + +--- + +## 3. The call chain the products flow through + +Two distinct consumers, different needs: + +``` +Phase 1 — RESOLUTION of A_C/τ (the bottleneck: ~97% of run time; MOD-τ products) + Resolution engine assembles a differential + └ module action → CTauAlgebra::multiply_basis_elements(out, coeff, R_deg, R_idx, s_deg, s_idx) + └ engine.product_indexed(...) → multiply_closed(a, b) [FULL A_C product] + └ keep only terms with τ-valuation 0 (drop τ-divisible) ← mod-τ + +Phase 2 — LIFT to A_C (correcting d²=0; FULL products with τ-powers) + MotivicResolution::correct → composite → engine.product_indexed(...) → multiply_closed +``` + +The fundamental kernel unit is the classical one: **a fixed operation +`a = Q(E₁)P(R₁)` × a general element `s` (a sum of `Q(E₂)P(R₂)` terms) → add +`a·s` to an F₂ output slice.** `a` repeats across many `s` in one assembly step — +that is the amortization a workgroup exploits (one workgroup per distinct `(E₁, +R₁)`; parallelize the term tests across `s`). + +> **Batching boundary — determine on the GPU box.** The classical kernel batches +> per Nassau `get_partial_matrix` (a proven-maximal, non-redundant unit; classical +> handoff §5). The motivic resolution currently uses the **standard** `Resolution` +> engine (`ext/src/motivic/mod.rs`, via `Resolution::new`), not Nassau, so the +> exact assembly boundary that maximizes same-`a` amortization must be located +> there (it is wherever a differential column is built from many algebra actions). +> Two options: (a) find/instrument the standard engine's equivalent of +> `get_partial_matrix`; or (b) port the motivic resolution to the Nassau +> signature algorithm and inherit its boundary and the classical §5 analysis +> wholesale. (b) is the cleaner long-term path and the one that reuses the most. + +> **CPU-side batch primitive already in place.** `MotivicMilnorAlgebra` now caches +> products as one dense [`ProductBlock`] per topological-degree pair `(t1, t2)` +> (`entry (idx1, idx2) = basis[t1][idx1] · basis[t2][idx2]`), with per-entry +> `OnceLock` fill. `MotivicMilnorAlgebra::fill_block(t1, t2)` computes an entire +> block in one (parallel) pass — this is the ready-made host handle for the GPU +> launch: hand the kernel a `(t1, t2)` block, get all `dim1 × dim2` structure +> constants back, store into the block's cells. +> +> **But do not fill whole blocks blindly.** Measured block *density* during a real +> resolution is only **~50%** (stem 30 and 40 both ≈ 50%): the standard engine +> requests roughly half of each touched block's `(idx1, idx2)` grid. So a GPU (or +> CPU) pass that fills the *entire* block does ~2× the necessary +> `multiply_closed` work — and since product compute is ~95% of resolution time, +> eager full-block fill is a **net CPU loss** (confirmed by measurement; that is +> why the CPU path stays lazy-per-entry). For the GPU this means either (i) pass +> the kernel only the *requested* `(idx1, idx2)` sublist of the block (keep the +> engine's demand-driven set), accepting a ragged batch; or (ii) fill the full +> block and eat the 2× overcompute, betting GPU throughput makes it a wash. Decide +> with an A/B on the GPU box; the ragged-batch path (i) is the safer default and +> mirrors what the lazy CPU cache already does. + +--- + +## 4. The math to port (`multiply_closed`, ρ = 0) + +`Q(E₁)P(R₁) · Q(E₂)P(R₂)` is a sum over **two** non-negative integer matrices +`X` (the ξ-part) and `Y` (the τ/Q-part), Kong–Lin Theorem 5.1 at ρ = 0. Sequences +are indexed from 0; **ξ-type** sequences (`R₁, R₂, S(X), R(X), T(X)`) drop index 0 +(ξ₀ = 1), **τ-type** sequences (`E, S(Y), T(Y)`) keep it (τ₀ ≠ 1). This ξ₀ = 1 / τ₀ +≠ 1 asymmetry is the whole subtlety — get it wrong and terms silently vanish. + +Matrix functionals (per matrix `M`, entries `x_{i,j}`, i = row, j = column): +- `S(M)_r = Σⱼ x_{r,j}` (row sum), `R(M)_j = Σᵢ 2ⁱ x_{i,j}` (column weighted sum), + `T(M)_r = Σ_{i+j=r} x_{i,j}` (anti-diagonal sum). +- `b(M) = Πᵣ multinomial(anti-diagonal r) mod 2` = **1 iff the entries on each + anti-diagonal are bitwise disjoint** (no carries). The CPU keeps a running XOR + per anti-diagonal and **rejects a column the instant an entry collides**, so + only `b = 1` matrices are ever visited and `b`/`T`/`S` are never recomputed + (`struct Acc`, `enum_x`/`enum_y`). Port this pruning — it is what makes the + enumeration cheap and it maps to per-anti-diagonal bit masks exactly like the + classical `masks[j]`. + +Enumeration (see `enum_x` → `on_x` → `enum_y`): +- **X** (`enum_x`): column j ≥ 1 ranges over `Σᵢ 2ⁱ x_{i,j} ≤ R₂[j]`; column 0 over + `x_{i,0} ≤ R₁[i]`; prune on the row budget `S(X)_r ≤ R₁[r]` (r ≥ 1) and on + b-collision. **For the mod-τ product the row budget is an equality + `S(X) = R₁` — i.e. the classical admissible matrix of `Sq(R₁)`.** +- **`on_x`**: `S′ = R₁ − S(X)` (ξ-part, index 0 = 0); τ-power `k = Σ(S′)`; the Y + column targets `R(Y)[j] = R₂[j] − R(X)[j]` (j ≥ 1) and the forced + `R(Y)[0] = Σ₂(E₁) + Σ₂(S′) − Σ_{j≥1} R(Y)[j]`; output P-part `T(X)` (drop index 0). +- **Y** (`enum_y`): each column an exact weighted sum `R(Y)[j]`, pruned on + b-collision. At a leaf the pair contributes iff the ρ = 0 degree equation + `popcount(E₁) + 2Σ(S′) = Σ(S(Y))` holds **and** `c(S(Y), S′) ≠ 0` (the τ-rewrite + coefficient `c_coeff`, Theorem 3.4 — pure integer/`binomial2` arithmetic, a + direct GPU port). Output monomial: `Q(E₂ + T(Y)) P(T(X))`, and it must be + square-free (`E₂ + T(Y) ∈ Seq₁`). + +**Coefficient and τ (the §0 fact 1 in detail).** Each surviving `(X, Y)` +contributes `b(X)·b(Y)·c(S(Y),S′)·binom(E₂+T(Y), E₂) = 1` — every factor is 0/1 and +we only reach the leaf when all are 1. So the coefficient of an output monomial is +`#{surviving (X,Y) landing on it} mod 2` — **an F₂ count = atomic XOR**, exactly +classical. The τ-power `k = Σ(S′)` is the *same* for every `(X,Y)` landing on a +given `z` (weight-homogeneity; `k = Σ(S′)`), so recover it host-side from the +weight table (`bidegree(...).1`, sign calibrated against `multiply_closed`) rather +than tracking it on device. **Mod-τ (resolution): keep only `z` with `k(z) = 0`** +— equivalently only enumerate `X` with `S(X) = R₁`. + +**Special case:** `a = 1` (E₁ = 0, R₁ = ∅) ⇒ output is `s` unchanged. Common +(identity operation). + +--- + +## 5. Batching / occupancy + +The classical §5 analysis (batch per assembly launch, one workgroup per distinct +operation, ~60% of term-work fills a 256-wide workgroup on same-operation +amortization at high stems, "widen the scope buys nothing") **carries over** — +the kernel unit is identical. Re-measure for the motivic workload once the +boundary of §3 is fixed; the motivic algebra is ~2× larger and the algebraic- +Novikov E₂ has more generators, so expect *denser* launches than classical at the +same stem (more `s` terms per operation), which only helps occupancy. + +CPU baseline to beat (this sandbox, `--features concurrent`, ~3–4 cores; use for +harness sanity, not as an absolute target): report `(40, 8)` resolution ≈ 2.1 s, +`(50, 10)` ≈ 18 s; the closed-form product is ~62 µs/call at stem 40, ~157 µs at +stem 50, and is ~100% of resolution time (linear algebra is ~3% — **the product is +the bottleneck, exactly as classical**). Regenerate with: + +``` +MOT_PROFILE=1 cargo run --release --no-default-features --features concurrent \ + --example resolve_motivic # prints resolution / products / lift / cohomology +``` + +`PRODUCT_HITS/MISSES/NANOS` (in `motivic/milnor.rs`) report the closed-form call +count and time; the `MOT_PROFILE` env gate is read in `MotivicResolution::new`. + +--- + +## 6. Data layout for the device + +Upload once per algebra (rebuild as `compute_basis` grows): + +- **Basis table + index.** The `A_C/τ` (= `A_C`) basis in a degree is a contiguous + sorted array of monomials `(e_mask: u32, R: [u32; MAX_XI])`, sorted by + `(e_mask, R)`. `index_of_ref` is a **binary search** over it — directly GPU- + portable (no hashing), unlike the classical `HashMap`. That is the simplest + device index and is enough to start. A hash-free *`seqno`-style* rank (classical + handoff §6) is a later optimization: the motivic rank factors as + `(rank of E among Seq₁ subsets at that weight) × (classical-seqno rank of R)` — + more arithmetic than classical but still hash-free; only build it if the binary + search shows up in the profile. +- **Weight table.** One `i32` per basis element = `bidegree(degree, idx).1`. The + host multiplies each present output index `z` by `τ^{k(z)}`, `k(z)` the weight + difference (magnitude `|w(a)+w(b)−w(z)|`; fix the sign once against + `multiply_closed`). For the mod-τ (resolution) kernel, keep only the `k(z) = 0` + slice (`w(z) = w(a)+w(b)`) and the output is pure F₂ — no weight multiply at all. +- **Per-launch inputs.** For each product: the operation `(E₁, R₁)` (fixed within a + workgroup) and `s` as a list of `(E₂, R₂)` terms (or indices into the shared + basis table), each padded to `MAX_XI` `u32`s for coalesced access. +- **Output.** F₂ vectors bit-packed into `Atomic` arrays, one per matrix row; + `atomic_xor(word, 1 << bit)` at `idx = index_of_ref(z)`. Match the crate's `fp` + p = 2 bit order (LSB-first within a limb — verify against `fp`). Multiple + `(X, Y)` pairs, possibly across threads, land on the same index and must cancel + mod 2 — **atomic XOR, not add** (classical gotcha §10 applies identically). + +--- + +## 7. Validation & benchmarking (on the GPU box) + +**Correctness oracle — already in the crate.** +`test_closed_form_matches_duality_oracle` (`motivic/milnor.rs` tests) checks +`multiply_closed` bit-for-bit against the independent duality product `multiply` +for **every** basis pair up to `deg(a)+deg(b) ≤ 18`, including the τ-carrying +(`Sq²·Sq² ∋ τQ₀Q₁`) and multi-Q cases. **Mirror it for the kernel:** for every +`(a, s)` in a degree range, compare the kernel's F₂ output slice (per output +index, i.e. the presence bits) against `multiply_closed(a, bⱼ)` reduced to +presence, and separately check the host-rebuilt τ-powers against +`multiply_closed`'s `Tau` coefficients. With the `cubecl-cpu` runtime this runs +without a GPU. Add a mod-τ variant comparing against +`CTauAlgebra::multiply_basis_elements`. + +**End-to-end A/B.** CPU baseline is the current default (`multiply_closed` behind +the `RwLock` cache). Feature-flag the GPU path in the resolution's product call, +keep device buffers resident across the resolution (upload only the degrees that +grew), and A/B `resolve_motivic` wall time at e.g. `(50,10)` and `(60,10)`. The +GPU must beat the CPU path *including* host↔device transfer. + +--- + +## 8. Suggested staging (each stage compiles + validates before the next) + +1. **Toolchain slice.** `cubecl` behind a `gpu` feature on `algebra` (classical + handoff §9.1 skeleton — a trivial F₂-XOR `#[cube]` kernel + host launch + test). +2. **Index on device.** Port `index_of_ref` as a `#[cube]` binary search over the + uploaded sorted basis; test it is the identity `index_of_ref(basisₖ) == k` over + a full degree. +3. **Single-`(E₁,R₁)` mod-τ multiply.** One workgroup: enumerate the classical + admissible matrices of `Sq(R₁)` (`S(X) = R₁`) — reuse the classical + `AdmissibleMatrix` port from `nassau_gpu` — then the Y (Q-part) enumeration and + the `c_coeff`/ρ leaf test, atomic-XOR the F₂ output. Validate one `(a, s)` + against `CTauAlgebra::multiply_basis_elements`. +4. **Full `A_C` multiply.** Relax `S(X) = R₁` to `S(X) ≤ R₁`, emit at all + τ-powers, rebuild `Tau` on the host from the weight table. Validate against + `multiply_closed` (§7 oracle). +5. **Batched launch + wire into the resolution + A/B** (§3 boundary, §5 occupancy, + §7 A/B). Keep buffers resident. + +--- + +## 9. Gotchas (motivic-specific; the classical §10 all still apply) + +- **ξ₀ = 1 but τ₀ ≠ 1.** ξ-sequences compare on indices ≥ 1, τ-sequences on ≥ 0. + This is the single most error-prone point — the naive Theorem 5.1 reading drops + the `Q₀·P(ξ₁) → Q₀P(ξ₁)` term because it treats `S(X)`'s index 0 as constrained. + The CPU reference and its oracle test pin it; keep the test. +- **τ never goes on the device.** Coefficients are F₂ presence; the τ-power is a + weight-table lookup on the host. Don't build `Tau` arithmetic into the kernel. +- **Two matrices, but Y is small.** X is the classical admissible matrix; Y + handles only the `τᵢ` interaction and is typically tiny (often trivial when + E₁ = E₂ = 0 and the product is pure-ξ). Enumerate X per workgroup (as classical); + the Y loop is a short inner sweep per `(X, term)`. +- **Mod-τ vs full.** The resolution wants `S(X) = R₁` (τ⁰) — a *restriction* of the + full enumeration, so a mod-τ-only kernel is strictly smaller. Ship it first (it + is the bottleneck); the full product is only for the Phase-2 lift. +- **Square-free output.** `E₂ + T(Y)` must be ∈ Seq₁ (each τᵢ at most once); a + collision kills the term. Q indices fit a `u32` mask. +- **`c_coeff` needs `r₀ = 0`.** `S′` has index 0 = 0 by construction; assert it. +- **Trailing-zero canonicalization.** Output P-parts are trimmed before indexing; + the sorted basis assumes trimmed monomials (same as classical §10). + +--- + +*Written from a sandbox without GPU access. The math, the CPU reference, and the +oracle are exercised on CPU (all `motivic` tests green); nothing here ran on a +GPU. Validate stage-by-stage on real hardware, starting from the classical +`nassau_gpu` admissible-matrix kernel — most of it is reused.* diff --git a/ext/crates/algebra/src/algebra/mod.rs b/ext/crates/algebra/src/algebra/mod.rs index 67baed2b04..44c0f4454a 100644 --- a/ext/crates/algebra/src/algebra/mod.rs +++ b/ext/crates/algebra/src/algebra/mod.rs @@ -18,6 +18,9 @@ pub use field::Field; pub mod milnor_algebra; pub use milnor_algebra::MilnorAlgebra; +pub mod motivic; +pub use motivic::{MotivicMilnorAlgebra, Tau}; + mod steenrod_algebra; pub use steenrod_algebra::{AlgebraType, SteenrodAlgebra}; diff --git a/ext/crates/algebra/src/algebra/motivic/milnor.rs b/ext/crates/algebra/src/algebra/motivic/milnor.rs new file mode 100644 index 0000000000..abd3f21326 --- /dev/null +++ b/ext/crates/algebra/src/algebra/motivic/milnor.rs @@ -0,0 +1,1772 @@ +//! The C-motivic (over $\mathbb{C}$, prime 2) dual Steenrod algebra $A_{**}$ and +//! the Steenrod algebra $A_C$ in the Milnor basis, over $\mathbb{F}_2[\tau]$. +//! +//! This follows Kong–Lin, *Product formulas for motivic Milnor basis* +//! (arXiv:2411.12890), specialized to the C-motivic point where $\rho = 0$ (the +//! motivic cohomology of a point over $\mathbb{C}$ is $\mathbb{F}_2[\tau]$). We +//! use the *conjugate* generators of that paper, so the coproduct and +//! Milnor-matrix product match Milnor's classical formulas on the polynomial +//! (`ξ`) part. +//! +//! The dual algebra is $A_{**} = \mathbb{F}_2[\tau][\xi_1, \xi_2, \dots] \otimes +//! E(\tau_0, \tau_1, \dots)$ with $\tau_i^2 = \tau \xi_{i+1}$ (the $\rho = 0$ +//! reduction of $\tau_i^2 = \tau\xi_{i+1} + \rho\tau_{i+1}$). Bidegrees: +//! $|\xi_i| = (2(2^i-1), 2^i-1)$, $|\tau_i| = (2^{i+1}-1, 2^i-1)$, +//! $|\tau| = (0, -1)$ in (stem, weight). +//! +//! Basis elements are $Q(E)P(R)$ with $E \in \mathrm{Seq}_1$ (entries in +//! $\{0,1\}$, encoded as a bitmask: bit $i$ is $\tau_i$) and $R \in \mathrm{Seq}$ +//! (the exponent vector of the $\xi_i$). +//! +//! The product in the Steenrod algebra $A_C$ is computed by **duality**: the +//! dual algebra $A_{**}$ is a commutative $\mathbb{F}_2[\tau]$-algebra whose +//! multiplication ([`dual_mul`]) reduces $\tau_i^2 = \tau\xi_{i+1}$ via +//! [`rewrite_tau`] (Kong–Lin Theorem 3.4), and its coproduct $\psi$ +//! ([`coproduct`], Kong–Lin §2.2) is a $\tau$-free algebra map. The product +//! $a \cdot b$ in $A_C$ is then read off from $\psi$: the coefficient of $z$ in +//! $a \cdot b$ is the coefficient of $\mathrm{mon}(a) \otimes \mathrm{mon}(b)$ in +//! $\psi(\mathrm{mon}(z))$ ([`multiply`]). This is a first, correctness-oriented +//! implementation; the closed-form product (Kong–Lin Theorem 5.1) can replace it +//! later for speed, validated against this one. +//! +//! Coefficients are homogeneous, so each one is a single power of $\tau$: a +//! [`Tau`] scalar (see [`super::tau`]). This is the $A_C$ product **engine**, +//! over $\mathbb{F}_2[\tau]$; the mod-$\tau$ reduction $A_C/\tau$ is presented to +//! the resolution engine as an ordinary $\mathbb{F}_2$-algebra in +//! [`super::ctau`]. +//! +//! Weight convention: the motivic weight of an algebra basis element is the +//! *negative* of the weight of the dual monomial it pairs with, so that products +//! are weight-homogeneous with $\tau$ (weight $-1$) absorbing the difference. + +use std::{ + cell::RefCell, + collections::BTreeMap, + rc::Rc, + sync::{ + Arc, OnceLock, RwLock, + atomic::{AtomicU64, Ordering}, + }, +}; + +/// Profiling counters for [`MotivicMilnorAlgebra::product_indexed`] (cache hits, +/// misses, and nanoseconds spent computing the closed-form product on misses). +pub static PRODUCT_HITS: AtomicU64 = AtomicU64::new(0); +pub static PRODUCT_MISSES: AtomicU64 = AtomicU64::new(0); +pub static PRODUCT_NANOS: AtomicU64 = AtomicU64::new(0); + +use fp::prime::Binomial; +use itertools::Itertools; +use maybe_rayon::prelude::*; +use once::OnceVec; +use rustc_hash::FxHashMap; + +use super::Tau; + +/// $\Sigma(R) = \sum_i r_i$. +fn sigma(r: &[u32]) -> u32 { + r.iter().sum() +} + +/// $\Sigma_2(R) = \sum_i r_i 2^i$. +fn sigma2(r: &[u32]) -> u64 { + r.iter().enumerate().map(|(i, &v)| (v as u64) << i).sum() +} + +/// The coefficient $c(S, R)$ of Notation 3.3 / Theorem 3.4, reduced mod 2. +/// +/// $c(S, R) = \prod_{n \ge 1} \binom{\lfloor \sum_{i=0}^{n-1} 2^{i-n}(s_i - r_i)\rfloor}{r_n}$ +/// when $r_0 = 0$, and $0$ otherwise. The inner floor is computed exactly with integer +/// (Euclidean) division: $\lfloor \sum_{i=0}^{n-1} 2^{i-n}(s_i-r_i)\rfloor +/// = \lfloor (\sum_{i=0}^{n-1} (s_i-r_i)2^i) / 2^n \rfloor$. Each factor is +/// [`Binomial::binomial2`], which is $0$ when the top is negative or below the bottom. +fn c_coeff(s: &[u32], r: &[u32]) -> u32 { + if r.first().copied().unwrap_or(0) > 0 { + return 0; + } + let get = |seq: &[u32], i: usize| seq.get(i).copied().unwrap_or(0) as i32; + let len = s.len().max(r.len()); + let mut prod = 1u32; + for n in 1..len { + let mut num: i32 = 0; + for i in 0..n { + num += (get(s, i) - get(r, i)) << i; + } + let floor = num.div_euclid(1i32 << n); + prod *= i32::binomial2(floor, get(r, n)) as u32; + if prod == 0 { + return 0; + } + } + prod +} + +/// A term $\tau^{\text{tau\_pow}} \cdot Q(E) P(R)$ (mod 2, C-motivic so no $\rho$), the output +/// of rewriting a monomial into the Milnor basis. `e_mask` encodes $E \in \mathrm{Seq}_1$. +#[derive(Debug, Clone, PartialEq, Eq)] +pub struct MotivicTerm { + pub tau_pow: u32, + pub e_mask: u32, + pub r: Vec, +} + +/// Rewrite $\tau(S) = \tau_0^{s_0}\tau_1^{s_1}\cdots$ (with `s` an arbitrary exponent vector) +/// into the C-motivic Milnor basis, using $\tau_i^2 = \tau\xi_{i+1}$ (Theorem 3.4 at +/// $\rho = 0$): +/// $$\tau(S) = \sum_{\substack{E \in \mathrm{Seq}_1,\ R \in \mathrm{Seq}\\ \Sigma_2(E)+\Sigma_2(R)=\Sigma_2(S)\\ \Sigma(E)+2\Sigma(R)=\Sigma(S)}} c(S,R)\, \tau^{\Sigma(R)}\, \tau(E)\xi(R),$$ +/// where the second constraint is the vanishing of the $\rho$-exponent +/// $\Sigma(S) - \Sigma(E) - 2\Sigma(R)$. +pub fn rewrite_tau(s: &[u32]) -> Vec { + let target2 = sigma2(s); + let target_sum = sigma(s); + // Highest index we must consider: any nonzero e_i or r_i contributes at least 2^i, so + // 2^i <= target2. + let max_idx = if target2 == 0 { + 0 + } else { + (63 - target2.leading_zeros()) as usize + }; + + let mut out = Vec::new(); + let mut e = vec![0u32; max_idx + 1]; + let mut r = vec![0u32; max_idx + 1]; + rewrite_tau_dfs(s, 0, max_idx, target2, target_sum, &mut e, &mut r, &mut out); + out +} + +#[allow(clippy::too_many_arguments)] +fn rewrite_tau_dfs( + s: &[u32], + idx: usize, + max_idx: usize, + rem2: u64, // remaining Σ2(E) + Σ2(R) budget + target_sum: u32, // Σ(S), for the ρ = 0 constraint Σ(E) + 2Σ(R) = Σ(S) + e: &mut [u32], + r: &mut [u32], + out: &mut Vec, +) { + if idx > max_idx { + if rem2 != 0 { + return; + } + // ρ-exponent must vanish (C-motivic): Σ(E) + 2Σ(R) = Σ(S). + if sigma(e) + 2 * sigma(r) != target_sum { + return; + } + let coeff = c_coeff(s, r); + if coeff == 0 { + return; + } + let mut e_mask = 0u32; + for (i, &ei) in e.iter().enumerate() { + if ei != 0 { + e_mask |= 1 << i; + } + } + let mut r_trimmed = r.to_vec(); + while let Some(&0) = r_trimmed.last() { + r_trimmed.pop(); + } + out.push(MotivicTerm { + tau_pow: sigma(r), + e_mask, + r: r_trimmed, + }); + return; + } + + let weight = 1u64 << idx; + // e_idx ∈ {0, 1} (E ∈ Seq_1), r_idx ≥ 0, with (e_idx + r_idx) * 2^idx ≤ rem2. + for e_idx in 0..=1u32 { + if (e_idx as u64) * weight > rem2 { + break; + } + let after_e = rem2 - (e_idx as u64) * weight; + let max_r = (after_e / weight) as u32; + for r_idx in 0..=max_r { + e[idx] = e_idx; + r[idx] = r_idx; + rewrite_tau_dfs( + s, + idx + 1, + max_idx, + after_e - (r_idx as u64) * weight, + target_sum, + e, + r, + out, + ); + } + } + e[idx] = 0; + r[idx] = 0; +} + +// --------------------------------------------------------------------------- +// The dual algebra A_** (used to compute the product in A_C by duality). +// --------------------------------------------------------------------------- + +/// An element of the dual algebra $A_{**}$: basis monomials $\tau(E)\xi(R)$ (`E` a square-free +/// bitmask, `R` a trimmed exponent vector) mapped to their $\mathbb{F}_2[\tau]$ coefficients. +/// +/// Everything we compute is homogeneous — $\psi$ and the products are graded — so every +/// coefficient is a single power of $\tau$: a [`Tau`] scalar. Zero coefficients are never stored, +/// so equality is canonical, and coefficient arithmetic is exactly [`Tau`]'s arithmetic — `mul` +/// adds valuations, `add` cancels equal powers mod 2 (unequal powers would be inhomogeneous and +/// cannot arise). +pub type DualElement = BTreeMap<(u32, Vec), Tau>; + +/// Elementwise sum of two exponent vectors, trimmed of trailing zeros. +fn vec_add(a: &[u32], b: &[u32]) -> Vec { + let n = a.len().max(b.len()); + let mut r: Vec = (0..n) + .map(|i| a.get(i).copied().unwrap_or(0) + b.get(i).copied().unwrap_or(0)) + .collect(); + while let Some(&0) = r.last() { + r.pop(); + } + r +} + +/// Add `coeff * key` into a sparse $\mathbb{F}_2[\tau]$-linear combination `acc`, dropping the +/// entry if it cancels to zero. Used for both [`DualElement`] and [`TensorElement`]; the +/// coefficient bookkeeping is entirely [`Tau`] arithmetic. +fn add_term(acc: &mut BTreeMap, key: K, coeff: Tau) { + use std::collections::btree_map::Entry; + if coeff.is_zero() { + return; + } + match acc.entry(key) { + Entry::Occupied(mut o) => { + let sum = *o.get() + coeff; + if sum.is_zero() { + o.remove(); + } else { + *o.get_mut() = sum; + } + } + Entry::Vacant(v) => { + v.insert(coeff); + } + } +} + +/// Multiply two basis monomials of $A_{**}$ and accumulate `coeff * (m1 * m2)` into `acc`. +/// +/// $A_{**}$ is commutative, so $\tau(E_1)\xi(R_1)\cdot\tau(E_2)\xi(R_2)$ is obtained by adding +/// exponents: the exterior parts form $S = E_1 + E_2$ (entries in $\{0,1,2\}$), which is +/// rewritten into the square-free basis via [`rewrite_tau`] ($\tau_i^2 = \tau\xi_{i+1}$), and +/// the resulting $\xi$ exponents are added to $R_1 + R_2$. +fn mul_monomials(m1: &(u32, Vec), m2: &(u32, Vec), coeff: Tau, acc: &mut DualElement) { + let (e1, e2) = (m1.0, m2.0); + let bits = u32::BITS - (e1 | e2).leading_zeros(); + let s: Vec = (0..bits) + .map(|i| ((e1 >> i) & 1) + ((e2 >> i) & 1)) + .collect(); + let r12 = vec_add(&m1.1, &m2.1); + for term in rewrite_tau(&s) { + let r = vec_add(&term.r, &r12); + add_term(acc, (term.e_mask, r), coeff * Tau::power(term.tau_pow)); + } +} + +/// The product of two elements of $A_{**}$. +pub fn dual_mul(a: &DualElement, b: &DualElement) -> DualElement { + let mut out = DualElement::new(); + for (m1, &c1) in a { + for (m2, &c2) in b { + mul_monomials(m1, m2, c1 * c2, &mut out); + } + } + out +} + +/// The generator $\tau_i \in A_{**}$. +pub fn tau_gen(i: usize) -> DualElement { + DualElement::from([((1u32 << i, vec![]), Tau::one())]) +} + +/// The generator $\xi_i \in A_{**}$ (for $i \ge 1$). +pub fn xi_gen(i: usize) -> DualElement { + let mut r = vec![0u32; i + 1]; + r[i] = 1; + DualElement::from([((0, r), Tau::one())]) +} + +/// The unit $1 \in A_{**}$. +pub fn dual_one() -> DualElement { + DualElement::from([((0, vec![]), Tau::one())]) +} + +// --------------------------------------------------------------------------- +// The coproduct ψ: A_** → A_** ⊗ A_**, an algebra map (Kong–Lin §2.2). +// --------------------------------------------------------------------------- + +/// An element of $A_{**} \otimes A_{**}$: pairs of basis monomials → $\mathbb{F}_2[\tau]$ +/// coefficients. Zero coefficients are never stored. +pub type TensorElement = BTreeMap<((u32, Vec), (u32, Vec)), Tau>; + +/// The unit $1 \otimes 1$. +fn tensor_one() -> TensorElement { + TensorElement::from([(((0, vec![]), (0, vec![])), Tau::one())]) +} + +/// Multiply in $A_{**} \otimes A_{**}$: $(a_L \otimes a_R)(b_L \otimes b_R) = (a_L b_L) \otimes (a_R b_R)$. +fn tensor_mul(t1: &TensorElement, t2: &TensorElement) -> TensorElement { + let mut out = TensorElement::new(); + for ((al, ar), &c1) in t1 { + for ((bl, br), &c2) in t2 { + let mut left = DualElement::new(); + mul_monomials(al, bl, Tau::one(), &mut left); + let mut right = DualElement::new(); + mul_monomials(ar, br, Tau::one(), &mut right); + let c = c1 * c2; + for (ml, &cl) in &left { + for (mr, &cr) in &right { + add_term(&mut out, (ml.clone(), mr.clone()), c * (cl * cr)); + } + } + } + } + out +} + +/// The monomial $\xi_j^p$ (with $\xi_0 = 1$). +fn xi_pow_mon(j: usize, p: u32) -> (u32, Vec) { + if j == 0 || p == 0 { + (0, vec![]) + } else { + let mut r = vec![0u32; j + 1]; + r[j] = p; + (0, r) + } +} + +/// $\psi(\tau_k) = 1 \otimes \tau_k + \sum_{i=0}^{k} \tau_i \otimes \xi_{k-i}^{2^i}$. +fn coprod_tau(k: usize) -> TensorElement { + let mut out = TensorElement::new(); + add_term(&mut out, ((0, vec![]), (1 << k, vec![])), Tau::one()); + for i in 0..=k { + add_term( + &mut out, + ((1 << i, vec![]), xi_pow_mon(k - i, 1 << i)), + Tau::one(), + ); + } + out +} + +/// $\psi(\xi_k) = \sum_{i=0}^{k} \xi_i \otimes \xi_{k-i}^{2^i}$ (with $\xi_0 = 1$). +fn coprod_xi(k: usize) -> TensorElement { + let mut out = TensorElement::new(); + for i in 0..=k { + add_term( + &mut out, + (xi_pow_mon(i, 1), xi_pow_mon(k - i, 1 << i)), + Tau::one(), + ); + } + out +} + +/// The coproduct of a basis monomial $\tau(E)\xi(R)$, computed as an algebra map: the product +/// of the coproducts of its generators. +fn coproduct_monomial(e_mask: u32, r: &[u32]) -> TensorElement { + let mut acc = tensor_one(); + for i in 0..u32::BITS { + if (e_mask >> i) & 1 != 0 { + acc = tensor_mul(&acc, &coprod_tau(i as usize)); + } + } + for (j, &rj) in r.iter().enumerate() { + for _ in 0..rj { + acc = tensor_mul(&acc, &coprod_xi(j)); + } + } + acc +} + +/// The coproduct of an arbitrary element of $A_{**}$ (extended $\mathbb{F}_2[\tau]$-linearly). +pub fn coproduct(elt: &DualElement) -> TensorElement { + let mut out = TensorElement::new(); + for ((e, r), &c) in elt { + for (key, cc) in coproduct_monomial(*e, r) { + add_term(&mut out, key, c * cc); + } + } + out +} + +// --------------------------------------------------------------------------- +// The antipode χ: A_** → A_** (conjugate generators ⟷ standard Milnor generators). +// --------------------------------------------------------------------------- + +/// $\xi_j^p \in A_{**}$ (with $\xi_0 = 1$). +fn xi_pow_elt(j: usize, p: u32) -> DualElement { + DualElement::from([(xi_pow_mon(j, p), Tau::one())]) +} + +/// $\chi(\xi_k)$, from the antipode axiom with this module's coproduct: +/// $\chi(\xi_k) = \sum_{i=0}^{k-1} \chi(\xi_i)\,\xi_{k-i}^{2^i}$, $\chi(\xi_0) = 1$. +fn chi_xi(k: usize) -> DualElement { + if k == 0 { + return dual_one(); + } + let mut acc = DualElement::new(); + for i in 0..k { + let term = dual_mul(&chi_xi(i), &xi_pow_elt(k - i, 1 << i)); + for (mon, c) in term { + add_term(&mut acc, mon, c); + } + } + acc +} + +/// $\chi(\tau_k)$, from the antipode axiom: +/// $\chi(\tau_k) = \tau_k + \sum_{i=0}^{k-1} \chi(\tau_i)\,\xi_{k-i}^{2^i}$. +fn chi_tau(k: usize) -> DualElement { + let mut acc = tau_gen(k); + for i in 0..k { + let term = dual_mul(&chi_tau(i), &xi_pow_elt(k - i, 1 << i)); + for (mon, c) in term { + add_term(&mut acc, mon, c); + } + } + acc +} + +/// The antipode $\chi$ of $A_{**}$ — an algebra map (since $A_{**}$ is commutative) extended +/// from `chi_tau`/`chi_xi` on the generators. It converts the paper's conjugate generators to +/// the standard Milnor/Voevodsky generators (`χ(ξ_2) = ξ_2 + ξ_1^3`, etc.) and is an +/// involution. +pub fn antipode(elt: &DualElement) -> DualElement { + let mut out = DualElement::new(); + for ((e, r), &c) in elt { + let mut m = dual_one(); + for i in 0..u32::BITS { + if (e >> i) & 1 != 0 { + m = dual_mul(&m, &chi_tau(i as usize)); + } + } + for (j, &rj) in r.iter().enumerate() { + for _ in 0..rj { + m = dual_mul(&m, &chi_xi(j)); + } + } + for (mon, cc) in m { + add_term(&mut out, mon, c * cc); + } + } + out +} + +// --------------------------------------------------------------------------- +// The product in the Steenrod algebra A_C, computed by dualizing ψ. +// --------------------------------------------------------------------------- + +/// The C-motivic (prime 2) bidegree `(t, w)` of the basis monomial `τ(E)ξ(R)`, in the +/// *paper's* R-convention (`r[j]` is the exponent of `ξ_j`, with `ξ_0 = 1` so `r[0]` is +/// ignored). `|τ_i| = (2^{i+1}-1, 2^i-1)`, `|ξ_j| = (2^{j+1}-2, 2^j-1)`. +fn paper_bidegree(e_mask: u32, r: &[u32]) -> (i32, i32) { + let (mut t, mut w) = (0i32, 0i32); + for i in 0..u32::BITS { + if (e_mask >> i) & 1 != 0 { + t += (1 << (i + 1)) - 1; + w += (1 << i) - 1; + } + } + for (j, &rj) in r.iter().enumerate().skip(1) { + t += rj as i32 * ((1 << (j + 1)) - 2); + w += rj as i32 * ((1 << j) - 1); + } + (t, w) +} + +/// All basis monomials `τ(E)ξ(R)` of a given topological degree `target`. +pub fn enum_basis(target: i32) -> Vec<(u32, Vec)> { + // Generators available up to `target`: τ_i of degree 2^{i+1}-1, ξ_j (j≥1) of degree + // 2^{j+1}-2. Represented as (is_tau, index, degree). + let mut gens: Vec<(bool, usize, i32)> = Vec::new(); + let mut i = 0; + while (1i32 << (i + 1)) - 1 <= target { + gens.push((true, i, (1 << (i + 1)) - 1)); + i += 1; + } + let mut j = 1; + while (1i32 << (j + 1)) - 2 <= target { + gens.push((false, j, (1 << (j + 1)) - 2)); + j += 1; + } + let r_len = gens + .iter() + .filter(|g| !g.0) + .map(|g| g.1 + 1) + .max() + .unwrap_or(0); + + let mut out = Vec::new(); + let mut e = 0u32; + let mut r = vec![0u32; r_len]; + enum_basis_dfs(&gens, 0, target, &mut e, &mut r, &mut out); + out +} + +fn enum_basis_dfs( + gens: &[(bool, usize, i32)], + idx: usize, + rem: i32, + e: &mut u32, + r: &mut [u32], + out: &mut Vec<(u32, Vec)>, +) { + if idx == gens.len() { + if rem == 0 { + let mut r = r.to_vec(); + while let Some(&0) = r.last() { + r.pop(); + } + out.push((*e, r)); + } + return; + } + let (is_tau, index, deg) = gens[idx]; + let max_mult = if is_tau { 1 } else { rem / deg }; + for m in 0..=max_mult { + if m * deg > rem { + break; + } + if is_tau && m == 1 { + *e |= 1 << index; + } else if !is_tau { + r[index] = m as u32; + } + enum_basis_dfs(gens, idx + 1, rem - m * deg, e, r, out); + if is_tau { + *e &= !(1 << index); + } else { + r[index] = 0; + } + } +} + +/// The product `a · b` in the C-motivic Steenrod algebra `A_C`, where `a`, `b` are Milnor basis +/// elements `Q(E)P(R)` given as the monomials `(E, R)` they are dual to. The result is a map +/// from basis monomials to their `𝔽₂[τ]` coefficients. +/// +/// Computed by duality: the coefficient of `z` in `a · b` is the coefficient of +/// `mon(a) ⊗ mon(b)` in `ψ(mon(z))`, summed over the basis `z` of the appropriate +/// topological degree. +pub fn multiply(a: &(u32, Vec), b: &(u32, Vec)) -> DualElement { + let key = (a.clone(), b.clone()); + let t = paper_bidegree(a.0, &a.1).0 + paper_bidegree(b.0, &b.1).0; + let mut out = DualElement::new(); + for z in enum_basis(t) { + if let Some(&c) = coproduct_monomial(z.0, &z.1).get(&key) + && !c.is_zero() + { + out.insert(z, c); + } + } + out +} + +// --------------------------------------------------------------------------- +// The closed-form product (Kong–Lin Theorem 5.1, C-motivic ρ = 0). +// --------------------------------------------------------------------------- +// +// Theorem 5.1 computes Q(E₁)P(R₁)·Q(E₂)P(R₂) as a sum over two Milnor-style +// matrices X, Y with the output monomial Q(E₂ + T(Y)) P(T(X)). The subtlety that +// makes the motivic case "much more complex than classical" (and that trips a +// naive reading) is the ξ₀ = 1 / τ₀ ≠ 1 asymmetry (Remark 4.3): **ξ-type** +// sequences (R₁, R₂, S(X), R(X), T(X)) are compared on indices ≥ 1 — index 0 is +// absorbed — while **τ-type** sequences (E, S(Y), T(Y)) keep index 0. +// +// The τ-rewriting piece of the formula — the coefficient c(S(Y), R₁−S(X)), the +// Σ₂ constraint, and the ρ = 0 filter — is exactly [`rewrite_tau`], which we +// reuse rather than re-derive. + +/// All columns `(v₀, v₁, …)` of non-negative integers with `Σᵢ 2ⁱ vᵢ == target`. +fn columns_eq(target: u32) -> Vec> { + if target == 0 { + return vec![vec![]]; + } + let mut out = Vec::new(); + let mut v0 = target % 2; + while v0 <= target { + for rest in columns_eq((target - v0) / 2) { + let mut col = vec![v0]; + col.extend(rest); + while let Some(&0) = col.last() { + col.pop(); + } + out.push(col); + } + v0 += 2; + } + out +} + +/// All columns with `Σᵢ 2ⁱ vᵢ ≤ bound` (distinct weighted sums, so no duplicates). +fn columns_le(bound: u32) -> Vec> { + (0..=bound).flat_map(columns_eq).collect() +} + +thread_local! { + /// Per-thread memo of [`columns_eq`] / [`columns_le`], keyed by target/bound. + /// The candidate lists are pure functions of one integer and are requested + /// over and over across products; caching them avoids re-enumerating and + /// re-allocating. (`Rc` is confined to a single product's call stack.) + static COLS_EQ: RefCell>>>> = RefCell::new(FxHashMap::default()); + static COLS_LE: RefCell>>>> = RefCell::new(FxHashMap::default()); +} + +fn columns_eq_cached(target: u32) -> Rc>> { + COLS_EQ.with(|c| Rc::clone(c.borrow_mut().entry(target).or_insert_with(|| Rc::new(columns_eq(target))))) +} + +fn columns_le_cached(bound: u32) -> Rc>> { + COLS_LE.with(|c| Rc::clone(c.borrow_mut().entry(bound).or_insert_with(|| Rc::new(columns_le(bound))))) +} + +/// The column-0 options for the `X` matrix: `x_{0,0} = 0`, and `x_{i,0} ≤ R₁[i]` +/// for `i ≥ 1` (bounded because it feeds the row sum `S(X)_i ≤ R₁[i]`). +fn col0_x_options(r1: &[u32]) -> Vec> { + let mut result = vec![vec![0u32]]; + for &bound in r1.iter().skip(1) { + let mut next = Vec::new(); + for base in &result { + for v in 0..=bound { + let mut c = base.clone(); + c.push(v); + next.push(c); + } + } + result = next; + } + for c in &mut result { + while let Some(&0) = c.last() { + c.pop(); + } + } + result +} + +/// Buffer size for row / column / anti-diagonal indices. All are bounded by +/// ~2·log₂(degree) + (number of ξ generators), far below this for any feasible +/// stem (a ξ or τ at index `i` already has degree ≥ 2^i). +const NB: usize = 64; + +fn trim(mut v: Vec) -> Vec { + while let Some(&0) = v.last() { + v.pop(); + } + v +} + +/// Immutable context for the closed-form product recursion (Theorem 5.1, ρ = 0). +struct Closed<'a> { + r1v: &'a [u32], + r2v: &'a [u32], + l: usize, + e1_mask: u32, + e2_mask: u32, + sigma2_e1: i64, + /// Per-column candidate columns for the `X` matrix. + x_cands: &'a [Rc>>], +} + +/// Incrementally maintained matrix state: row sums `S`, a running XOR per +/// anti-diagonal (`or`, which equals the OR because we only keep disjoint entries) +/// for the `b`-coefficient pruning, and the anti-diagonal sums `T`. +struct Acc { + rows: [u32; NB], + or: [u32; NB], + sum: [u32; NB], +} + +impl Acc { + fn zero() -> Self { + Self { + rows: [0; NB], + or: [0; NB], + sum: [0; NB], + } + } +} + +impl<'a> Closed<'a> { + /// Enumerate the `X` matrix column by column. A column is rejected — pruning + /// the whole subtree — if a row sum would exceed `R₁` (rows ≥ 1; row 0 is + /// absorbed by ξ₀ = 1) or if any entry collides on its anti-diagonal (which + /// would make `b(X) = 0`). So only `b = 1` matrices are visited, and `b`/`T`/`S` + /// never need recomputing. The accumulator `acc` and `xcols` are reused. + fn enum_x(&self, j: usize, xcols: &mut Vec<&'a [u32]>, acc: &mut Acc, out: &mut DualElement) { + if j == self.l { + self.on_x(xcols, acc, out); + return; + } + for cand in self.x_cands[j].iter() { + let fits = cand.iter().enumerate().all(|(i, &v)| { + v == 0 + || (!(i >= 1 && acc.rows[i] + v > self.r1v.get(i).copied().unwrap_or(0)) + && acc.or[i + j] & v == 0) + }); + if !fits { + continue; + } + for (i, &v) in cand.iter().enumerate() { + acc.rows[i] += v; + acc.or[i + j] ^= v; + acc.sum[i + j] += v; + } + xcols.push(cand); + self.enum_x(j + 1, xcols, acc, out); + xcols.pop(); + for (i, &v) in cand.iter().enumerate() { + acc.rows[i] -= v; + acc.or[i + j] ^= v; + acc.sum[i + j] -= v; + } + } + } + + /// A complete `X` (already `b = 1`): form `S′`, the `Y` column targets, and the + /// output P-part `T(X) = acc.sum`, then enumerate `Y`. + fn on_x(&self, xcols: &[&'a [u32]], acc: &Acc, out: &mut DualElement) { + // S′ = R₁ − S(X) on indices ≥ 1 (S′[0] = 0); τ-power = Σ(S′). + let sprime: Vec = (0..self.l) + .map(|j| if j == 0 { 0 } else { self.r1v[j] - acc.rows[j] }) + .collect(); + let sigma2_sprime: i64 = sprime.iter().enumerate().map(|(j, &v)| (v as i64) << j).sum(); + let sigma_sprime: u32 = sprime.iter().sum(); + + // RY targets: column j ≥ 1 of Y has weighted sum R₂[j] − R(X)[j]. + let mut ry: Vec = vec![0; self.l]; + for (j, ry_j) in ry.iter_mut().enumerate().skip(1) { + let rxj: u32 = xcols[j].iter().enumerate().map(|(i, &v)| v << i).sum(); + *ry_j = self.r2v[j] as i64 - rxj as i64; + if *ry_j < 0 { + return; + } + } + // Σ₂(S(Y)) = Σ₂(E₁) + Σ₂(S′) is forced and equals Σ_{j≥0} R(Y)[j], so column + // 0 of Y has a determined weighted sum. + let ry0 = (self.sigma2_e1 + sigma2_sprime) - ry[1..].iter().sum::(); + if ry0 < 0 { + return; + } + + // The output P-part T(X) (index 0 dropped: ξ₀ = 1). + let out_r = trim(std::iter::once(0).chain((1..NB).map(|d| acc.sum[d])).collect()); + + let mut y_cands: Vec>>> = Vec::with_capacity(self.l); + y_cands.push(columns_eq_cached(ry0 as u32)); + for &t in ry.iter().skip(1) { + y_cands.push(columns_eq_cached(t as u32)); + } + + let mut ycols: Vec<&[u32]> = Vec::with_capacity(self.l); + let mut yacc = Acc::zero(); + enum_y( + &y_cands, + 0, + &mut ycols, + &mut yacc, + self.e1_mask, + self.e2_mask, + &sprime, + sigma_sprime, + &out_r, + out, + ); + } +} + +/// Enumerate the `Y` matrix column by column (each column an exact weighted sum), +/// pruning on anti-diagonal collisions (`b(Y) = 0`) and accumulating each matching +/// contribution. +#[allow(clippy::too_many_arguments)] +fn enum_y<'b>( + y_cands: &'b [Rc>>], + j: usize, + ycols: &mut Vec<&'b [u32]>, + acc: &mut Acc, + e1_mask: u32, + e2_mask: u32, + sprime: &[u32], + sigma_sprime: u32, + out_r: &[u32], + out: &mut DualElement, +) { + if j == y_cands.len() { + // The rewriting τ(S(Y)) contains the term with exterior part E₁ and ξ-part + // S′ iff the ρ = 0 degree equation Σ(E₁) + 2Σ(S′) = Σ(S(Y)) holds (with + // Σ(E₁) = popcount) and c(S(Y), S′) ≠ 0. (Σ₂(E₁) + Σ₂(S′) = Σ₂(S(Y)) is + // forced by the Y-column construction.) That term contributes at τ^{Σ(S′)}. + let sigma_sy: u32 = acc.rows.iter().sum(); + if e1_mask.count_ones() + 2 * sigma_sprime != sigma_sy { + return; + } + let sy_len = (0..NB).rev().find(|&i| acc.rows[i] > 0).map_or(0, |i| i + 1); + if c_coeff(&acc.rows[..sy_len], sprime) == 0 { + return; + } + // E_out = E₂ + T(Y) must be square-free (Seq₁). Q indices fit in a u32. + let mut out_e_mask = 0u32; + for i in 0..32 { + let val = ((e2_mask >> i) & 1) + acc.sum[i]; + if val > 1 { + return; + } + if val == 1 { + out_e_mask |= 1 << i; + } + } + add_term(out, (out_e_mask, out_r.to_vec()), Tau::power(sigma_sprime)); + return; + } + for cand in y_cands[j].iter() { + if !cand.iter().enumerate().all(|(i, &v)| acc.or[i + j] & v == 0) { + continue; + } + for (i, &v) in cand.iter().enumerate() { + acc.rows[i] += v; + acc.or[i + j] ^= v; + acc.sum[i + j] += v; + } + ycols.push(cand); + enum_y( + y_cands, j + 1, ycols, acc, e1_mask, e2_mask, sprime, sigma_sprime, out_r, out, + ); + ycols.pop(); + for (i, &v) in cand.iter().enumerate() { + acc.rows[i] -= v; + acc.or[i + j] ^= v; + acc.sum[i + j] -= v; + } + } +} + +/// The product `a · b` in `A_C` via Kong–Lin Theorem 5.1 (ρ = 0). Same contract +/// as [`multiply`] (the duality oracle it is validated against). +pub fn multiply_closed(a: &(u32, Vec), b: &(u32, Vec)) -> DualElement { + let (e1_mask, r1) = (a.0, &a.1); + let (e2_mask, r2) = (b.0, &b.1); + let l = r1.len().max(r2.len()).max(1); + let r1v: Vec = (0..l).map(|j| r1.get(j).copied().unwrap_or(0)).collect(); + let r2v: Vec = (0..l).map(|j| r2.get(j).copied().unwrap_or(0)).collect(); + + // Per-column candidate columns for X: column 0 bounded by R₁ rows, columns + // j ≥ 1 by R₂[j] weighted. + let mut x_cands: Vec>>> = Vec::with_capacity(l); + x_cands.push(Rc::new(col0_x_options(&r1v))); + for &bound in r2v.iter().skip(1) { + x_cands.push(columns_le_cached(bound)); + } + + let ctx = Closed { + r1v: &r1v, + r2v: &r2v, + l, + e1_mask, + e2_mask, + sigma2_e1: e1_mask as i64, // Σ (bit i)·2ⁱ = e1_mask + x_cands: &x_cands, + }; + let mut out = DualElement::new(); + let mut xcols: Vec<&[u32]> = Vec::with_capacity(l); + let mut acc = Acc::zero(); + ctx.enum_x(0, &mut xcols, &mut acc, &mut out); + out +} + +// --------------------------------------------------------------------------- +// The C-motivic Steenrod algebra A_C as a free 𝔽₂[τ]-module on the Milnor basis. +// --------------------------------------------------------------------------- + +/// The C-motivic (prime 2) Steenrod algebra $A_C$, presented as a free +/// $\mathbb{F}_2[\tau]$-module on the Milnor basis $\{Q(E)P(R)\}$ with lazy +/// per-(topological-)degree basis indexing. +/// +/// This is the $A_C$ product **engine** over $\mathbb{F}_2[\tau]$: the dual-based product is +/// exposed through [`product_indexed`](MotivicMilnorAlgebra::product_indexed), whose coefficients +/// are [`Tau`] scalars. It is deliberately not an [`Algebra`](crate::algebra::Algebra) +/// implementation, because that trait is over $\mathbb{F}_p$; the mod-$\tau$ reduction, which *is* +/// such an algebra, lives in [`super::ctau`], and the honest $\mathbb{F}_2[\tau]$ resolution is +/// built by lifting against this engine (Phase 2). +/// +/// Weight convention: [`bidegree`](MotivicMilnorAlgebra::bidegree) returns `(t, w)` where `t` is +/// the topological degree and `w` is the motivic weight in the presentation where $\tau$ has +/// weight $-1$ and products are weight-homogeneous (i.e. `w = -(dual monomial weight)`; see the +/// module-level note). +#[derive(Default)] +pub struct MotivicMilnorAlgebra { + /// `basis[t]` is the $\mathbb{F}_2[\tau]$-basis in topological degree `t`, sorted for stable + /// indexing. + basis: OnceVec)>>, + /// Memoized basis-element products, one dense [`ProductBlock`] per pair of + /// topological degrees `(t1, t2)`. The duality product is expensive and a + /// resolution asks for the same structure constants repeatedly, so we cache + /// them (the role the classical Milnor algebra's `cache-multiplication` table + /// plays). Blocking by degree pair — rather than a flat `(t1, idx1, t2, idx2)` + /// map — makes the cache the natural *batch unit* (see [`ProductBlock`]): the + /// per-degree-pair registry is small and read-mostly, and each block's entries + /// fill through independent [`OnceLock`]s, so a concurrent resolution reads + /// hits lock-free instead of serializing on one global product lock. + blocks: RwLock>>, +} + +/// A dense block of basis-element products for one pair of topological degrees +/// `(t1, t2)`. Entry `(idx1, idx2)`, at flat position `idx1 * dim2 + idx2`, is the +/// product of the `idx1`-th basis element in degree `t1` with the `idx2`-th in +/// degree `t2` — the `(Tau, index-in-degree-(t1+t2))` list [`multiply_closed`] +/// returns. +/// +/// This is the natural **batch unit** for the product: every structure constant +/// the resolution needs for a given degree pair lives in one block. Entries fill +/// lazily (each an independent [`OnceLock`], so reads are lock-free and a +/// concurrent resolution never serializes on a global product lock), and +/// [`MotivicMilnorAlgebra::fill_block`] computes the whole block at once — the +/// shape a GPU kernel would fill in a single launch. +pub struct ProductBlock { + /// The number of basis elements in degree `t2` (the block's row stride). + dim2: usize, + entries: Vec>>, +} + +impl MotivicMilnorAlgebra { + pub fn new() -> Self { + Self { + basis: OnceVec::new(), + blocks: RwLock::new(FxHashMap::default()), + } + } + + /// Compute and cache the basis in every topological degree up to and including `degree`. + /// Idempotent and cheap to re-call. + pub fn compute_basis(&self, degree: i32) { + for t in self.basis.len() as i32..=degree { + let mut b = enum_basis(t); + b.sort(); + self.basis.push(b); + } + } + + /// The $\mathbb{F}_2[\tau]$-rank of `A_C` in topological degree `degree`. + pub fn dimension(&self, degree: i32) -> usize { + if degree < 0 { + return 0; + } + self.basis[degree as usize].len() + } + + /// The `idx`-th basis monomial `(E, R)` in degree `degree`. + pub fn basis_element(&self, degree: i32, idx: usize) -> &(u32, Vec) { + &self.basis[degree as usize][idx] + } + + /// The index of a basis monomial in its degree, if present. + pub fn index_of(&self, degree: i32, elt: &(u32, Vec)) -> Option { + self.basis[degree as usize].binary_search(elt).ok() + } + + /// The `(topological degree, motivic weight)` of a basis element. + pub fn bidegree(&self, degree: i32, idx: usize) -> (i32, i32) { + let (e, r) = &self.basis[degree as usize][idx]; + let (t, w) = paper_bidegree(*e, r); + (t, -w) + } + + /// Parse a basis element written as `basis_element_to_string` prints it — + /// `1`, `Q_i`, `P(r_0, r_1, …)`, or a space-separated product `Q_i … P(R)` — + /// into its `(topological degree, index)`. The inverse of + /// [`Self::basis_element_to_string`]; returns `None` on a malformed string or a + /// monomial absent from the basis. This is what lets `.json` module descriptors + /// be written over the motivic Steenrod algebra. + pub fn basis_element_from_string(&self, elt: &str) -> Option<(i32, usize)> { + let mut e_mask: u32 = 0; + let mut r: Vec = Vec::new(); + // The `P(…)` block (if any) is last and may contain spaces after commas, so + // split it off before whitespace-tokenizing the `Q_i` factors. + let (q_part, p_part) = match elt.split_once("P(") { + Some((q, rest)) => (q, Some(rest.strip_suffix(')')?)), + None => (elt, None), + }; + for token in q_part.split_whitespace() { + if token == "1" { + continue; // the unit factor + } + let i: u32 = token.strip_prefix("Q_")?.parse().ok()?; + e_mask |= 1 << i; + } + if let Some(p) = p_part { + r = p + .split(',') + .map(|x| x.trim().parse::()) + .collect::>() + .ok()?; + } + + let degree = paper_bidegree(e_mask, &r).0; + if degree < 0 { + return None; + } + self.compute_basis(degree); + // Match up to trailing zeros in the ξ-exponent vector (its canonical length + // is degree-dependent). + let trim = |v: &[u32]| { + let mut v = v.to_vec(); + while v.last() == Some(&0) { + v.pop(); + } + v + }; + let want = trim(&r); + (0..self.dimension(degree)) + .find(|&idx| { + let (e2, r2) = self.basis_element(degree, idx); + *e2 == e_mask && trim(r2) == want + }) + .map(|idx| (degree, idx)) + } + + /// A display string for a basis element (`Q_i … P(R)`). + pub fn basis_element_to_string(&self, degree: i32, idx: usize) -> String { + let (e_mask, r) = self.basis_element(degree, idx); + let mut parts = Vec::new(); + for i in fp::prime::iter::BitflagIterator::set_bit_iterator(*e_mask as u64) { + parts.push(format!("Q_{i}")); + } + if !r.is_empty() { + parts.push(format!("P({})", r.iter().format(", "))); + } + if parts.is_empty() { + "1".to_string() + } else { + parts.join(" ") + } + } + + /// The dense [`ProductBlock`] for the degree pair `(t1, t2)`, created (with all + /// entries empty) on first request and shared thereafter. Creation is the only + /// step that touches the block registry's write lock; entry computation happens + /// lock-free through the block's [`OnceLock`]s. + fn block(&self, t1: i32, t2: i32) -> Arc { + if let Some(block) = self.blocks.read().unwrap().get(&(t1, t2)) { + return Arc::clone(block); + } + // Ensure the operand and output bases exist before sizing/indexing. + self.compute_basis(t1); + self.compute_basis(t2); + self.compute_basis(t1 + t2); + let dim1 = self.dimension(t1); + let dim2 = self.dimension(t2); + let mut w = self.blocks.write().unwrap(); + Arc::clone(w.entry((t1, t2)).or_insert_with(|| { + Arc::new(ProductBlock { + dim2, + entries: (0..dim1 * dim2).map(|_| OnceLock::new()).collect(), + }) + })) + } + + /// The closed-form product (Kong–Lin Theorem 5.1) of the two basis elements at + /// `(t1, idx1)` and `(t2, idx2)`, mapped into the degree-`(t1+t2)` basis. The + /// bases must already be computed (the caller through [`Self::block`] ensures + /// this). + fn compute_product(&self, t1: i32, idx1: usize, t2: i32, idx2: usize) -> Vec<(Tau, usize)> { + let a = &self.basis[t1 as usize][idx1]; + let b = &self.basis[t2 as usize][idx2]; + let t = t1 + t2; + multiply_closed(a, b) + .into_iter() + .map(|(z, c)| { + ( + c, + self.index_of(t, &z) + .expect("product landed outside the basis"), + ) + }) + .collect() + } + + /// The product of two basis elements, as an $\mathbb{F}_2[\tau]$-linear combination of basis + /// elements in degree `t1 + t2`: a list of `(coefficient, index)` pairs. The product of two + /// homogeneous basis elements is weight-homogeneous, so each coefficient is a single power of + /// $\tau$ — a [`Tau`] scalar. + pub fn product_indexed( + &self, + t1: i32, + idx1: usize, + t2: i32, + idx2: usize, + ) -> Vec<(Tau, usize)> { + self.product_indexed_with(t1, idx1, t2, idx2, <[_]>::to_vec) + } + + /// Run `f` on the cached product of the two basis elements **by reference**, + /// without cloning the term list. A hot consumer (the Phase 2 lift's + /// `composite`, which walks millions of cached products) should prefer this to + /// [`Self::product_indexed`] to avoid a `Vec` clone per structure constant. + pub fn product_indexed_with( + &self, + t1: i32, + idx1: usize, + t2: i32, + idx2: usize, + f: impl FnOnce(&[(Tau, usize)]) -> R, + ) -> R { + let block = self.block(t1, t2); + let cell = &block.entries[idx1 * block.dim2 + idx2]; + if let Some(cached) = cell.get() { + PRODUCT_HITS.fetch_add(1, Ordering::Relaxed); + return f(cached); + } + // The closed-form product (Kong–Lin Theorem 5.1), validated exhaustively + // against the duality oracle `multiply`. + let terms = cell.get_or_init(|| { + let start = std::time::Instant::now(); + let result = self.compute_product(t1, idx1, t2, idx2); + PRODUCT_MISSES.fetch_add(1, Ordering::Relaxed); + PRODUCT_NANOS.fetch_add(start.elapsed().as_nanos() as u64, Ordering::Relaxed); + result + }); + f(terms) + } + + /// Compute *every* entry of the `(t1, t2)` product block at once, filling it in + /// parallel (under the `concurrent` feature) — the batch shape a GPU kernel + /// fills in a single launch. Idempotent: entries already present are kept, and + /// a second call is a cheap no-op. This is the batching boundary the GPU port + /// hooks into; on CPU it lets a caller warm a whole degree pair up front + /// instead of paying a lock-free-but-serial miss per structure constant. + pub fn fill_block(&self, t1: i32, t2: i32) { + let block = self.block(t1, t2); + let dim2 = block.dim2; + if dim2 == 0 || block.entries.is_empty() { + return; + } + let dim1 = block.entries.len() / dim2; + (0..dim1).into_maybe_par_iter().for_each(|idx1| { + for idx2 in 0..dim2 { + block.entries[idx1 * dim2 + idx2].get_or_init(|| { + let start = std::time::Instant::now(); + let result = self.compute_product(t1, idx1, t2, idx2); + PRODUCT_MISSES.fetch_add(1, Ordering::Relaxed); + PRODUCT_NANOS.fetch_add(start.elapsed().as_nanos() as u64, Ordering::Relaxed); + result + }); + } + }); + } +} + +impl std::fmt::Display for MotivicMilnorAlgebra { + fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result { + write!(f, "MotivicMilnorAlgebra(C, p=2)") + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn test_algebra_basis_and_multiply() { + let alg = MotivicMilnorAlgebra::new(); + alg.compute_basis(8); + + // Degree 0 is the unit; degrees 1 and 2 are 1-dimensional (Q_0 and P(ξ_1)). + assert_eq!(alg.dimension(0), 1); + assert_eq!(alg.basis_element(0, 0), &(0u32, vec![])); + assert_eq!(alg.dimension(1), 1); + assert_eq!(alg.basis_element(1, 0), &(0b1u32, vec![])); // Q_0 + assert_eq!(alg.dimension(2), 1); + assert_eq!(alg.basis_element(2, 0), &(0u32, vec![0, 1])); // P(ξ_1) + + // bidegree: Q_0 is (1, 0), P(ξ_1) is (2, -1) in this presentation. + assert_eq!(alg.bidegree(1, 0), (1, 0)); + assert_eq!(alg.bidegree(2, 0), (2, -1)); + + // Q_0 · P(ξ_1) = Q_1 + Q_0 P(ξ_1), reconstructed from indices. + let terms: DualElement = alg + .product_indexed(1, 0, 2, 0) + .into_iter() + .map(|(c, idx)| (alg.basis_element(3, idx).clone(), c)) + .collect(); + assert_eq!( + terms, + DualElement::from([((0b10, vec![]), Tau::one()), ((0b1, vec![0, 1]), Tau::one())]) + ); + } + + #[test] + fn test_algebra_multiply_matches_raw_and_is_homogeneous() { + let alg = MotivicMilnorAlgebra::new(); + alg.compute_basis(12); + for t1 in 0..=6 { + for idx1 in 0..alg.dimension(t1) { + for t2 in 0..=6 { + for idx2 in 0..alg.dimension(t2) { + let a = alg.basis_element(t1, idx1).clone(); + let b = alg.basis_element(t2, idx2).clone(); + // Indexed product agrees with the raw monomial product. + let indexed: DualElement = alg + .product_indexed(t1, idx1, t2, idx2) + .into_iter() + .map(|(c, idx)| (alg.basis_element(t1 + t2, idx).clone(), c)) + .collect(); + assert_eq!(indexed, multiply(&a, &b)); + + // Weight-homogeneous: w_out - (τ-power) = w_a + w_b. + let (wa, wb) = (alg.bidegree(t1, idx1).1, alg.bidegree(t2, idx2).1); + for (c, idx) in alg.product_indexed(t1, idx1, t2, idx2) { + let w = alg.bidegree(t1 + t2, idx).1; + assert_eq!(w - c.valuation().unwrap() as i32, wa + wb); + } + } + } + } + } + } + + #[test] + fn test_dual_mul_relations() { + // τ_0^2 = τ ξ_1 (the defining relation at ρ = 0). + assert_eq!( + dual_mul(&tau_gen(0), &tau_gen(0)), + DualElement::from([((0, vec![0, 1]), Tau::power(1))]) // τ^1 · ξ_1 + ); + // ξ_1^2 is just the monomial ξ_1^2. + assert_eq!( + dual_mul(&xi_gen(1), &xi_gen(1)), + DualElement::from([((0, vec![0, 2]), Tau::one())]) + ); + // Distinct τ's commute and stay square-free: τ_0 τ_1. + assert_eq!( + dual_mul(&tau_gen(0), &tau_gen(1)), + DualElement::from([((0b11, vec![]), Tau::one())]) + ); + // τ_1^2 = τ ξ_2. + assert_eq!( + dual_mul(&tau_gen(1), &tau_gen(1)), + DualElement::from([((0, vec![0, 0, 1]), Tau::power(1))]) + ); + // Multiplication by the unit is the identity. + assert_eq!(dual_mul(&dual_one(), &tau_gen(2)), tau_gen(2)); + } + + #[test] + fn test_pure_xi_matches_classical_milnor() { + // Motivically the even part is NOT closed: P(R_1)·P(R_2) can leak into Q-terms + // carrying τ (e.g. Sq^2·Sq^2 ∋ τ Q_0 Q_1). But the pure-ξ *outputs* (E = 0) are τ-free + // and must agree with the ordinary mod-2 Milnor product. Compare those against the + // codebase's `MilnorAlgebra` (paper R = [0, p_part...]; drop the ξ_0 slot to convert). + use fp::{prime::TWO, vector::FpVector}; + + use crate::algebra::{ + Algebra, + milnor_algebra::{MilnorAlgebra, MilnorBasisElement}, + }; + + let alg = MilnorAlgebra::new(TWO, false); + let trim = |mut v: Vec| { + while let Some(&0) = v.last() { + v.pop(); + } + v + }; + + let classical = |p1: &[u32], p2: &[u32]| -> Vec> { + let mk = |p: &[u32]| { + let mut m = MilnorBasisElement { + q_part: 0, + p_part: p.to_vec(), + degree: 0, + }; + m.compute_degree(TWO); + m + }; + let (m1, m2) = (mk(p1), mk(p2)); + let deg = m1.degree + m2.degree; + alg.compute_basis(deg); + let mut res = FpVector::new(TWO, alg.dimension(deg)); + alg.multiply(res.as_slice_mut(), 1, &m1, &m2); + let mut out: Vec> = res + .iter_nonzero() + .map(|(idx, _)| trim(alg.basis_element_from_index(deg, idx).p_part.clone())) + .collect(); + out.sort(); + out + }; + + let motivic = |p1: &[u32], p2: &[u32]| -> Vec> { + let to_paper = |pp: &[u32]| { + let mut r = vec![0u32]; + r.extend_from_slice(pp); + trim(r) + }; + let mut out: Vec> = multiply(&(0, to_paper(p1)), &(0, to_paper(p2))) + .into_iter() + .filter(|((e, _), _)| *e == 0) // keep the pure-ξ outputs + .map(|((_, r), c)| { + assert_eq!(c, Tau::one(), "pure-ξ output carried a τ power"); + trim(r.get(1..).map(<[u32]>::to_vec).unwrap_or_default()) + }) + .collect(); + out.sort(); + out + }; + + // The paper uses *conjugate* generators, so P(paper R) differs from the codebase's + // standard Milnor basis by the antipode χ for ξ_2 and higher (e.g. χ(ξ_2)=ξ_2+ξ_1^3). + // They agree on the self-conjugate ξ_1-power sector (ξ_1 is primitive), so restrict the + // comparison to outputs that are pure powers of ξ_1 (a `p_part` of length ≤ 1). + let xi1_only = |mut v: Vec>| { + v.retain(|pp| pp.len() <= 1); + v + }; + for (p1, p2) in [ + (vec![1u32], vec![1u32]), + (vec![2], vec![1]), + (vec![1], vec![2]), + (vec![3], vec![1]), + (vec![2, 1], vec![1]), + (vec![0, 1], vec![0, 1]), + (vec![2], vec![2]), + (vec![3], vec![3]), + ] { + assert_eq!( + xi1_only(motivic(&p1, &p2)), + xi1_only(classical(&p1, &p2)), + "ξ_1-sector mismatch for {p1:?} * {p2:?}" + ); + } + } + + #[test] + fn test_full_xi_matches_classical_via_antipode() { + // Full (all-ξ) cross-check against the codebase's classical Milnor product, using the + // antipode χ to reconcile the paper's conjugate generators with the codebase's standard + // ones. Identity checked (pure-ξ, τ-free, so over 𝔽₂): expressing both inputs and + // outputs of the motivic conjugate product in the standard basis via χ reproduces the + // codebase product. `A[μ][W] = coeff of conj-mon μ in std ξ(W) = [χ(mon_W)]_μ`. + use fp::{prime::TWO, vector::FpVector}; + + use crate::algebra::{ + Algebra, + milnor_algebra::{MilnorAlgebra, MilnorBasisElement}, + }; + + let alg = MilnorAlgebra::new(TWO, false); + let trim = |mut v: Vec| { + while let Some(&0) = v.last() { + v.pop(); + } + v + }; + let pp_to_paper = |pp: &[u32]| { + let mut r = vec![0u32]; + r.extend_from_slice(pp); + trim(r) + }; + let paper_to_pp = |r: &[u32]| trim(r.get(1..).map(<[u32]>::to_vec).unwrap_or_default()); + let mdeg = |pp: &[u32]| paper_bidegree(0, &pp_to_paper(pp)).0; + let xi_pps = |t: i32| -> Vec> { + enum_basis(t) + .into_iter() + .filter(|(e, _)| *e == 0) + .map(|(_, r)| paper_to_pp(&r)) + .collect() + }; + let a_coeff = |mu_pp: &[u32], w_pp: &[u32]| -> u32 { + let ap = antipode(&DualElement::from([((0u32, pp_to_paper(w_pp)), Tau::one())])); + // The antipode of a pure-ξ element is pure-ξ and τ-free, so a present monomial has + // coefficient τ^0. + ap.get(&(0u32, pp_to_paper(mu_pp))) + .map_or(0, |&c| u32::from(c == Tau::one())) + }; + let cb_mul = |a_pp: &[u32], b_pp: &[u32], acc: &mut BTreeMap, u32>| { + let mk = |p: &[u32]| { + let mut m = MilnorBasisElement { + q_part: 0, + p_part: p.to_vec(), + degree: 0, + }; + m.compute_degree(TWO); + m + }; + let (m1, m2) = (mk(a_pp), mk(b_pp)); + let deg = m1.degree + m2.degree; + alg.compute_basis(deg); + let mut res = FpVector::new(TWO, alg.dimension(deg)); + alg.multiply(res.as_slice_mut(), 1, &m1, &m2); + for (idx, _) in res.iter_nonzero() { + let pp = trim(alg.basis_element_from_index(deg, idx).p_part.clone()); + *acc.entry(pp).or_insert(0) ^= 1; + } + }; + + for (r1, r2) in [ + (vec![2u32], vec![1u32]), + (vec![1], vec![2]), + (vec![0, 1], vec![1]), + (vec![2], vec![2]), + (vec![3], vec![1]), + (vec![1, 1], vec![1]), + (vec![4], vec![1]), + (vec![2, 1], vec![2]), + ] { + // LHS: motivic conjugate product (E=0), converted to std via χ. + let conj: Vec> = multiply(&(0, pp_to_paper(&r1)), &(0, pp_to_paper(&r2))) + .into_iter() + .filter(|((e, _), _)| *e == 0) + .map(|((_, r), _)| paper_to_pp(&r)) + .collect(); + let mut lhs: BTreeMap, u32> = BTreeMap::new(); + for w in xi_pps(mdeg(&r1) + mdeg(&r2)) { + let s = conj.iter().fold(0u32, |acc, mu| acc ^ a_coeff(mu, &w)); + if s != 0 { + lhs.insert(w, s); + } + } + + // RHS: convert inputs to std via χ, multiply in the codebase. + let in1: Vec> = xi_pps(mdeg(&r1)) + .into_iter() + .filter(|w| a_coeff(&r1, w) == 1) + .collect(); + let in2: Vec> = xi_pps(mdeg(&r2)) + .into_iter() + .filter(|w| a_coeff(&r2, w) == 1) + .collect(); + let mut rhs: BTreeMap, u32> = BTreeMap::new(); + for w1 in &in1 { + for w2 in &in2 { + cb_mul(w1, w2, &mut rhs); + } + } + rhs.retain(|_, v| *v != 0); + + assert_eq!(lhs, rhs, "full ξ cross-check failed for {r1:?} * {r2:?}"); + } + } + + #[test] + fn test_product_associative() { + // Associativity (a·b)·c = a·(b·c) is convention-independent and a strong global check. + // Extend the basis×basis product to element×basis, F_2[τ]-linearly (add_term handles the + // coefficient sum + zero cancellation via Tau arithmetic). + fn mul_elt_basis(x: &DualElement, c: &(u32, Vec)) -> DualElement { + let mut out = DualElement::new(); + for (z, &cz) in x { + for (w, cw) in multiply(z, c) { + add_term(&mut out, w, cz * cw); + } + } + out + } + + let basis: Vec<(u32, Vec)> = (0..=5).flat_map(enum_basis).collect(); + for a in &basis { + for b in &basis { + let ab = multiply(a, b); + for c in &basis { + let lhs = mul_elt_basis(&ab, c); + let bc = multiply(b, c); + // a · (b·c): multiply basis `a` on the left of each term of bc. + let mut rhs = DualElement::new(); + for (z, &cz) in &bc { + for (w, cw) in multiply(a, z) { + add_term(&mut rhs, w, cz * cw); + } + } + assert_eq!(lhs, rhs, "associativity failed at {a:?},{b:?},{c:?}"); + } + } + } + } + + #[test] + fn test_product_q0_p_xi1() { + // The case my reading of Kong–Lin Theorem 5.1 got wrong. By duality: + // Q_0 · P(ξ_1) = Q_1 + Q_0 P(ξ_1) (both with coefficient τ^0). + let q0 = (0b1, vec![]); + let p_xi1 = (0, vec![0, 1]); + assert_eq!( + multiply(&q0, &p_xi1), + DualElement::from([ + ((0b10, vec![]), Tau::one()), // Q_1 + ((0b1, vec![0, 1]), Tau::one()), // Q_0 P(ξ_1) + ]) + ); + } + + #[test] + fn test_product_unit_and_squares() { + // 1 · x = x. + let x = (0b101, vec![0, 2]); + assert_eq!( + multiply(&(0, vec![]), &x), + DualElement::from([(x.clone(), Tau::one())]) + ); + // Q_i^2 = 0. (Q_0 = Sq^1 is the motivic Bockstein; P(ξ_1) = Sq^2 does NOT square to + // zero — its square is τ Q_0 Q_1 + …, a genuine motivic feature.) + for i in 0..3 { + assert!( + multiply(&(1 << i, vec![]), &(1 << i, vec![])).is_empty(), + "Q_{i}^2 ≠ 0" + ); + } + } + + #[test] + fn test_product_weight_homogeneous_and_tau_appears() { + // Every product is weight-homogeneous: term z occurs with a single power τ^k, and by + // weight-preservation of ψ, k = (w_a + w_b) - w_z ≥ 0 in *dual*-monomial weights (the + // algebra weight of a basis element is the negative of its dual monomial's weight, so + // τ, of weight -1, absorbs the difference). Also confirm τ genuinely enters at least + // one product (i.e. the τ_i^2 = τξ_{i+1} relation fires). + let basis: Vec<(u32, Vec)> = (0..=6).flat_map(enum_basis).collect(); + let mut saw_tau = false; + for a in &basis { + for b in &basis { + let (wa, wb) = (paper_bidegree(a.0, &a.1).1, paper_bidegree(b.0, &b.1).1); + for (z, c) in multiply(a, b) { + // Each coefficient is a single τ power by construction (a homogeneous Tau). + let k = c.valuation().expect("stored coefficients are never zero") as i32; + let wz = paper_bidegree(z.0, &z.1).1; + assert_eq!(wa + wb - wz, k, "weight mismatch: {a:?}*{b:?} → z={z:?}"); + assert!(k >= 0); + saw_tau |= k > 0; + } + } + } + assert!( + saw_tau, + "no product produced a τ coefficient — τ_i^2 relation never exercised" + ); + } + + #[test] + fn test_antipode() { + // χ(ξ_1) = ξ_1 (primitive); χ(ξ_2) = ξ_2 + ξ_1^3; χ(τ_0) = τ_0. + assert_eq!(antipode(&xi_gen(1)), xi_gen(1)); + assert_eq!( + antipode(&xi_gen(2)), + DualElement::from([ + ((0, vec![0, 0, 1]), Tau::one()), + ((0, vec![0, 3]), Tau::one()) + ]) + ); + assert_eq!(antipode(&tau_gen(0)), tau_gen(0)); + + // χ is an involution: χ(χ(x)) = x. + for x in [ + xi_gen(1), + xi_gen(2), + xi_gen(3), + tau_gen(0), + tau_gen(1), + tau_gen(2), + dual_mul(&tau_gen(0), &xi_gen(2)), + ] { + assert_eq!(antipode(&antipode(&x)), x, "χ not an involution on {x:?}"); + } + + // Antipode axiom: m(χ ⊗ id)ψ(x) = η ε(x); for positive-degree x this is 0. + let counit_via_antipode = |x: &DualElement| -> DualElement { + let mut out = DualElement::new(); + for ((l, r), c) in coproduct(x) { + let prod = dual_mul( + &antipode(&DualElement::from([(l, Tau::one())])), + &DualElement::from([(r, Tau::one())]), + ); + for (mon, cc) in prod { + add_term(&mut out, mon, c * cc); + } + } + out + }; + for x in [ + xi_gen(2), + xi_gen(3), + tau_gen(1), + tau_gen(2), + dual_mul(&tau_gen(0), &tau_gen(1)), + ] { + assert!( + counit_via_antipode(&x).is_empty(), + "antipode axiom failed on {x:?}" + ); + } + } + + #[test] + fn test_coproduct_generators() { + // ψ(τ_1) = 1⊗τ_1 + τ_0⊗ξ_1 + τ_1⊗1 (Kong–Lin §2.2). + assert_eq!( + coproduct(&tau_gen(1)), + TensorElement::from([ + (((0, vec![]), (0b10, vec![])), Tau::one()), // 1 ⊗ τ_1 + (((0b1, vec![]), (0, vec![0, 1])), Tau::one()), // τ_0 ⊗ ξ_1 + (((0b10, vec![]), (0, vec![])), Tau::one()), // τ_1 ⊗ 1 + ]) + ); + // ψ(ξ_2) = 1⊗ξ_2 + ξ_1⊗ξ_1^2 + ξ_2⊗1 (matches Milnor). + assert_eq!( + coproduct(&xi_gen(2)), + TensorElement::from([ + (((0, vec![]), (0, vec![0, 0, 1])), Tau::one()), // 1 ⊗ ξ_2 + (((0, vec![0, 1]), (0, vec![0, 2])), Tau::one()), // ξ_1 ⊗ ξ_1^2 + (((0, vec![0, 0, 1]), (0, vec![])), Tau::one()), // ξ_2 ⊗ 1 + ]) + ); + } + + #[test] + fn test_coproduct_is_algebra_map() { + // ψ(xy) = ψ(x)ψ(y); in particular this checks that the τ_i^2 = τξ_{i+1} reduction is + // compatible with the (τ-free) coproduct. + for (x, y) in [ + (tau_gen(0), tau_gen(0)), + (tau_gen(1), tau_gen(1)), + (tau_gen(0), tau_gen(1)), + (xi_gen(1), xi_gen(1)), + (tau_gen(0), xi_gen(1)), + (tau_gen(1), xi_gen(2)), + ] { + let lhs = coproduct(&dual_mul(&x, &y)); + let rhs = tensor_mul(&coproduct(&x), &coproduct(&y)); + assert_eq!(lhs, rhs, "ψ not an algebra map on {x:?} * {y:?}"); + } + } + + #[test] + fn test_rewrite_tau_square_free_is_identity() { + // A square-free τ(E) is already a basis element: τ_0 τ_2 → Q(E)P(0), τ^0. + let terms = rewrite_tau(&[1, 0, 1]); + assert_eq!( + terms, + vec![MotivicTerm { + tau_pow: 0, + e_mask: 0b101, + r: vec![], + }] + ); + } + + #[test] + fn test_rewrite_tau_example_3_1() { + // Kong–Lin Example 3.1, specialized to ρ = 0 (C-motivic): + // τ_0^2 τ_1 = τ · τ_1 ξ_1 + // Encoded: S = (2, 1); expected single term τ^1 · Q({1}) P((0,1)). + let terms = rewrite_tau(&[2, 1]); + assert_eq!( + terms, + vec![MotivicTerm { + tau_pow: 1, + e_mask: 0b10, // τ_1 + r: vec![0, 1], // ξ_1^1 + }], + "τ_0^2 τ_1 = τ τ_1 ξ_1 (ρ=0)" + ); + + // τ_0^4 = τ^2 ξ_1^2 + // Encoded: S = (4); expected single term τ^2 · Q(∅) P((0,2)). + let terms = rewrite_tau(&[4]); + assert_eq!( + terms, + vec![MotivicTerm { + tau_pow: 2, + e_mask: 0, + r: vec![0, 2], // ξ_1^2 + }], + "τ_0^4 = τ^2 ξ_1^2 (ρ=0)" + ); + } + + #[test] + fn test_closed_form_small_cases() { + // The case a naive reading of Theorem 5.1 gets wrong (the ξ₀ = 1 index + // absorption): Q_0 · P(ξ_1) = Q_1 + Q_0 P(ξ_1). + let q0 = (0b1, vec![]); + let p_xi1 = (0, vec![0, 1]); + assert_eq!(multiply_closed(&q0, &p_xi1), multiply(&q0, &p_xi1)); + + // The τ-generating case: P(ξ_1)² has a τ Q_0 Q_1 term. + assert_eq!( + multiply_closed(&p_xi1, &p_xi1), + multiply(&p_xi1, &p_xi1) + ); + + // Q_i² = 0. + for i in 0..3 { + assert_eq!( + multiply_closed(&(1 << i, vec![]), &(1 << i, vec![])), + multiply(&(1 << i, vec![]), &(1 << i, vec![])) + ); + } + } + + #[test] + fn test_closed_form_matches_duality_oracle() { + // Exhaustive fuzz: the closed form (Theorem 5.1) agrees with the duality + // oracle on every ordered pair of basis elements with deg(a)+deg(b) ≤ 18, + // including the τ-carrying and multi-Q cases. (The bound is kept modest + // because the *oracle* is slow; the closed form is not.) + let basis: Vec<(i32, (u32, Vec))> = (0..=14) + .flat_map(|t| enum_basis(t).into_iter().map(move |m| (t, m))) + .collect(); + for (ta, a) in &basis { + for (tb, b) in &basis { + if ta + tb > 18 { + continue; + } + assert_eq!( + multiply_closed(a, b), + multiply(a, b), + "closed form ≠ oracle for {a:?} · {b:?}" + ); + } + } + } + + #[test] + fn test_fill_block_matches_lazy_products() { + // The batch path (`fill_block`, which a GPU kernel would drive) must agree + // entry-for-entry with the lazy per-product path. Fill whole degree-pair + // blocks up front on one algebra, compute the same products lazily on + // another, and compare every structure constant. + let batched = MotivicMilnorAlgebra::new(); + let lazy = MotivicMilnorAlgebra::new(); + batched.compute_basis(12); + lazy.compute_basis(12); + for t1 in 0..=6 { + for t2 in 0..=6 { + batched.fill_block(t1, t2); + for idx1 in 0..batched.dimension(t1) { + for idx2 in 0..batched.dimension(t2) { + assert_eq!( + batched.product_indexed(t1, idx1, t2, idx2), + lazy.product_indexed(t1, idx1, t2, idx2), + "batch ≠ lazy at ({t1},{idx1})·({t2},{idx2})" + ); + } + } + } + } + } + + #[test] + fn test_rewrite_tau_unit() { + assert_eq!( + rewrite_tau(&[]), + vec![MotivicTerm { + tau_pow: 0, + e_mask: 0, + r: vec![], + }] + ); + } +} diff --git a/ext/crates/algebra/src/algebra/motivic/mod.rs b/ext/crates/algebra/src/algebra/motivic/mod.rs new file mode 100644 index 0000000000..04c4dc9949 --- /dev/null +++ b/ext/crates/algebra/src/algebra/motivic/mod.rs @@ -0,0 +1,28 @@ +//! The C-motivic (over $\mathbb{C}$, prime 2) Steenrod algebra and its mod-$\tau$ +//! reduction. +//! +//! This layer implements the *deformation* view of the C-motivic Adams $E_2$: the +//! C-motivic dual Steenrod algebra $A_C$ over $\mathbb{F}_2[\tau]$, its mod-$\tau$ +//! reduction $A_C/\tau$ (a connected finite-type $\mathbb{F}_2$-algebra), and the +//! coefficient ring $\mathbb{F}_2[\tau]$ itself. +//! +//! This module provides the foundation layer (see `MOTIVIC_PLAN.md`): +//! +//! - [`Tau`] — the coefficient ring $\mathbb{F}_2[\tau]$ as a small homogeneous +//! scalar ([`tau`]); the whole $\tau$-tower is carried here rather than threaded +//! through the resolution engine. +//! - [`MotivicMilnorAlgebra`] — $A_C$, a free $\mathbb{F}_2[\tau]$-module on the +//! Milnor basis, with the Kong–Lin product ([`milnor`]). This is the product +//! **engine**; it is deliberately not an [`Algebra`](crate::algebra::Algebra) +//! (that trait is over $\mathbb{F}_p$). +//! +//! The mod-$\tau$ reduction $A_C/\tau$ — the connected finite-type +//! $\mathbb{F}_2$-algebra that the existing resolution engine resolves to yield the +//! algebraic Novikov $E_2$ — is presented as an [`Algebra`](crate::algebra::Algebra) +//! in a follow-up on top of this engine. + +pub mod tau; +pub use tau::Tau; + +pub mod milnor; +pub use milnor::MotivicMilnorAlgebra; diff --git a/ext/crates/algebra/src/algebra/motivic/tau.rs b/ext/crates/algebra/src/algebra/motivic/tau.rs new file mode 100644 index 0000000000..685262de24 --- /dev/null +++ b/ext/crates/algebra/src/algebra/motivic/tau.rs @@ -0,0 +1,143 @@ +//! The coefficient ring $\mathbb{F}_2[\tau]$, represented by $\tau$-valuations. +//! +//! Every computation over the C-motivic Steenrod algebra that we perform is +//! **homogeneous** — differentials, products, and coproducts all preserve the +//! motivic weight — and a homogeneous element of $\mathbb{F}_2[\tau]$ in a fixed +//! bidegree is exactly $\tau^k$ for a single $k$ (or $0$). So we never need a +//! polynomial: a coefficient is one integer, its $\tau$-valuation, with a +//! sentinel for zero. That is [`Tau`]. +//! +//! This is a small, self-contained scalar type for the motivic layer only. It is +//! deliberately **not** an [`fp`](fp) replacement threaded through the resolution +//! engine: the engine stays $\mathbb{F}_p$-only (see the `motivic` module docs), +//! and `Tau` is used solely by the motivic algebra and the Phase 2 lift. + +/// A homogeneous $\mathbb{F}_2[\tau]$ coefficient: either $0$ or $\tau^k$ for a +/// unique $k \ge 0$, stored as the $\tau$-valuation with a zero sentinel. +/// +/// Over $\mathbb{F}_2$ a nonzero homogeneous coefficient has unit $\mathbb{F}_2$ +/// part, so the valuation is the whole story: `None` is $0$, `Some(k)` is +/// $\tau^k$. Multiplication adds valuations; addition of two coefficients in the +/// *same* bidegree can only ever combine equal powers (homogeneity), so it +/// cancels mod 2. +#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash, Default)] +pub struct Tau(Option); + +impl Tau { + /// The zero coefficient. + pub const ZERO: Self = Self(None); + + /// The unit $\tau^0 = 1$. + pub const ONE: Self = Self(Some(0)); + + /// The zero coefficient (method form, for symmetry with [`Tau::one`]). + pub fn zero() -> Self { + Self::ZERO + } + + /// The unit $1 = \tau^0$. + pub fn one() -> Self { + Self::ONE + } + + /// The coefficient $\tau^k$. + pub fn power(k: u32) -> Self { + Self(Some(k)) + } + + /// Whether this coefficient is zero. + pub fn is_zero(self) -> bool { + self.0.is_none() + } + + /// The $\tau$-valuation: `Some(k)` if this is $\tau^k$, `None` if zero. + pub fn valuation(self) -> Option { + self.0 + } + + /// The $n$-th power $(\tau^k)^n = \tau^{kn}$, with the convention $x^0 = 1$. + pub fn pow(self, n: u32) -> Self { + if n == 0 { + Self::ONE + } else { + Self(self.0.map(|k| k * n)) + } + } + + /// Multiply by $\tau^k$, i.e. raise the valuation by `k` (zero stays zero). + pub fn shift(self, k: u32) -> Self { + Self(self.0.map(|v| v + k)) + } +} + +impl std::ops::Mul for Tau { + type Output = Self; + + /// Product: $\tau^i \cdot \tau^j = \tau^{i+j}$; zero absorbs. + // Multiplying τ-powers genuinely *adds* their valuations, so the `+` here is correct. + #[allow(clippy::suspicious_arithmetic_impl)] + fn mul(self, other: Self) -> Self { + match (self.0, other.0) { + (Some(i), Some(j)) => Self(Some(i + j)), + _ => Self::ZERO, + } + } +} + +impl std::ops::Add for Tau { + type Output = Self; + + /// Sum of two homogeneous coefficients in one bidegree. By homogeneity the + /// two nonzero powers must be equal, so $\tau^k + \tau^k = 0$ over + /// $\mathbb{F}_2$; a zero summand is absorbed. Adding *distinct* nonzero + /// powers is inhomogeneous and cannot arise (checked in debug builds). + fn add(self, other: Self) -> Self { + match (self.0, other.0) { + (None, x) | (x, None) => Self(x), + (Some(i), Some(j)) => { + debug_assert_eq!(i, j, "adding inhomogeneous tau powers τ^{i} + τ^{j}"); + Self::ZERO + } + } + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn test_zero_one() { + assert!(Tau::zero().is_zero()); + assert!(!Tau::one().is_zero()); + assert_eq!(Tau::one().valuation(), Some(0)); + assert_eq!(Tau::power(3).valuation(), Some(3)); + assert_eq!(Tau::zero().valuation(), None); + } + + #[test] + fn test_mul_adds_valuations() { + assert_eq!(Tau::power(2) * Tau::power(3), Tau::power(5)); + assert_eq!(Tau::one() * Tau::power(4), Tau::power(4)); + assert_eq!(Tau::zero() * Tau::power(4), Tau::zero()); + assert_eq!(Tau::power(4) * Tau::zero(), Tau::zero()); + } + + #[test] + fn test_add_cancels_equal_powers() { + assert_eq!(Tau::power(2) + Tau::power(2), Tau::zero()); + assert_eq!(Tau::power(2) + Tau::zero(), Tau::power(2)); + assert_eq!(Tau::zero() + Tau::power(5), Tau::power(5)); + assert_eq!(Tau::zero() + Tau::zero(), Tau::zero()); + } + + #[test] + fn test_pow_and_shift() { + assert_eq!(Tau::power(3).pow(2), Tau::power(6)); + assert_eq!(Tau::power(3).pow(0), Tau::one()); + assert_eq!(Tau::zero().pow(0), Tau::one()); + assert_eq!(Tau::zero().pow(3), Tau::zero()); + assert_eq!(Tau::power(2).shift(3), Tau::power(5)); + assert_eq!(Tau::zero().shift(3), Tau::zero()); + } +} From cd83f678ed839bb5a147f15726068d8f39e36065 Mon Sep 17 00:00:00 2001 From: Claude Date: Sun, 12 Jul 2026 06:12:30 +0000 Subject: [PATCH 2/4] motivic: drop rustdoc links to the (follow-up) ctau module The two [`super::ctau`] intra-doc links referenced a module that isn't part of this foundation PR, breaking the docs build. Reword to plain prose. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013ePtYD7Bt4iPeCtmqtqvZE --- ext/crates/algebra/src/algebra/motivic/milnor.rs | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/ext/crates/algebra/src/algebra/motivic/milnor.rs b/ext/crates/algebra/src/algebra/motivic/milnor.rs index abd3f21326..a5fce5b290 100644 --- a/ext/crates/algebra/src/algebra/motivic/milnor.rs +++ b/ext/crates/algebra/src/algebra/motivic/milnor.rs @@ -32,8 +32,8 @@ //! Coefficients are homogeneous, so each one is a single power of $\tau$: a //! [`Tau`] scalar (see [`super::tau`]). This is the $A_C$ product **engine**, //! over $\mathbb{F}_2[\tau]$; the mod-$\tau$ reduction $A_C/\tau$ is presented to -//! the resolution engine as an ordinary $\mathbb{F}_2$-algebra in -//! [`super::ctau`]. +//! the resolution engine as an ordinary $\mathbb{F}_2$-algebra in a follow-up +//! layer built on top of this engine. //! //! Weight convention: the motivic weight of an algebra basis element is the //! *negative* of the weight of the dual monomial it pairs with, so that products @@ -892,8 +892,8 @@ pub fn multiply_closed(a: &(u32, Vec), b: &(u32, Vec)) -> DualElement /// exposed through [`product_indexed`](MotivicMilnorAlgebra::product_indexed), whose coefficients /// are [`Tau`] scalars. It is deliberately not an [`Algebra`](crate::algebra::Algebra) /// implementation, because that trait is over $\mathbb{F}_p$; the mod-$\tau$ reduction, which *is* -/// such an algebra, lives in [`super::ctau`], and the honest $\mathbb{F}_2[\tau]$ resolution is -/// built by lifting against this engine (Phase 2). +/// such an algebra, is a follow-up layer on top of this engine, and the honest +/// $\mathbb{F}_2[\tau]$ resolution is built by lifting against this engine (Phase 2). /// /// Weight convention: [`bidegree`](MotivicMilnorAlgebra::bidegree) returns `(t, w)` where `t` is /// the topological degree and `w` is the motivic weight in the presentation where $\tau$ has From d6b27bd5e7f7806730f3595c7b7b2620cf4fbaa5 Mon Sep 17 00:00:00 2001 From: Claude Date: Sun, 12 Jul 2026 06:32:29 +0000 Subject: [PATCH 3/4] motivic: drop the plan doc, trim the GPU handoff to durable notes MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Remove MOTIVIC_PLAN.md (a working design doc, not repo-durable) and cut the GPU handoff from a branch-specific handoff brief down to the two engine facts that actually inform a kernel port (F₂ atomic-XOR output with a host-side τ-power lookup; the mod-τ product is odd-primary-Milnor structured) plus where the batch boundary is (`fill_block` / `ProductBlock`). Drop the prose reference to the removed plan doc. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013ePtYD7Bt4iPeCtmqtqvZE --- MOTIVIC_PLAN.md | 275 ------------- ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md | 371 ++---------------- ext/crates/algebra/src/algebra/motivic/mod.rs | 2 +- 3 files changed, 43 insertions(+), 605 deletions(-) delete mode 100644 MOTIVIC_PLAN.md diff --git a/MOTIVIC_PLAN.md b/MOTIVIC_PLAN.md deleted file mode 100644 index cccc1ec22d..0000000000 --- a/MOTIVIC_PLAN.md +++ /dev/null @@ -1,275 +0,0 @@ -# C-motivic Ext: end-to-end implementation plan - -Computing the **C-motivic Adams E₂ page** — `Ext` over the C-motivic mod-2 -Steenrod algebra, `p = 2`, over `ℂ` — as a first-class computation that reuses -the existing `ext` engine without disturbing the classical path. - -> Base: this branch starts fresh from `master`. The earlier exploratory work -> (a generic-ground-ring engine + `FpTau` `Solvable`) is on -> `73de90ee4671068b3a018640e993ad1a61f0d272` and is recoverable for porting the -> validated *algebra* (`motivic_milnor.rs`), which is the only piece we reuse. - -> **Status (implemented).** Phases 0–4 are done. Phase 0: `algebra::motivic` -> (`Tau`, `MotivicMilnorAlgebra`, `CTauAlgebra`). Phase 1: -> `examples/resolve_motivic_ctau.rs` + `tests/motivic_ctau.rs`. Phases 2–3: -> `ext::motivic` (weight-homogeneous lift, `H(δ)`). Phase 4: -> `examples/resolve_motivic.rs` + the `examples/benchmarks/motivic-S_2` golden -> fixture. **All three anchors pass**: invert `τ` reproduces classical `Ext_A` -> rank-for-rank; `τ = 0` gives the algebraic Novikov ranks; keeping `τ` exhibits -> the `h₁`-tower `τ`-torsion (`h₁⁴`, `h₁⁵`). The classical path is untouched and -> the 69 resolution benchmarks stay bit-identical. One deviation from §2.5: the -> lift and `H(δ)` live in `ext::motivic` (they consume the resolution's -> quasi-inverses, populated by the `ext` engine) rather than the `algebra` crate; -> the pure `Tau` linear algebra could still be factored down later. - ---- - -## 0. Success criteria (the three anchors) - -The motivic E₂ is pinned by three views of one object, so a correct -implementation must reproduce all three: - -1. **invert τ → classical Adams E₂.** Regress against - `ext/examples/benchmarks/num_gens-S_2`. -2. **kill τ (τ = 0) → algebraic Novikov E₂.** The generator counts of the - resolution; equivalently the motivic Adams E₂ of `Cτ` (Gheorghe–Wang–Xu). -3. **keep τ → motivic Adams E₂ with τ-torsion.** The full bigraded answer, - including genuine τ-torsion classes (e.g. `h₁⁵` at `(s,t,w)`), regressed - against a published C-motivic chart (Isaksen–Wang–Xu). - -**Guardrail:** the classical path stays **bit-identical** — all 69 resolution -benchmarks unchanged. We achieve this by *not modifying the `ext` engine at all* -(see §2). - ---- - -## 1. The mathematics we are implementing - -Three algebras, not two. Let `A_C` be the C-motivic dual Steenrod algebra -(Voevodsky): - -``` -A_C = F₂[τ][ξ₁, ξ₂, …, τ₀, τ₁, …] / (τᵢ² = τ·ξᵢ₊₁), |τ| = (0, −1) (stem, weight) -``` - -- **`A_C`** — over `F₂[τ]`, bigraded by (stem, weight). Resolving the sphere - `F₂[τ]` over it gives the **motivic Adams E₂**. -- **`A_C/τ`** — set `τ = 0`: `A_C/τ ≅ F₂[ξᵢ] ⊗ Λ(τᵢ)`. This is a *connected, - finite-type `F₂`-algebra* (every generator has positive stem once τ is gone), - and it is **strictly bigger than the classical `A_* = F₂[ξᵢ]`**. Its Ext is - the **algebraic Novikov E₂**. -- **`A`** — the classical Steenrod algebra. A *further* collapse of `A_C` - (invert τ, or run the algebraic Novikov SS off `A_C/τ`). This is what `ext` - resolves today. - -### Why this makes motivic first-class *and* fast - -The minimal resolution over `A_C` is a genuine object over the motivic Steenrod -algebra. We compute it the way you compute over any τ-adically filtered ring — -**associated graded first, then deform** — which here means: - -- **Leading term:** resolve `A_C/τ` over the field `F₂`. This is a *motivic* - algebra (not classical `A`), resolved by the ordinary field engine. Fast. -- **Deform:** lift the differentials to `A_C` over `F₂[τ]`; the lifted `d∘d` is - zero only mod τ, and a small τ-power (Gaussian-elimination) correction - restores `d² = 0`. This is the honest `A_C` minimal resolution. - -Because the mod-τ object is `A_C/τ` (motivic), this is **not** "compute -classical and correct" — both objects are motivic. First-class semantics, and -the heavy combinatorics live over a field. (This is Guozhen Wang's Adams–Novikov -algorithm, module side; see §7.) - -### Representation: one `F₂[τ]`-matrix per stem, valuation entries - -A homogeneous element of `F₂[τ]` in a fixed bidegree is exactly `τᵏ` for a -unique `k` — store it as the integer `k` (`τ`-valuation), with a sentinel for 0. -Homogeneity fixes every differential entry: generator `j` at weight `wⱼ` maps to -generator `i` at weight `wᵢ` via `τ^{wᵢ − wⱼ}`. So the resolution differential at -each stem is **one matrix over `F₂[τ]`, the same size as the classical `F₂` -matrix**, with integer entries — *not* a sequence of matrices per weight. The -τ-tower is stored once. Linear algebra over it is `F₂` elimination with graded-PID -pivoting (pick the lowest-valuation entry in a column) plus integer bookkeeping — -a small constant over `F₂`, never `×(#weights)`, never polynomial arithmetic. - -### The E₂ is a cohomology, not a generator count - -`Ext_{A_C}(F₂[τ], F₂[τ]) = H(Hom_{A_C}(F_•, F₂[τ]), δ)`, where -`Hom_{A_C}(F_s, F₂[τ]) = Vₛ*` is a **free `F₂[τ]`-module** (rank = # generators) -and `δ` is the **identity-operation part** of the resolution differentials — -apply the augmentation `A_C → F₂[τ]`, which kills positive Steenrod operations -and leaves the τ-power on the unit. Every entry of `δ` is `τᵏ` (or 0), and -`δ² = 0`. The three anchors are three operations on this one complex: - -| operation on τ | result | -|---|---| -| invert τ | classical Adams E₂ | -| set τ = 0 | `δ ≡ 0` ⟹ generator counts = algebraic Novikov E₂ | -| keep τ, take `H(δ)` | motivic Adams E₂, with `F₂[τ]/τᵏ` torsion | - -The first (exploratory) implementation printed the middle row — generator counts -— and never applied `δ`. The real work is `H(δ)`, organised efficiently by the -**τ-Bockstein spectral sequence / Curtis table** (Wang §6): leading-term -cycle/tag table, `a` = permanent cycle, `a→b` = differential whose *length* is -the τ-filtration difference (= torsion order). This keeps the heavy pass over -`F₂` while producing honest `F₂[τ]`-cohomology. - ---- - -## 2. Architecture decisions (committed) - -1. **Do not genericize the `ext` engine.** Motivic factors through the *field* - engine via `A_C/τ`, so the engine stays `F₂`-only and classical stays - bit-identical by construction. (This supersedes / closes PR #263 — the - ground-ring genericization is unnecessary under this factoring.) -2. **Bigrading = internal degree `t` (engine) + weight carried alongside.** - `A_C/τ` is single-`t`-graded with weight a bounded per-basis label; `A_C` - carries weight in its `F₂[τ]` coefficients (the τ-tower). No ℤ²-graded engine. -3. **The deformation lives in `algebra`/the motivic layer, never in `ext`.** -4. **Module picture now; comodules deferred.** BP/Adams–Novikov would add a - *parallel* comodule engine later (§8) — additive, not a rewrite. -5. **Crate placement.** New code in the `algebra` crate (algebras, coefficient - rep, lift) plus, on the `ext` side, the cohomology (`H(δ)`) as a - generalization of the existing `ExtAlgebra` bridge (§6), and a thin - driver/example. The `ext` *resolution engine* and all classical modules are - untouched; `ExtAlgebra` is *extended*, with the classical `δ = 0` path - preserved. - ---- - -## 3. Phase 0 — Coefficient ring and the two algebras - -**Goal:** the `F₂[τ]` valuation ring and the algebras `A_C` (over `F₂[τ]`) and -`A_C/τ` (over `F₂`), with no resolution yet. - -- `algebra/src/algebra/motivic/tau.rs` — `Tau`: homogeneous `F₂[τ]` as a - `(coeff ∈ F₂, valuation ∈ u32)` with a zero sentinel. Ops: add (same - valuation), multiply (add valuations), `valuation`, `shift`, `pow`. Vector / - matrix helpers over `Tau` for the motivic layer only (small, self-contained; - **not** an `fp` replacement threaded through `ext`). -- `algebra/src/algebra/motivic/milnor.rs` — port `MotivicMilnorAlgebra` from - `73de90ee46` (validated: bigrading, Milnor product, coproduct, antipode, - cross-checked against classical). Strip its dependence on the old - `Ring`/`Solvable` tower; it needs only `Tau` arithmetic and the (stem, weight) - bigrading. -- `A_C/τ`: expose the mod-τ reduction as an ordinary `F₂` `Algebra` - (implements the existing `algebra::Algebra` trait), basis by internal degree - `t`, each basis element tagged with its weight. This is the object `ext` - resolves in Phase 1. - -**Tests (port + extend):** product vs classical Milnor after inverting τ; -associativity; coproduct is an algebra map; antipode; `A_C/τ` basis dimensions -per (t, weight) against hand computation for small t. - ---- - -## 4. Phase 1 — Resolve `A_C/τ` over `F₂` → algebraic Novikov E₂ - -**Goal:** first end-to-end numbers, using the **unmodified** `ext` engine. - -- Present `A_C/τ` to `ext` as an `Algebra`; resolve the trivial module `F₂` over - it with the existing minimal-resolution engine. Generators per (s, t) = - algebraic Novikov E₂ ranks; carry the weight tag per generator (needed by the - lift in Phase 2). -- `ext/examples/resolve_motivic_ctau.rs` — a `resolve`-style driver printing the - `A_C/τ` chart. - -**Milestone value:** a real, publishable computation (algebraic Novikov E₂ / -motivic Ext of `Cτ`) with essentially zero engine risk, and it validates the -`A_C/τ` algebra structure before the harder τ work. - -**Tests:** generator counts contain the classical ranks (each classical class -survives); spot-check low-degree algebraic Novikov ranks against Wang/known -values; `d² = 0` in the `F₂` resolution (engine already guarantees this). - ---- - -## 5. Phase 2 — Lift to `A_C` over `F₂[τ]` → the motivic minimal resolution - -**Goal:** the honest minimal resolution over `A_C`, in the valuation rep. - -- Take the Phase-1 `A_C/τ` resolution as the **model**: the generator skeleton - (ranks, weights, and the mod-τ differentials) is fixed. -- Lift each differential entry from `A_C/τ` to `A_C` (choose the `F₂[τ]` lift of - the mod-τ operation). Compute `d∘d`; it is `≡ 0 mod τ`, so it equals τ·(error). - Solve the small `F₂[τ]` linear system (graded-PID elimination on valuation - matrices) that adjusts the lift so `d² = 0` exactly — Wang §4–5. Optionally use - the "project to cogenerators" size reduction later as an optimization. -- Store differentials as per-stem `Tau` matrices (§1 representation). - -**Tests:** generator counts identical to Phase 1 (the lift preserves the -skeleton); `d² = 0` over `F₂[τ]` on every generator through a fixed range; -reducing the lifted differentials mod τ reproduces Phase 1 exactly. - ---- - -## 6. Phase 3 — `H(δ)` in `ExtAlgebra` → the motivic Adams E₂ - -`ExtAlgebra` (`ext/src/ext_algebra/`) is already the intended **resolution → Ext -bridge**: it wraps a resolution and presents `Ext(M, k)` as a bigraded algebra. -Today that bridge is *trivial* — `dimension(b) = number_of_gens_in_bidegree(b)` -— because over a field the Hom-complex differential vanishes (`δ = 0`), so -`Ext = generators`, read off. The motivic case is the **same bridge finally doing -the real work**: `Ext = H(δ)`, with the classical field case (`δ = 0`) recovering -today's behavior. Phase 3 is therefore a *generalization of `ExtAlgebra`*, not a -new extractor. - -**Goal:** teach `ExtAlgebra` that "Ext from a resolution" means taking cohomology, -then read the motivic E₂ off it. - -- Give `ExtAlgebra` a δ-aware Ext: extract `δ` (the identity-operation, τ-power - part of the resolution differentials — apply the augmentation `A_C → F₂[τ]`), a - complex of free `F₂[τ]`-modules with `δ² = 0`, and compute `Ext = H(δ)` as - f.g. `F₂[τ]`-modules (free ⊕ `F₂[τ]/τᵏ`). The current `dimension`/`basis` - become the `δ = 0` specialization. -- Organize the cohomology by the τ-Bockstein Curtis table (Wang §6): leading-term - cycle/tag, `diff_length` = τ-filtration difference = torsion order. This lands - next to the existing `secondary` (d₂) submodule — the τ-BSS / algebraic-Novikov - differentials are *more differentials on Ext*, the same family `ExtAlgebra` - already hosts — so products, d₂, and Massey products keep operating on the - motivic Ext through one surface. -- The two derived views fall out: τ = 0 (drop δ) = algebraic Novikov (Phase 1); - invert τ (rank of δ over `F₂(τ)`) = classical. - -**Tests — the three anchors:** -1. invert τ reproduces `num_gens-S_2` exactly, *including killing τ-torsion* - (`h₁⁵` disappears). -2. τ = 0 reproduces the Phase-1 algebraic Novikov ranks. -3. the full E₂ matches a published C-motivic chart over a fixed range: the - `h₁`-tower is τ-torsion-free/periodic as expected, and named τ-torsion classes - (`h₁⁵` etc.) appear with the correct `F₂[τ]/τᵏ` order. - ---- - -## 7. Phase 4 — Presentation, fixtures, docs - -- `ext/examples/resolve_motivic.rs` — print the bigraded motivic E₂ (ASCII, or - emit a SeqSee/JSON chart with weight and τ-torsion annotations). -- Golden fixtures for the motivic chart (small range) under - `ext/examples/benchmarks/` so the anchors regress in CI. -- Doc pass: a `motivic` module doc explaining the three algebras, the three - pages, and where the deformation lives; keep the model identifier out of any - committed text. - ---- - -## 8. Out of scope (future) - -- **BP / Adams–Novikov.** Wang's algorithm is natively about `BP∗BP`-comodules; - motivic is the `A_C/τ ↔ P` special case. Reaching BP means adding a *parallel* - comodule / Hopf-algebroid resolver (cofree resolutions, cobar `δ`) that shares - the coefficient layer (extended to p-adic `Z/pᵏ`) and the Curtis/SS extractor. - Additive; the module engine here is untouched. -- **McTague's Steenrod-algebra data structure.** A faster cycle/tag reduction - for the `F₂` skeleton (Phase 1). Slot it behind the leading-term-reduction - interface if/when Phase 1 becomes the bottleneck; not on the critical path. - ---- - -## 9. Reference - -- G. Wang, *Computations of the Adams–Novikov E₂-term* (2020) — the algorithm - (`MinimalResolution.pdf` in `joeybf/minimalresolution`; §4–5 model-and-lift, - §6 Curtis table). -- Gheorghe–Wang–Xu, *The special fiber of the motivic deformation…* - (arXiv:1809.09290) — algebraic Novikov E₂ = motivic Adams E₂ of `Cτ`. -- Isaksen–Wang–Xu — C-motivic charts for the anchor regressions. diff --git a/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md b/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md index 3d54b83ba4..ad78bc7622 100644 --- a/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md +++ b/ext/crates/algebra/MOTIVIC_GPU_HANDOFF.md @@ -1,329 +1,42 @@ -# Handoff: CubeCL kernel for the C-motivic Milnor multiply (ℂ, p = 2) - -A self-contained brief for building a GPU kernel that offloads the product in the -C-motivic Steenrod algebra `A_C` (and its mod-τ reduction `A_C/τ`). It mirrors the -classical Nassau handoff (`GPU_KERNEL_HANDOFF.md`, branch `nassau_gpu`) and reuses -as much of that machinery as the math allows. Written from a **CPU-only container -that cannot run GPU code**, so everything a GPU box needs — the exact CPU -reference to port, the output encoding, the validation oracle — is captured here. - -The motivic layer lives on branch **`claude/awesome-davinci-khkr8x`** of -`JoeyBF/sseq` (`ext/` workspace), in `ext/crates/algebra/src/algebra/motivic/`. -Pull that first. Read the classical `GPU_KERNEL_HANDOFF.md` alongside this — the -device plumbing (CubeCL toolchain, batching boundary, atomic-XOR output, hash-free -indexing, staging, gotchas) is shared and not repeated in full here. - ---- - -## 0. TL;DR — what to build, and the two facts that make it tractable - -Build a CubeCL kernel that computes, in one launch, the batch of `Q(E₁)P(R₁) · s` -products the motivic resolution issues at a differential-assembly step (the same -"fixed operation × general element" unit as classical; §3). Port the CPU -**two-matrix closed form** `multiply_closed` (Kong–Lin Theorem 5.1, ρ = 0) in -`motivic/milnor.rs`. Two facts collapse most of the apparent extra complexity: - -1. **The device output stays F₂ atomic-XOR — the τ never touches it.** Every - surviving `(X, Y)` pair contributes coefficient exactly `1`, and by weight- - homogeneity a fixed output monomial `z` carries a *single, determined* τ-power - `k(z)` — the weight difference between `a·b` and `z` (magnitude - `|w(a)+w(b)−w(z)|`; `bidegree(...).1` uses a negated-weight convention, so pin - the exact expression against `multiply_closed`, which the §7 oracle checks). So - the kernel computes **F₂ presence** of each output index (atomic XOR, identical - to classical); the τ-power is a host-side weight lookup, not part of the device - accumulator. **You do not need `Tau` arithmetic on the GPU.** - -2. **The resolution's products (the bottleneck) are the mod-τ product, which is - odd-primary-Milnor-structured at p = 2 — so the classical admissible-matrix - kernel ports almost directly.** `A_C/τ ≅ 𝔽₂[ξᵢ] ⊗ E(τᵢ)` has the odd-primary - dual's shape (with `2ⁱ` powers, not `pⁱ`), and mod τ the X-matrix constraint - collapses from `S(X) ≤ R₁` to `S(X) = R₁` exactly — **the classical admissible - matrix**. The only additions over the classical kernel are the second matrix - `Y` (the `τᵢ` / Q-part interaction) and `2ⁱ` weights. The *full* `A_C` product - (with τ, needed only by the Phase-2 lift) needs both matrices and the τ-power. - -So: start from the classical `nassau_gpu` kernel, keep its admissible-matrix X -enumeration and atomic-XOR output, add the Y enumeration and the weight table. -Validate bit-for-bit against `multiply_closed`. - ---- - -## 1. Environment reality - -Same as the classical handoff §1 — everything there applies verbatim (CubeCL -`0.10` + `wgpu` compiles headless; `cubecl-cpu` needs a network fetch that was -403-blocked in the sandbox but works on a normal box and gives a GPU-less -correctness backend; pick `wgpu`/`cuda`/`hip` per hardware; ignore the sandbox's -`-C target-cpu=native` SIGILL). **Nothing in this doc was executed on a GPU.** - ---- - -## 2. Where things live (branch `claude/awesome-davinci-khkr8x`) - -Key file: **`ext/crates/algebra/src/algebra/motivic/milnor.rs`** - -| Item | Symbol | Line (approx.) | -|---|---|---| -| **CPU reference to port** | `multiply_closed` (returns `DualElement`; the fuzz-validated product) | ~851 | -| X enumeration (ξ-matrix) | `Closed::enum_x` | ~702 | -| X→Y bridge (S′, RY targets, T(X)) | `Closed::on_x` | ~734 | -| Y enumeration (τ/Q-matrix) | `enum_y` | ~788 | -| Incremental b / T / S state | `struct Acc` | ~ (near `Closed`) | -| τ-rewrite coefficient `c(S, R)` | `c_coeff` | ~82 | -| Column enumerators | `columns_eq`, `columns_le`, `col0_x_options` | ~588, ~609, ~632 | -| Buffer bound (rows/anti-diagonals) | `const NB` (= 64) | ~656 | -| Product used by the engine | `MotivicMilnorAlgebra::product_indexed` | (returns `Vec<(Tau, usize)>`) | -| Basis index (binary search) | `index_of`, `index_of_ref` | ~945 | -| **Weight of a basis element** | `MotivicMilnorAlgebra::bidegree` (`.1` = weight) | ~950 | -| Product cache (replace on GPU) | `products: RwLock>` | struct field | -| CPU profiling counters | `PRODUCT_HITS/MISSES/NANOS` (env `MOT_PROFILE`) | ~52 | - -Mod-τ reduction (**the resolution algebra**): -**`ext/crates/algebra/src/algebra/motivic/ctau.rs`** - -| Item | Symbol | Line | -|---|---|---| -| `A_C/τ` product = τ⁰ part of `A_C` | `CTauAlgebra::multiply_basis_elements` | ~98 | -| Weight label per basis element | `CTauAlgebra::weight` | ~64 | -| Underlying `A_C` engine | `CTauAlgebra::engine` | ~55 | - -Coefficient ring (host-side only; **not** needed on device): `motivic/tau.rs` -(`Tau` = a τ-valuation, `Option`). - -Full-product consumer (Phase-2 lift): **`ext/src/motivic/mod.rs`** — -`MotivicResolution::lift`/`correct`/`composite` call `engine().product_indexed` -and *do* use the τ-powers. Driver + profiling: `ext/examples/resolve_motivic.rs`. - -Validation oracle (already in the crate): `multiply` (duality product) in -`motivic/milnor.rs`, and the test `test_closed_form_matches_duality_oracle`. - ---- - -## 3. The call chain the products flow through - -Two distinct consumers, different needs: - -``` -Phase 1 — RESOLUTION of A_C/τ (the bottleneck: ~97% of run time; MOD-τ products) - Resolution engine assembles a differential - └ module action → CTauAlgebra::multiply_basis_elements(out, coeff, R_deg, R_idx, s_deg, s_idx) - └ engine.product_indexed(...) → multiply_closed(a, b) [FULL A_C product] - └ keep only terms with τ-valuation 0 (drop τ-divisible) ← mod-τ - -Phase 2 — LIFT to A_C (correcting d²=0; FULL products with τ-powers) - MotivicResolution::correct → composite → engine.product_indexed(...) → multiply_closed -``` - -The fundamental kernel unit is the classical one: **a fixed operation -`a = Q(E₁)P(R₁)` × a general element `s` (a sum of `Q(E₂)P(R₂)` terms) → add -`a·s` to an F₂ output slice.** `a` repeats across many `s` in one assembly step — -that is the amortization a workgroup exploits (one workgroup per distinct `(E₁, -R₁)`; parallelize the term tests across `s`). - -> **Batching boundary — determine on the GPU box.** The classical kernel batches -> per Nassau `get_partial_matrix` (a proven-maximal, non-redundant unit; classical -> handoff §5). The motivic resolution currently uses the **standard** `Resolution` -> engine (`ext/src/motivic/mod.rs`, via `Resolution::new`), not Nassau, so the -> exact assembly boundary that maximizes same-`a` amortization must be located -> there (it is wherever a differential column is built from many algebra actions). -> Two options: (a) find/instrument the standard engine's equivalent of -> `get_partial_matrix`; or (b) port the motivic resolution to the Nassau -> signature algorithm and inherit its boundary and the classical §5 analysis -> wholesale. (b) is the cleaner long-term path and the one that reuses the most. - -> **CPU-side batch primitive already in place.** `MotivicMilnorAlgebra` now caches -> products as one dense [`ProductBlock`] per topological-degree pair `(t1, t2)` -> (`entry (idx1, idx2) = basis[t1][idx1] · basis[t2][idx2]`), with per-entry -> `OnceLock` fill. `MotivicMilnorAlgebra::fill_block(t1, t2)` computes an entire -> block in one (parallel) pass — this is the ready-made host handle for the GPU -> launch: hand the kernel a `(t1, t2)` block, get all `dim1 × dim2` structure -> constants back, store into the block's cells. -> -> **But do not fill whole blocks blindly.** Measured block *density* during a real -> resolution is only **~50%** (stem 30 and 40 both ≈ 50%): the standard engine -> requests roughly half of each touched block's `(idx1, idx2)` grid. So a GPU (or -> CPU) pass that fills the *entire* block does ~2× the necessary -> `multiply_closed` work — and since product compute is ~95% of resolution time, -> eager full-block fill is a **net CPU loss** (confirmed by measurement; that is -> why the CPU path stays lazy-per-entry). For the GPU this means either (i) pass -> the kernel only the *requested* `(idx1, idx2)` sublist of the block (keep the -> engine's demand-driven set), accepting a ragged batch; or (ii) fill the full -> block and eat the 2× overcompute, betting GPU throughput makes it a wash. Decide -> with an A/B on the GPU box; the ragged-batch path (i) is the safer default and -> mirrors what the lazy CPU cache already does. - ---- - -## 4. The math to port (`multiply_closed`, ρ = 0) - -`Q(E₁)P(R₁) · Q(E₂)P(R₂)` is a sum over **two** non-negative integer matrices -`X` (the ξ-part) and `Y` (the τ/Q-part), Kong–Lin Theorem 5.1 at ρ = 0. Sequences -are indexed from 0; **ξ-type** sequences (`R₁, R₂, S(X), R(X), T(X)`) drop index 0 -(ξ₀ = 1), **τ-type** sequences (`E, S(Y), T(Y)`) keep it (τ₀ ≠ 1). This ξ₀ = 1 / τ₀ -≠ 1 asymmetry is the whole subtlety — get it wrong and terms silently vanish. - -Matrix functionals (per matrix `M`, entries `x_{i,j}`, i = row, j = column): -- `S(M)_r = Σⱼ x_{r,j}` (row sum), `R(M)_j = Σᵢ 2ⁱ x_{i,j}` (column weighted sum), - `T(M)_r = Σ_{i+j=r} x_{i,j}` (anti-diagonal sum). -- `b(M) = Πᵣ multinomial(anti-diagonal r) mod 2` = **1 iff the entries on each - anti-diagonal are bitwise disjoint** (no carries). The CPU keeps a running XOR - per anti-diagonal and **rejects a column the instant an entry collides**, so - only `b = 1` matrices are ever visited and `b`/`T`/`S` are never recomputed - (`struct Acc`, `enum_x`/`enum_y`). Port this pruning — it is what makes the - enumeration cheap and it maps to per-anti-diagonal bit masks exactly like the - classical `masks[j]`. - -Enumeration (see `enum_x` → `on_x` → `enum_y`): -- **X** (`enum_x`): column j ≥ 1 ranges over `Σᵢ 2ⁱ x_{i,j} ≤ R₂[j]`; column 0 over - `x_{i,0} ≤ R₁[i]`; prune on the row budget `S(X)_r ≤ R₁[r]` (r ≥ 1) and on - b-collision. **For the mod-τ product the row budget is an equality - `S(X) = R₁` — i.e. the classical admissible matrix of `Sq(R₁)`.** -- **`on_x`**: `S′ = R₁ − S(X)` (ξ-part, index 0 = 0); τ-power `k = Σ(S′)`; the Y - column targets `R(Y)[j] = R₂[j] − R(X)[j]` (j ≥ 1) and the forced - `R(Y)[0] = Σ₂(E₁) + Σ₂(S′) − Σ_{j≥1} R(Y)[j]`; output P-part `T(X)` (drop index 0). -- **Y** (`enum_y`): each column an exact weighted sum `R(Y)[j]`, pruned on - b-collision. At a leaf the pair contributes iff the ρ = 0 degree equation - `popcount(E₁) + 2Σ(S′) = Σ(S(Y))` holds **and** `c(S(Y), S′) ≠ 0` (the τ-rewrite - coefficient `c_coeff`, Theorem 3.4 — pure integer/`binomial2` arithmetic, a - direct GPU port). Output monomial: `Q(E₂ + T(Y)) P(T(X))`, and it must be - square-free (`E₂ + T(Y) ∈ Seq₁`). - -**Coefficient and τ (the §0 fact 1 in detail).** Each surviving `(X, Y)` -contributes `b(X)·b(Y)·c(S(Y),S′)·binom(E₂+T(Y), E₂) = 1` — every factor is 0/1 and -we only reach the leaf when all are 1. So the coefficient of an output monomial is -`#{surviving (X,Y) landing on it} mod 2` — **an F₂ count = atomic XOR**, exactly -classical. The τ-power `k = Σ(S′)` is the *same* for every `(X,Y)` landing on a -given `z` (weight-homogeneity; `k = Σ(S′)`), so recover it host-side from the -weight table (`bidegree(...).1`, sign calibrated against `multiply_closed`) rather -than tracking it on device. **Mod-τ (resolution): keep only `z` with `k(z) = 0`** -— equivalently only enumerate `X` with `S(X) = R₁`. - -**Special case:** `a = 1` (E₁ = 0, R₁ = ∅) ⇒ output is `s` unchanged. Common -(identity operation). - ---- - -## 5. Batching / occupancy - -The classical §5 analysis (batch per assembly launch, one workgroup per distinct -operation, ~60% of term-work fills a 256-wide workgroup on same-operation -amortization at high stems, "widen the scope buys nothing") **carries over** — -the kernel unit is identical. Re-measure for the motivic workload once the -boundary of §3 is fixed; the motivic algebra is ~2× larger and the algebraic- -Novikov E₂ has more generators, so expect *denser* launches than classical at the -same stem (more `s` terms per operation), which only helps occupancy. - -CPU baseline to beat (this sandbox, `--features concurrent`, ~3–4 cores; use for -harness sanity, not as an absolute target): report `(40, 8)` resolution ≈ 2.1 s, -`(50, 10)` ≈ 18 s; the closed-form product is ~62 µs/call at stem 40, ~157 µs at -stem 50, and is ~100% of resolution time (linear algebra is ~3% — **the product is -the bottleneck, exactly as classical**). Regenerate with: - -``` -MOT_PROFILE=1 cargo run --release --no-default-features --features concurrent \ - --example resolve_motivic # prints resolution / products / lift / cohomology -``` - -`PRODUCT_HITS/MISSES/NANOS` (in `motivic/milnor.rs`) report the closed-form call -count and time; the `MOT_PROFILE` env gate is read in `MotivicResolution::new`. - ---- - -## 6. Data layout for the device - -Upload once per algebra (rebuild as `compute_basis` grows): - -- **Basis table + index.** The `A_C/τ` (= `A_C`) basis in a degree is a contiguous - sorted array of monomials `(e_mask: u32, R: [u32; MAX_XI])`, sorted by - `(e_mask, R)`. `index_of_ref` is a **binary search** over it — directly GPU- - portable (no hashing), unlike the classical `HashMap`. That is the simplest - device index and is enough to start. A hash-free *`seqno`-style* rank (classical - handoff §6) is a later optimization: the motivic rank factors as - `(rank of E among Seq₁ subsets at that weight) × (classical-seqno rank of R)` — - more arithmetic than classical but still hash-free; only build it if the binary - search shows up in the profile. -- **Weight table.** One `i32` per basis element = `bidegree(degree, idx).1`. The - host multiplies each present output index `z` by `τ^{k(z)}`, `k(z)` the weight - difference (magnitude `|w(a)+w(b)−w(z)|`; fix the sign once against - `multiply_closed`). For the mod-τ (resolution) kernel, keep only the `k(z) = 0` - slice (`w(z) = w(a)+w(b)`) and the output is pure F₂ — no weight multiply at all. -- **Per-launch inputs.** For each product: the operation `(E₁, R₁)` (fixed within a - workgroup) and `s` as a list of `(E₂, R₂)` terms (or indices into the shared - basis table), each padded to `MAX_XI` `u32`s for coalesced access. -- **Output.** F₂ vectors bit-packed into `Atomic` arrays, one per matrix row; - `atomic_xor(word, 1 << bit)` at `idx = index_of_ref(z)`. Match the crate's `fp` - p = 2 bit order (LSB-first within a limb — verify against `fp`). Multiple - `(X, Y)` pairs, possibly across threads, land on the same index and must cancel - mod 2 — **atomic XOR, not add** (classical gotcha §10 applies identically). - ---- - -## 7. Validation & benchmarking (on the GPU box) - -**Correctness oracle — already in the crate.** -`test_closed_form_matches_duality_oracle` (`motivic/milnor.rs` tests) checks -`multiply_closed` bit-for-bit against the independent duality product `multiply` -for **every** basis pair up to `deg(a)+deg(b) ≤ 18`, including the τ-carrying -(`Sq²·Sq² ∋ τQ₀Q₁`) and multi-Q cases. **Mirror it for the kernel:** for every -`(a, s)` in a degree range, compare the kernel's F₂ output slice (per output -index, i.e. the presence bits) against `multiply_closed(a, bⱼ)` reduced to -presence, and separately check the host-rebuilt τ-powers against -`multiply_closed`'s `Tau` coefficients. With the `cubecl-cpu` runtime this runs -without a GPU. Add a mod-τ variant comparing against -`CTauAlgebra::multiply_basis_elements`. - -**End-to-end A/B.** CPU baseline is the current default (`multiply_closed` behind -the `RwLock` cache). Feature-flag the GPU path in the resolution's product call, -keep device buffers resident across the resolution (upload only the degrees that -grew), and A/B `resolve_motivic` wall time at e.g. `(50,10)` and `(60,10)`. The -GPU must beat the CPU path *including* host↔device transfer. - ---- - -## 8. Suggested staging (each stage compiles + validates before the next) - -1. **Toolchain slice.** `cubecl` behind a `gpu` feature on `algebra` (classical - handoff §9.1 skeleton — a trivial F₂-XOR `#[cube]` kernel + host launch + test). -2. **Index on device.** Port `index_of_ref` as a `#[cube]` binary search over the - uploaded sorted basis; test it is the identity `index_of_ref(basisₖ) == k` over - a full degree. -3. **Single-`(E₁,R₁)` mod-τ multiply.** One workgroup: enumerate the classical - admissible matrices of `Sq(R₁)` (`S(X) = R₁`) — reuse the classical - `AdmissibleMatrix` port from `nassau_gpu` — then the Y (Q-part) enumeration and - the `c_coeff`/ρ leaf test, atomic-XOR the F₂ output. Validate one `(a, s)` - against `CTauAlgebra::multiply_basis_elements`. -4. **Full `A_C` multiply.** Relax `S(X) = R₁` to `S(X) ≤ R₁`, emit at all - τ-powers, rebuild `Tau` on the host from the weight table. Validate against - `multiply_closed` (§7 oracle). -5. **Batched launch + wire into the resolution + A/B** (§3 boundary, §5 occupancy, - §7 A/B). Keep buffers resident. - ---- - -## 9. Gotchas (motivic-specific; the classical §10 all still apply) - -- **ξ₀ = 1 but τ₀ ≠ 1.** ξ-sequences compare on indices ≥ 1, τ-sequences on ≥ 0. - This is the single most error-prone point — the naive Theorem 5.1 reading drops - the `Q₀·P(ξ₁) → Q₀P(ξ₁)` term because it treats `S(X)`'s index 0 as constrained. - The CPU reference and its oracle test pin it; keep the test. -- **τ never goes on the device.** Coefficients are F₂ presence; the τ-power is a - weight-table lookup on the host. Don't build `Tau` arithmetic into the kernel. -- **Two matrices, but Y is small.** X is the classical admissible matrix; Y - handles only the `τᵢ` interaction and is typically tiny (often trivial when - E₁ = E₂ = 0 and the product is pure-ξ). Enumerate X per workgroup (as classical); - the Y loop is a short inner sweep per `(X, term)`. -- **Mod-τ vs full.** The resolution wants `S(X) = R₁` (τ⁰) — a *restriction* of the - full enumeration, so a mod-τ-only kernel is strictly smaller. Ship it first (it - is the bottleneck); the full product is only for the Phase-2 lift. -- **Square-free output.** `E₂ + T(Y)` must be ∈ Seq₁ (each τᵢ at most once); a - collision kills the term. Q indices fit a `u32` mask. -- **`c_coeff` needs `r₀ = 0`.** `S′` has index 0 = 0 by construction; assert it. -- **Trailing-zero canonicalization.** Output P-parts are trimmed before indexing; - the sorted basis assumes trimmed monomials (same as classical §10). - ---- - -*Written from a sandbox without GPU access. The math, the CPU reference, and the -oracle are exercised on CPU (all `motivic` tests green); nothing here ran on a -GPU. Validate stage-by-stage on real hardware, starting from the classical -`nassau_gpu` admissible-matrix kernel — most of it is reused.* +# GPU port notes: the C-motivic Milnor multiply (ℂ, p = 2) + +The C-motivic product (`multiply_closed`, Kong–Lin Theorem 5.1 at ρ = 0, in +`motivic/milnor.rs`) is the arithmetic bottleneck of the deformation pipeline and +the natural GPU-offload target. This note records the two facts that make a kernel +tractable and where the batch boundary is; the device plumbing (CubeCL toolchain, +atomic-XOR output, hash-free indexing) mirrors the classical Nassau kernel and is +not repeated here. + +## The batch unit + +`MotivicMilnorAlgebra::fill_block(t1, t2)` computes an entire `ProductBlock` — every +structure constant `Q(E₁)P(R₁) · Q(E₂)P(R₂)` for a pair of topological degrees — in +one shot, filling independent `OnceLock` cells in parallel. That block is the shape a +GPU kernel would fill in a single launch: the resolution asks for the same +degree-pair block repeatedly, so the per-degree-pair block is the batching boundary. + +## Two facts that collapse the apparent complexity + +1. **The device output stays F₂ atomic-XOR — τ never touches it.** Every surviving + matrix pair `(X, Y)` contributes coefficient exactly `1`, and by weight-homogeneity + a fixed output monomial `z` carries a single, determined τ-power — the weight + difference between `a·b` and `z`. So the kernel computes only **F₂ presence** of + each output index (identical to the classical kernel); the τ-power is a host-side + weight lookup, not part of the device accumulator. No `Tau` arithmetic on the GPU. + (Pin the exact τ-power expression against `multiply_closed` — `bidegree(...)` uses a + negated-weight convention.) + +2. **The bottleneck products are the mod-τ product, which is odd-primary-Milnor + structured at p = 2.** `A_C/τ ≅ 𝔽₂[ξᵢ] ⊗ E(τᵢ)` has the odd-primary dual's shape + (with `2ⁱ` powers, not `pⁱ`), and mod τ the X-matrix constraint collapses from + `S(X) ≤ R₁` to `S(X) = R₁` — exactly the classical admissible matrix. The only + additions over the classical kernel are the second matrix `Y` (the `τᵢ`/Q-part + interaction) and the `2ⁱ` weights. The full `A_C` product (with τ, needed only by + the Phase-2 lift) needs both matrices plus the τ-power. + +## Validation + +Validate any kernel bit-for-bit against `multiply_closed`, which is itself checked +exhaustively against the duality oracle `multiply` in the crate's tests. Start from +the classical admissible-matrix X enumeration, add the Y enumeration and the weight +table. diff --git a/ext/crates/algebra/src/algebra/motivic/mod.rs b/ext/crates/algebra/src/algebra/motivic/mod.rs index 04c4dc9949..94f1fa1bd3 100644 --- a/ext/crates/algebra/src/algebra/motivic/mod.rs +++ b/ext/crates/algebra/src/algebra/motivic/mod.rs @@ -6,7 +6,7 @@ //! reduction $A_C/\tau$ (a connected finite-type $\mathbb{F}_2$-algebra), and the //! coefficient ring $\mathbb{F}_2[\tau]$ itself. //! -//! This module provides the foundation layer (see `MOTIVIC_PLAN.md`): +//! This module provides the foundation layer for that computation: //! //! - [`Tau`] — the coefficient ring $\mathbb{F}_2[\tau]$ as a small homogeneous //! scalar ([`tau`]); the whole $\tau$-tower is carried here rather than threaded From 032123a133f76eeba82467f7fbcac1d3354cf8ed Mon Sep 17 00:00:00 2001 From: Claude Date: Sun, 12 Jul 2026 07:34:32 +0000 Subject: [PATCH 4/4] =?UTF-8?q?motivic:=20fix=20CI=20=E2=80=94=20redundant?= =?UTF-8?q?=20rustdoc=20link=20+=20nightly=20rustfmt?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit - `tau.rs`: `[`fp`](fp)` is a redundant explicit link (the label `fp` already resolves to the crate), which `-D rustdoc::redundant-explicit-links` (implied by the docs job's `-D warnings`) rejects. Drop the explicit target. - Run nightly `cargo fmt` for the `lint` job: the repo's rustfmt.toml enables unstable options (reorder_impl_items, group_imports, format_strings) that stable rustfmt silently skips, so the motivic files were never formatted the way CI checks. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013ePtYD7Bt4iPeCtmqtqvZE --- .../algebra/src/algebra/motivic/milnor.rs | 77 +++++++++++++------ ext/crates/algebra/src/algebra/motivic/tau.rs | 7 +- 2 files changed, 58 insertions(+), 26 deletions(-) diff --git a/ext/crates/algebra/src/algebra/motivic/milnor.rs b/ext/crates/algebra/src/algebra/motivic/milnor.rs index a5fce5b290..cf3d792261 100644 --- a/ext/crates/algebra/src/algebra/motivic/milnor.rs +++ b/ext/crates/algebra/src/algebra/motivic/milnor.rs @@ -621,11 +621,23 @@ thread_local! { } fn columns_eq_cached(target: u32) -> Rc>> { - COLS_EQ.with(|c| Rc::clone(c.borrow_mut().entry(target).or_insert_with(|| Rc::new(columns_eq(target))))) + COLS_EQ.with(|c| { + Rc::clone( + c.borrow_mut() + .entry(target) + .or_insert_with(|| Rc::new(columns_eq(target))), + ) + }) } fn columns_le_cached(bound: u32) -> Rc>> { - COLS_LE.with(|c| Rc::clone(c.borrow_mut().entry(bound).or_insert_with(|| Rc::new(columns_le(bound))))) + COLS_LE.with(|c| { + Rc::clone( + c.borrow_mut() + .entry(bound) + .or_insert_with(|| Rc::new(columns_le(bound))), + ) + }) } /// The column-0 options for the `X` matrix: `x_{0,0} = 0`, and `x_{i,0} ≤ R₁[i]` @@ -737,7 +749,11 @@ impl<'a> Closed<'a> { let sprime: Vec = (0..self.l) .map(|j| if j == 0 { 0 } else { self.r1v[j] - acc.rows[j] }) .collect(); - let sigma2_sprime: i64 = sprime.iter().enumerate().map(|(j, &v)| (v as i64) << j).sum(); + let sigma2_sprime: i64 = sprime + .iter() + .enumerate() + .map(|(j, &v)| (v as i64) << j) + .sum(); let sigma_sprime: u32 = sprime.iter().sum(); // RY targets: column j ≥ 1 of Y has weighted sum R₂[j] − R(X)[j]. @@ -757,7 +773,11 @@ impl<'a> Closed<'a> { } // The output P-part T(X) (index 0 dropped: ξ₀ = 1). - let out_r = trim(std::iter::once(0).chain((1..NB).map(|d| acc.sum[d])).collect()); + let out_r = trim( + std::iter::once(0) + .chain((1..NB).map(|d| acc.sum[d])) + .collect(), + ); let mut y_cands: Vec>>> = Vec::with_capacity(self.l); y_cands.push(columns_eq_cached(ry0 as u32)); @@ -807,7 +827,10 @@ fn enum_y<'b>( if e1_mask.count_ones() + 2 * sigma_sprime != sigma_sy { return; } - let sy_len = (0..NB).rev().find(|&i| acc.rows[i] > 0).map_or(0, |i| i + 1); + let sy_len = (0..NB) + .rev() + .find(|&i| acc.rows[i] > 0) + .map_or(0, |i| i + 1); if c_coeff(&acc.rows[..sy_len], sprime) == 0 { return; } @@ -826,7 +849,11 @@ fn enum_y<'b>( return; } for cand in y_cands[j].iter() { - if !cand.iter().enumerate().all(|(i, &v)| acc.or[i + j] & v == 0) { + if !cand + .iter() + .enumerate() + .all(|(i, &v)| acc.or[i + j] & v == 0) + { continue; } for (i, &v) in cand.iter().enumerate() { @@ -836,7 +863,16 @@ fn enum_y<'b>( } ycols.push(cand); enum_y( - y_cands, j + 1, ycols, acc, e1_mask, e2_mask, sprime, sigma_sprime, out_r, out, + y_cands, + j + 1, + ycols, + acc, + e1_mask, + e2_mask, + sprime, + sigma_sprime, + out_r, + out, ); ycols.pop(); for (i, &v) in cand.iter().enumerate() { @@ -1094,13 +1130,7 @@ impl MotivicMilnorAlgebra { /// elements in degree `t1 + t2`: a list of `(coefficient, index)` pairs. The product of two /// homogeneous basis elements is weight-homogeneous, so each coefficient is a single power of /// $\tau$ — a [`Tau`] scalar. - pub fn product_indexed( - &self, - t1: i32, - idx1: usize, - t2: i32, - idx2: usize, - ) -> Vec<(Tau, usize)> { + pub fn product_indexed(&self, t1: i32, idx1: usize, t2: i32, idx2: usize) -> Vec<(Tau, usize)> { self.product_indexed_with(t1, idx1, t2, idx2, <[_]>::to_vec) } @@ -1196,7 +1226,10 @@ mod tests { .collect(); assert_eq!( terms, - DualElement::from([((0b10, vec![]), Tau::one()), ((0b1, vec![0, 1]), Tau::one())]) + DualElement::from([ + ((0b10, vec![]), Tau::one()), + ((0b1, vec![0, 1]), Tau::one()) + ]) ); } @@ -1380,7 +1413,10 @@ mod tests { .collect() }; let a_coeff = |mu_pp: &[u32], w_pp: &[u32]| -> u32 { - let ap = antipode(&DualElement::from([((0u32, pp_to_paper(w_pp)), Tau::one())])); + let ap = antipode(&DualElement::from([( + (0u32, pp_to_paper(w_pp)), + Tau::one(), + )])); // The antipode of a pure-ξ element is pure-ξ and τ-free, so a present monomial has // coefficient τ^0. ap.get(&(0u32, pp_to_paper(mu_pp))) @@ -1617,9 +1653,9 @@ mod tests { assert_eq!( coproduct(&xi_gen(2)), TensorElement::from([ - (((0, vec![]), (0, vec![0, 0, 1])), Tau::one()), // 1 ⊗ ξ_2 + (((0, vec![]), (0, vec![0, 0, 1])), Tau::one()), // 1 ⊗ ξ_2 (((0, vec![0, 1]), (0, vec![0, 2])), Tau::one()), // ξ_1 ⊗ ξ_1^2 - (((0, vec![0, 0, 1]), (0, vec![])), Tau::one()), // ξ_2 ⊗ 1 + (((0, vec![0, 0, 1]), (0, vec![])), Tau::one()), // ξ_2 ⊗ 1 ]) ); } @@ -1695,10 +1731,7 @@ mod tests { assert_eq!(multiply_closed(&q0, &p_xi1), multiply(&q0, &p_xi1)); // The τ-generating case: P(ξ_1)² has a τ Q_0 Q_1 term. - assert_eq!( - multiply_closed(&p_xi1, &p_xi1), - multiply(&p_xi1, &p_xi1) - ); + assert_eq!(multiply_closed(&p_xi1, &p_xi1), multiply(&p_xi1, &p_xi1)); // Q_i² = 0. for i in 0..3 { diff --git a/ext/crates/algebra/src/algebra/motivic/tau.rs b/ext/crates/algebra/src/algebra/motivic/tau.rs index 685262de24..1a6e53c526 100644 --- a/ext/crates/algebra/src/algebra/motivic/tau.rs +++ b/ext/crates/algebra/src/algebra/motivic/tau.rs @@ -8,7 +8,7 @@ //! sentinel for zero. That is [`Tau`]. //! //! This is a small, self-contained scalar type for the motivic layer only. It is -//! deliberately **not** an [`fp`](fp) replacement threaded through the resolution +//! deliberately **not** an [`fp`] replacement threaded through the resolution //! engine: the engine stays $\mathbb{F}_p$-only (see the `motivic` module docs), //! and `Tau` is used solely by the motivic algebra and the Phase 2 lift. @@ -24,11 +24,10 @@ pub struct Tau(Option); impl Tau { - /// The zero coefficient. - pub const ZERO: Self = Self(None); - /// The unit $\tau^0 = 1$. pub const ONE: Self = Self(Some(0)); + /// The zero coefficient. + pub const ZERO: Self = Self(None); /// The zero coefficient (method form, for symmetry with [`Tau::one`]). pub fn zero() -> Self {