diff --git a/.cargo/config-avx512.toml b/.cargo/config-avx512.toml new file mode 100644 index 000000000..721e59fbb --- /dev/null +++ b/.cargo/config-avx512.toml @@ -0,0 +1,18 @@ +# AVX-512 build. Mirrors ndarray/.cargo/config-avx512.toml (`sapphirerapids`). +# cargo build --config .cargo/config-avx512.toml ... +# The `simd.rs` dispatch picks the `simd_avx512` backend (native __m512d). +# +# `sapphirerapids` is a SUPERSET of most AVX-512 silicon in the fleet. It +# implies x86-64-v4 PLUS avx512_bf16 / avx512_fp16 / AMX. On Sapphire / +# Emerald / Granite Rapids that is the right pick. On Cascade Lake (family 6 +# model 0x55 — what this session's host is) or Ice Lake, it is NOT: AMX is +# inline-asm and runtime-gated by `amx_available()`, so it stays inert, but +# anything compiled under `cfg(target_feature = "avx512bf16")` will SIGILL on +# silicon that lacks it. There, use one of: +# CARGO_BUILD_RUSTFLAGS='-Ctarget-cpu=x86-64-v4' # what the probe was run with +# cargo build --config .cargo/config-native.toml # host CPUID, non-portable +# Rule of thumb: `sapphirerapids` for a KNOWN SPR/EMR/GNR deployment target; +# `x86-64-v4` for portable AVX-512; `native` for this-machine-only. +[build] +[target.'cfg(target_arch = "x86_64")'] +rustflags = ["-Ctarget-cpu=sapphirerapids"] diff --git a/.cargo/config-native.toml b/.cargo/config-native.toml new file mode 100644 index 000000000..bf9ed2396 --- /dev/null +++ b/.cargo/config-native.toml @@ -0,0 +1,6 @@ +# Build-machine-tuned. rustc resolves the host CPUID at compile time; the +# binary is NOT portable off this machine. Mirrors ndarray/.cargo/config-native.toml. +# cargo build --config .cargo/config-native.toml ... +[build] +[target.'cfg(target_arch = "x86_64")'] +rustflags = ["-Ctarget-cpu=native"] diff --git a/.cargo/config.toml b/.cargo/config.toml new file mode 100644 index 000000000..6dd5d98a7 --- /dev/null +++ b/.cargo/config.toml @@ -0,0 +1,34 @@ +# lance-graph cargo config — mirrors ndarray's `.cargo/config.toml` policy +# (operator ruling, 2026-09-02): GitHub needs v3; everything else v4 or native. +# +# ndarray is the MANDATORY execution substrate for every crate here, and its +# `simd.rs` dispatch is COMPILE-TIME (`cfg(target_feature = "avx512f")` / +# `"avx2"`). Whatever `target-cpu` this build runs at decides which backend +# `ndarray::simd::{F64x8, F32x16, I16x32}` resolves to. `.cargo/config.toml` +# is read from the INVOKING directory, not the dependency's — so ndarray's +# own config does NOT apply when a lance-graph crate is built from here. +# Without this file a local build compiles the substrate at generic x86-64 +# (SSE2), i.e. the scalar arm. +# +# Measured on the Goursat wavefront (crates/jc/examples/goursat_substrate_probe, +# 4097-point paths, release): under v4 `F64x8::mul_add` is 8.6x faster than +# the flat scalar loop; at generic x86-64 the SAME code is 4.5x SLOWER than +# scalar, because the polyfill's scalar arm is emulating eight lanes through +# arrays. A build with no target-cpu is a regression, not a no-op. +# +# Default here is v3 (AVX2) for the same reason ndarray chose it: portable +# across all x86_64 since Haswell, and what CI runs (`.github/workflows/*.yml` +# set `-C target-cpu=x86-64-v3` explicitly). Anything below v3 SIGILLs on the +# `simd_avx2` intrinsics. +# +# AVX-512 (v4): cargo build --config .cargo/config-avx512.toml ... +# or CARGO_BUILD_RUSTFLAGS='-Ctarget-cpu=x86-64-v4' +# Build machine: cargo build --config .cargo/config-native.toml ... +# +# The two crypto cfgs ndarray's config also carries (`curve25519_dalek_backend`, +# `poly1305_force_soft`) are deliberately NOT repeated: neither crate is in +# this workspace's Cargo.lock (verified 2026-09-02, 0 hits), so there is no +# AVX2 crypto surface here to compile out. Re-check if `encryption` or +# `ed25519-dalek` ever enters the tree. +[target.'cfg(target_arch = "x86_64")'] +rustflags = ["-Ctarget-cpu=x86-64-v3"] diff --git a/.claude/board/EPIPHANIES.md b/.claude/board/EPIPHANIES.md index 869bbd8fc..b49fce5d4 100644 --- a/.claude/board/EPIPHANIES.md +++ b/.claude/board/EPIPHANIES.md @@ -1,3 +1,56 @@ +## 2026-09-02 — E-PILLAR-11-PUBLISHED-BOUND-NEEDS-ITS-OWN-NUMERIC-GUARD-1 + +**Status:** FINDING (measured) — closes the Hambly-Lyons edition-constant arc. +**Confidence:** every number below was run, not inferred; the disable arm was +executed and restored. + +The published constant correction (#1133/`1a3c7294`, arXiv `e` → +Annals `2e`) doubles the certificate depth. Measured consequences: + +| | old arXiv `e` | corrected Annals `2e` | +|---|---|---| +| `c` | `e·ln(1+√2)` = 2.3958 | `2e·ln(1+√2)` = 4.79155 | +| depth at `L=3` | 7 | **14** | +| coefficients (`d=2`) | 255 | **32 767** | +| `pillar_passes`, debug | 25.98 s | **26.17 s** | + +**A 128x coefficient increase costs ~0.2 s inside a 26 s run — inside noise.** +The exhaustive arm's cost is the enumeration of reduced words of length <= 3 +(484 of them), not the signature depth, so raising the cutoff moves nothing. +The other three tests are 0.07-0.10 s each, the `L=8` figure-of-eight arm +included: the correction only LOOSENED its bound (`<= 19` -> `<= 38`), making +it cheaper to satisfy. The 26 s is the PDE depth-infinity leg (8 Goursat solves +on 4609-point paths in debug) — established by elimination, since the constant +flip did not move it and the other tests are ~0.1 s, and corroborated by the +release figure of 0.18 s for a 4097-point solve. Not directly instrumented. + +**The load-bearing half: the behavioural tests do NOT guard the constant.** +Run with the wrong (arXiv `e`) constant, `pillar_passes` **still passes**. What +fails is `theorem2_depth_is_the_paper_floor`: + +``` +assertion `left == right` failed: L=1 + left: 2 floor(e·ln(1+√2)) + right: 4 floor(2e·ln(1+√2)) +``` + +That is the correct division of labour, not a gap. The theorem gives a +SUFFICIENT depth; the fixtures separate at depth 3, far inside either bound, so +a looser constant cannot fail a behaviour test that clears it comfortably. + +**Consequence, and the reason this is banked:** `theorem2_depth_is_the_paper_floor` +is an independent PROVENANCE guard — it is the only thing standing between the +tree and a silent regression to the arXiv edition. Beside a green +`pillar_passes` it reads like a tautology restating a constant, which is exactly +how a future tidy-up would justify deleting it. **It must not be removed as +redundant.** Numbering provenance is pinned at the definition site +(`hambly_lyons.rs:365`): Annals Theorem 5 = arXiv v2 Theorem 2; Annals +Theorem 2 is NOT the lattice cutoff. + +**Generalisation:** when a constant enters a SUFFICIENT bound, behavioural tests +downstream of it are structurally incapable of pinning it — any bound at least +as large passes. Such a constant needs its own numeric guard, and that guard's +apparent triviality is a property of correct design, not a smell. ## 2026-09-02 — E-THE-DTO-LADDER-IS-THE-ALU-BUS-AND-WAS-ALREADY-RULED-1 — a session re-asked a question `.claude/v3/` had answered on 2026-07-02, and proposed a probe on the wrong axis before reading it **Status:** CORRECTION + FENCE (same day, same branch, before merge). **Confidence:** High. diff --git a/.claude/board/TECH_DEBT.md b/.claude/board/TECH_DEBT.md index ac81df884..91e87705a 100644 --- a/.claude/board/TECH_DEBT.md +++ b/.claude/board/TECH_DEBT.md @@ -1,3 +1,162 @@ +## TD-PILLAR11-SCIENTIFIC-LOOPS-BYPASS-NDARRAY-SIMD-1 (2026-09-02) — OPEN + +**The debt is not missing SIMD support. It is scientific code bypassing the +already-complete `ndarray` execution substrate.** (Operator ruling, 2026-09-02.) + +No violation was committed — the Pillar-11 legs contain zero raw intrinsics +(`core::arch` / `_mm*` / `target_feature`: 0 hits in `crates/sigker/src` and +`crates/jc/src`). They are scalar `f64`, hand-written beside a substrate that +already owns the machine vocabulary. + +**The law, corrected from an earlier weaker reading.** There is no +"SIMD op exists -> polyfill, else -> scalar fallback" decision at consumer +level. `ndarray::simd` hides the architecture choice behind ONE bit-exact typed +surface, and **scalar is one BACKEND of that surface**, not a consumer-authored +alternative: `F64x8` resolves to AVX-512 (`__m512d`), AVX2 (`f64x4` x2), NEON +(`[float64x2_t; 4]`), wasm32+simd128 (`[v128; 4]`), and `scalar::F64x8` for +"other non-x86 targets ... full scalar fallback" (`simd.rs` dispatch arms, +verified 2026-09-02). So: + +**Method-level parity CONFIRMED (operator, 2026-09-02):** *"every backend +implements every method with parity coverage."* The dispatch check above only +established that each backend EXPORTS the type; this closes the stronger +question a composing consumer actually depends on — that a method reached +through the typed surface is implemented, and bit-exact, on every arm. It is +what makes "compose from polyfill methods" a safe instruction rather than a +per-method gamble, and it removes the last conditional from the law below. + +``` +named high-level algorithm exists -> call it +it does not exist -> COMPOSE it from ndarray polyfill methods +never -> a consumer-local scalar arithmetic path +never -> consumer-local intrinsics +``` + + JC tells us what is mathematically true. + ndarray tells the machine how to execute it. + Everything else is composition, never a second arithmetic implementation. + +**For the Goursat kernel the arithmetic is entirely available; the only +unsolved part is algorithmic SCHEDULING.** `signature_kernel_pde` +(`sigker/src/kernel.rs:106-131`) computes + +``` +k[i+1][j+1] = k[i+1][j] + k[i][j+1] - k[i][j] + c_ij·k[i][j] + = mul_add(c_ij, diag, left + up - diag) +``` + +`mul_add`, `from_slice`, `copy_to_slice`, `reduce_sum`, `select` all ship on +`F64x8` today. What blocks it is that `k[i+1][j+1]` reads `k[i+1][j]` — a +strict serial recurrence along `j` in row-major order. Independence lives on +the **anti-diagonal**; the transformation needed is +`row-major serial recurrence -> anti-diagonal / rolling-wavefront formulation`. + +**A1 determines A2's shape — the architectural reason not to jump from +`Vec` straight to SIMD.** In a flat row-major buffer the anti-diagonal is +STRIDED (stride `m-1`), so a naive wavefront needs gather. **Three rolling +anti-diagonal buffers make every wavefront load contiguous and can eliminate +the gather entirely.** The storage decision therefore fixes which lane ops A2 +needs at all; it is not a warm-up measurement. + +**W1.5 falsifier (do not run from a banking session):** + +| arm | shape | relation | +|---|---|---| +| A0 | current `Vec>`, row-major | reference | +| A1 | flat / rolling storage, SAME recurrence and order | **A0 = A1 exactly** — only storage changes | +| A2 | rolling anti-diagonal traversal via `ndarray::simd::method()` | A1 <-> A2 gets a PREDECLARED solver tolerance: traversal changes evaluation order | + +The A1<->A2 tolerance is stated in advance, never discovered after the fact. + +**Inventory — the scalar surfaces in the same two crates.** Straightforwardly +canonical (`reduce_sum` over `mul_add`): `LogSignature::{dot, cosine}`, +`RandomizedSignature::{dot, cosine}`, `linear_path_kernel_closed_form`. +Accumulation-shaped: `signature_truncated`, `log_signature_truncated`, +`RandomizedSignature::encode`, `hydrate_signature`, `signature_kernel`, +`signature_kernel_normalized`. Combinatorial rather than arithmetic, and +plausibly staying scalar: `shuffle_product`, `enumerate_lyndon_words`, +`witt_component` / `witt_dimension`. + +**Not shipped:** `signature_pde_sweep` and `shuffle_product_lift` return zero +hits in `ndarray/src/` — they remain the W1.5 catalogue's shopping list. The +missing thing is traversal/composition, NOT arithmetic; whether a composition +is later promoted to a named `signature_pde_sweep()` is purely an API/reuse +question, not a blocker. + +**Pay by:** W1.5, gated on `jc Pillar 11` (green for the lattice leg, +`7751581f`). **Not W5** — W5 is workload-pressure machinery and stays HOLD +(4609 vs the 11 585-point 1 GiB threshold) regardless of any speedup here. + +**A0/A1/A2 RUN (2026-09-02, same day — operator: "you didn't try the 25-26 +seconds with ndarray yet").** `crates/jc/examples/goursat_substrate_probe.rs`, +release, 4097-point paths (16.8M cells). Every A2 lane op is +`ndarray::simd::F64x8::mul_add` — the body is three FMAs +(`t = mul_add(1,left,up)`, `u = mul_add(-1,diag,t)`, `new = mul_add(c,diag,u)`); +no `Add`/`Sub` operator was needed, and none was minted. + +| target-cpu | backend | A0 | A1 | A2 | A1/A2 | A0=A1 | \|A1-A2\|/A1 | +|---|---|---|---|---|---|---|---| +| generic x86-64 (no target-cpu) | scalar arm | 0.182 s | 0.169 s | **0.773 s** | **0.22x** | bit-exact | 2.431e-13 | +| x86-64-v3 (AVX2) | `f64x4` x2 | 0.163 s | 0.150 s | 0.0165 s | **9.12x** | bit-exact | 2.431e-13 | +| x86-64-v4 (AVX-512) | `__m512d` | 0.171 s | 0.157 s | 0.0182 s | **8.62x** | bit-exact | 2.431e-13 | + +Four findings, each falsifiable and each measured: + +1. **Storage was NOT the wall; the recurrence was.** A1/A0 = 1.09x. The + hypothesis banked above — that a flat buffer alone would close most of the + gap — is **FALSIFIED**. The probe was built to answer that and it did. +2. **A build with no `target-cpu` is a REGRESSION, not a no-op.** At generic x86-64 the + polyfill's scalar arm runs the wavefront 4.5x SLOWER than the flat scalar + loop: eight lanes emulated through arrays plus the wavefront bookkeeping. + Until this commit lance-graph had NO `.cargo/config.toml`, so every local + build landed there. `.cargo/config.toml` now pins v3 (CI already did via + `RUSTFLAGS` in `.github/workflows/*.yml`), with `config-avx512.toml` / + `config-native.toml` mirroring ndarray's. +3. **Bit-exact across backends, as confirmed.** `|A1-A2|/A1` is identical to the + last digit on scalar, AVX2 and AVX-512 at every size. The delta itself + (~1e-13 at 4097) is the fused-vs-separate rounding of `c·diag` + accumulated over 16.8M cells — A2 is the MORE accurate arm. +4. **AVX2 ~ AVX-512 here** (9.1x vs 8.6x): the wavefront is latency-bound on + the diagonal recurrence, not width-bound. Widening lanes buys nothing until + the dependency chain is restructured; that is a scheduling question, not a + substrate one. + +**Correction (2026-09-02, same session):** an earlier revision of this entry +labelled the no-`target-cpu` row `"386"` and attributed the phrase to the +operator. That was a misreading of `x86-64-v4` in a terse message. The +measurement is unchanged; only the label was wrong and is removed above. The +commit message on `5df2d785` still carries it and is not rewritten (pushed). +Also recorded: `.cargo/config-avx512.toml`'s `sapphirerapids` is a SUPERSET of +Cascade Lake / Ice Lake silicon — this session's host is Cascade Lake (family 6 +model 0x55; `amx_report()`: `cpu_model()=OtherX86`, `expects_amx=false`). On +such hosts `x86-64-v4` or `native` is the correct pick; the probe was run with +`x86-64-v4`, so its numbers stand. + +**SPR vs EMR — resolved (2026-09-02, operator refinement + ndarray git history):** +not the enablement, the DETECTION. Pre-PR-#217 (`src/simd_caps.rs` @ `bdf243cc`, +2026-06-13) detected AMX by CPUID feature bits alone (`amx_tile`/`amx_int8`/ +`amx_bf16`/`amx_fp16`, EDX bits 24/25/22) — no XCR0 gate, no model table, and +the `arch_prctl` issued on syscall 157, so it always failed. PR #217 +(`e563fdcd`, 2026-06-14) replaced it with the four-gate detector (CPUID + +OSXSAVE + XCR0 + `arch_prctl` 158) PLUS the CPUID model table (`CpuModel`: +SPR 0x8F / EMR 0xCF / GNR 0xAD,0xAE / SRF 0xAF), added to tell "no silicon" +from "not OS-enabled". On SPR the old detector said *present* while nothing +ever executed (every tile test early-returned — Gotcha 9); on EMR the new +detector said *present AND enabled* and tiles ran. The `arch_prctl` grant is +the same on both; what differed was the detection code, and the change landed +on EMR silicon. `amx-enablement-and-kernel.md` §2 says this ("EMR was simply +the host where gate 4 got fixed first"); the operator's "detected differently" +is the same fact from the outside. Minor inconsistency noticed, not chased: +`cpu_ops.rs:186` says "Linux 5.19+", the doc says "5.16+". + +Storage detail that held: with `dy` stored REVERSED, the anti-diagonal walk +is forward in `i`, so k-buffers, `dx` and `dy` are all contiguous slices — +**no gather**, exactly as predicted by "A1 determines A2's shape". + +`jc` now depends on `ndarray` as a plain, non-optional `[dependencies]` +entry; its "zero external deps in production" header is retired (operator: +ndarray is mandatory everywhere). What stays standalone is the PROOF. + ## TD-GHOST-TIER-NAME-COLLISION-1 (2026-09-02) — OPEN, doc-only `crates/lance-graph-contract/src/counterfactual.rs` calls the −6 minority-pole diff --git a/crates/jc/Cargo.lock b/crates/jc/Cargo.lock index 9be3afa8b..550bd302f 100644 --- a/crates/jc/Cargo.lock +++ b/crates/jc/Cargo.lock @@ -11,6 +11,12 @@ dependencies = [ "memchr", ] +[[package]] +name = "autocfg" +version = "1.5.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f2032f911046de80f0a198e0901378627c33f59ea0ac00e363d481118bd70a53" + [[package]] name = "equivalent" version = "1.0.2" @@ -57,6 +63,7 @@ name = "jc" version = "0.1.0" dependencies = [ "lance-graph-contract", + "ndarray", "neural-debug", "p64", "sigker", @@ -77,12 +84,38 @@ version = "0.2.16" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "b6d2cec3eae94f9f509c767b45932f1ada8350c4bdb85af2fcab4a3c14807981" +[[package]] +name = "matrixmultiply" +version = "0.3.11" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "3f607c237553f086e7043417a51df26b2eb899d3caff94e6a67592ff992fedc7" +dependencies = [ + "autocfg", + "rawpointer", +] + [[package]] name = "memchr" version = "2.8.3" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "cf8baf1c55e62ffcace7a9f06f4bd9cd3f0c4beb022d3b367256b91b87513d98" +[[package]] +name = "ndarray" +version = "0.17.2" +dependencies = [ + "fractal", + "matrixmultiply", + "num-complex", + "num-integer", + "num-traits", + "p64", + "paste", + "portable-atomic", + "portable-atomic-util", + "rawpointer", +] + [[package]] name = "neural-debug" version = "0.1.0" @@ -93,6 +126,33 @@ dependencies = [ "walkdir", ] +[[package]] +name = "num-complex" +version = "0.4.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "73f88a1307638156682bada9d7604135552957b7818057dcef22705b4d509495" +dependencies = [ + "num-traits", +] + +[[package]] +name = "num-integer" +version = "0.1.47" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "7ce2d95d4b3734dc35aa2f45e1aa22cd416814592a4f9d9205e11affd5b8e10b" +dependencies = [ + "num-traits", +] + +[[package]] +name = "num-traits" +version = "0.2.19" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "071dfc062690e90b734c0b2273ce72ad0ffa95f0c74596bc250dcfd960262841" +dependencies = [ + "autocfg", +] + [[package]] name = "p64" version = "0.1.0" @@ -100,6 +160,27 @@ dependencies = [ "fractal", ] +[[package]] +name = "paste" +version = "1.0.15" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "57c0d7b74b563b49d38dae00a0c37d4d6de9b432382b2892f0574ddcae73fd0a" + +[[package]] +name = "portable-atomic" +version = "1.15.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "05c8b63e8d9609db387f0324918f81d68fe27748f084ef092fb35954d0539a85" + +[[package]] +name = "portable-atomic-util" +version = "0.2.7" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "c2a106d1259c23fac8e543272398ae0e3c0b8d33c88ed73d0cc71b0f1d902618" +dependencies = [ + "portable-atomic", +] + [[package]] name = "proc-macro2" version = "1.0.106" @@ -118,6 +199,12 @@ dependencies = [ "proc-macro2", ] +[[package]] +name = "rawpointer" +version = "0.2.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "60a357793950651c4ed0f3f52338f53b2f809f32d83a07f72909fa13e4c6c1e3" + [[package]] name = "regex" version = "1.13.1" diff --git a/crates/jc/Cargo.toml b/crates/jc/Cargo.toml index 639f3a9d4..8657b4d21 100644 --- a/crates/jc/Cargo.toml +++ b/crates/jc/Cargo.toml @@ -5,34 +5,47 @@ edition = "2021" description = "Jirak-Cartan: five-pillar proof-in-code for binary-Hamming causal field computation" license = "Apache-2.0" -# Zero EXTERNAL deps in production — standalone, like deepnsm and bgz17. -# The proof is the proof regardless of SIMD path. -# Workspace-sibling path-deps are opt-in via features (see [features] below). - -# Optional workspace-sibling deps — gated by feature flags to preserve -# the default zero-dep build. `cargo build` (default) gives a fully -# standalone JC; `cargo build --features hambly-lyons` activates Pillar 11 -# by pulling in the sigker workspace sibling. +# ndarray is MANDATORY — here and everywhere in lance-graph (operator ruling, +# 2026-09-02). It is the execution substrate: every numeric path is written +# against `ndarray::simd::{F64x8, F32x16, I16x32}::*`, whose scalar / AVX2 / +# AVX-512 / NEON / WASM backends are bit-exact behind one typed surface. JC +# certifies what is mathematically true; ndarray tells the machine how to +# execute it. A consumer-local scalar arithmetic path is not a fallback, it +# is a second arithmetic implementation, and is forbidden. +# +# The former "zero external deps in production" constitution is retired by +# that ruling. What remains standalone is the PROOF: no other crate's +# behaviour is a premise of any pillar. Workspace siblings that carry +# domain content (sigker) stay feature-gated so a pillar can be run alone. [dependencies] +ndarray = { path = "../../../ndarray" } sigker = { path = "../sigker", optional = true } [features] -# Default build is zero-dep — honors the standalone constitution. +# Default build: ndarray only (mandatory). Pillar 11 opts in below. default = [] # Activates Pillar 11 (Hambly-Lyons signature uniqueness) by pulling in the # sigker workspace sibling. See `src/hambly_lyons.rs` for the probe + the # DEFERRED fallback used when the feature is off. hambly-lyons = ["dep:sigker"] -# Dev-only deps for cross-crate bridge examples (production stays zero-dep). +# Dev-only deps for cross-crate bridge examples. [dev-dependencies] lance-graph-contract = { path = "../lance-graph-contract" } # `substrate_compare` example only — the on-demand P64-vs-V3 comparison drives # the REAL Palette64 (compute) + NodeGuid (address) and records into the -# neural-debug runtime registry. Dev-deps, so production JC stays zero-dep. +# neural-debug runtime registry. Dev-deps: examples only, not the library. p64 = { path = "../../../ndarray/crates/p64" } neural-debug = { path = "../neural-debug" } +[[example]] +# A0/A1/A2 Goursat probe (TD-PILLAR11-SCIENTIFIC-LOOPS-BYPASS-NDARRAY-SIMD-1). +# Imports `sigker::signature_kernel_pde` as the A0 reference, so it only +# exists when Pillar 11 is on — CI builds examples with default features and +# an unconditional `use sigker` there is a hard E0432, not a warning. +name = "goursat_substrate_probe" +required-features = ["hambly-lyons"] + [[example]] name = "prove_it" diff --git a/crates/jc/examples/goursat_substrate_probe.rs b/crates/jc/examples/goursat_substrate_probe.rs new file mode 100644 index 000000000..e47840f8c --- /dev/null +++ b/crates/jc/examples/goursat_substrate_probe.rs @@ -0,0 +1,184 @@ +//! A0/A1/A2 — the Goursat solve on the canonical ndarray substrate. +//! +//! The falsifier recorded in `TD-PILLAR11-SCIENTIFIC-LOOPS-BYPASS-NDARRAY-SIMD-1`: +//! +//! | arm | storage | traversal | arithmetic | +//! |---|---|---|---| +//! | A0 | `Vec>` | row-major | scalar (the shipped `signature_kernel_pde`) | +//! | A1 | flat `Vec` | row-major, SAME order | scalar | +//! | A2 | three rolling anti-diagonals | wavefront | `ndarray::simd::F64x8::mul_add` | +//! +//! Contract: A0 = A1 **bit-exact** (only storage changes). A1 <-> A2 under a +//! predeclared tolerance: the wavefront's `mul_add(c, diag, u)` fuses the +//! `c·diag + u` rounding the scalar path performs as two roundings, so the two +//! arms legitimately differ in the last bits. No consumer intrinsics; every +//! lane op is `ndarray::simd::method()`. +//! +//! Run: `cargo run --release --manifest-path crates/jc/Cargo.toml \ +//! --features hambly-lyons --example goursat_substrate_probe` + +use ndarray::simd::F64x8; +use sigker::signature_kernel_pde; +use std::time::Instant; + +const LANES: usize = 8; + +fn path(n: usize) -> Vec> { + (0..=n) + .map(|i| { + let t = i as f64 / n as f64; + vec![ + t + 0.05 * (260.0 * t).cos(), + 0.5 * t + 0.05 * (260.0 * t).sin(), + ] + }) + .collect() +} + +/// Increments dx_i = x[i+1] - x[i], one contiguous block per step. +fn increments(x: &[Vec]) -> (Vec, usize) { + let dim = x[0].len(); + let mut out = Vec::with_capacity((x.len() - 1) * dim); + for w in x.windows(2) { + out.extend(w[1].iter().zip(&w[0]).map(|(next, prev)| next - prev)); + } + (out, dim) +} + +/// A1 — flat row-major storage, the shipped recurrence in the shipped order. +/// `k[i+1][j+1] = k[i+1][j] + k[i][j+1] - k[i][j] + c_ij * k[i][j]`, evaluated +/// left-to-right exactly as `signature_kernel_pde` writes it. +fn goursat_flat(x: &[Vec], y: &[Vec]) -> f64 { + let (n, m) = (x.len(), y.len()); + let (dx, dim) = increments(x); + let (dy, _) = increments(y); + let mut k = vec![1.0f64; n * m]; + for i in 0..n - 1 { + let dxi = &dx[i * dim..(i + 1) * dim]; + for j in 0..m - 1 { + let dyj = &dy[j * dim..(j + 1) * dim]; + let c: f64 = (0..dim).map(|a| dxi[a] * dyj[a]).sum(); + let (left, up, diag) = (k[(i + 1) * m + j], k[i * m + j + 1], k[i * m + j]); + k[(i + 1) * m + j + 1] = left + up - diag + c * diag; + } + } + k[n * m - 1] +} + +/// A2 — wavefront over three rolling anti-diagonal buffers, every lane op +/// `ndarray::simd::F64x8::mul_add`. +/// +/// On diagonal `d` (cells with `i + j = d`), indexed by row `i`: +/// left = k[i][j-1] = prev1[i] up = k[i-1][j] = prev1[i-1] +/// diag = k[i-1][j-1] = prev2[i-1] c = +/// With `dy` stored REVERSED, `dy[j-1] = dyr[m-1-d+i]` walks FORWARD in `i`, +/// so every operand — k-buffers, dx, dy — is a contiguous slice. No gather. +/// +/// Body, three FMAs (±1.0 multipliers are exact, so the first two round +/// exactly like `+`/`-`; only the last fuses what A1 rounds twice): +/// t = mul_add( 1, left, up) u = mul_add(-1, diag, t) new = mul_add(c, diag, u) +fn goursat_wavefront(x: &[Vec], y: &[Vec]) -> f64 { + let (n, m) = (x.len(), y.len()); + let (dx, dim) = increments(x); + let (dy, _) = increments(y); + assert_eq!(dim, 2, "probe fixes dim = 2 (the W5 path shape)"); + // Split by component; reverse dy so the anti-diagonal walk is forward. + let dx0: Vec = dx.iter().step_by(2).copied().collect(); + let dx1: Vec = dx.iter().skip(1).step_by(2).copied().collect(); + let mut dyr0: Vec = dy.iter().step_by(2).copied().collect(); + let mut dyr1: Vec = dy.iter().skip(1).step_by(2).copied().collect(); + dyr0.reverse(); + dyr1.reverse(); + + let mut prev2 = vec![1.0f64; n]; + let mut prev1 = vec![1.0f64; n]; + let mut cur = vec![1.0f64; n]; + let (one, neg_one, zero) = (F64x8::splat(1.0), F64x8::splat(-1.0), F64x8::splat(0.0)); + let mut out = [0.0f64; LANES]; + + for d in 2..(n + m - 1) { + // Boundaries on this diagonal: k[0][d] and k[d][0] are 1. + if d < m { + cur[0] = 1.0; + } + if d < n { + cur[d] = 1.0; + } + // Interior rows: i >= 1, j = d - i >= 1, i <= n-1, j <= m-1. + let lo = 1usize.max(d.saturating_sub(m - 1)); + let hi = (d - 1).min(n - 1); + if lo > hi { + std::mem::swap(&mut prev2, &mut prev1); + std::mem::swap(&mut prev1, &mut cur); + continue; + } + // dyr index for row i is (m-1-d)+i; the difference may be negative but the sum is not. + let base = (m - 1).wrapping_sub(d); + let mut i = lo; + while i + LANES <= hi + 1 { + let left = F64x8::from_slice(&prev1[i..i + LANES]); + let up = F64x8::from_slice(&prev1[i - 1..i - 1 + LANES]); + let diag = F64x8::from_slice(&prev2[i - 1..i - 1 + LANES]); + let a0 = F64x8::from_slice(&dx0[i - 1..i - 1 + LANES]); + let a1 = F64x8::from_slice(&dx1[i - 1..i - 1 + LANES]); + let r = base.wrapping_add(i); // == m-1-d+i, in range for interior rows + let b0 = F64x8::from_slice(&dyr0[r..r + LANES]); + let b1 = F64x8::from_slice(&dyr1[r..r + LANES]); + let c = a1.mul_add(b1, a0.mul_add(b0, zero)); + let t = one.mul_add(left, up); + let u = neg_one.mul_add(diag, t); + c.mul_add(diag, u).copy_to_slice(&mut out); + cur[i..i + LANES].copy_from_slice(&out); + i += LANES; + } + // Scalar tail: same three-FMA arithmetic, so A2 is internally uniform. + while i <= hi { + let r = base.wrapping_add(i); + let c = dx1[i - 1].mul_add(dyr1[r], dx0[i - 1] * dyr0[r]); + let t = 1.0f64.mul_add(prev1[i], prev1[i - 1]); + let u = (-1.0f64).mul_add(prev2[i - 1], t); + cur[i] = c.mul_add(prev2[i - 1], u); + i += 1; + } + std::mem::swap(&mut prev2, &mut prev1); + std::mem::swap(&mut prev1, &mut cur); + } + // Final cell k[n-1][m-1] is on diagonal n+m-2, the last computed, now in prev1. + prev1[n - 1] +} + +fn main() { + println!( + "{:>6} {:>10} {:>10} {:>10} {:>9} {:>9} {:>12} {:>12}", + "len", "A0 secs", "A1 secs", "A2 secs", "A0/A1", "A1/A2", "A0==A1", "|A1-A2|/A1" + ); + for &n in &[256usize, 1024, 2048, 4096] { + let (x, y) = (path(n), path(n)); + let t = Instant::now(); + let a0 = signature_kernel_pde(&x, &y); + let s0 = t.elapsed().as_secs_f64(); + let t = Instant::now(); + let a1 = goursat_flat(&x, &y); + let s1 = t.elapsed().as_secs_f64(); + let t = Instant::now(); + let a2 = goursat_wavefront(&x, &y); + let s2 = t.elapsed().as_secs_f64(); + let exact = a0.to_bits() == a1.to_bits(); + let rel = ((a1 - a2) / a1).abs(); + println!( + "{:>6} {s0:>10.4} {s1:>10.4} {s2:>10.4} {:>8.2}x {:>8.2}x {:>12} {rel:>12.3e}", + n + 1, + s0 / s1, + s1 / s2, + if exact { "bit-exact" } else { "DIFFERS" } + ); + assert!( + exact, + "A0 != A1 at n={n}: {a0:e} vs {a1:e} — storage change altered the result" + ); + } + println!("\nA0 = A1 bit-exact at every size (only storage changed)."); + println!( + "A1 <-> A2 differ by fused vs. separate rounding of c·diag — the predeclared tolerance." + ); +}