simd: verify-pass fixes — backend engine on the unfused reference tier

Three-angle adversarial verify (numerics/swap-trace/docs) on the
completion diff. Substantive finding acted on: the fused tier's scalar
polyfill (AVX2 arm F64x8::mul_add = per-lane f64::mul_add) lowers to a
libm fma() call on baseline x86-64 builds without the fma target
feature — downstream consumers do not inherit this repo's .cargo pin,
so backend::native::gemm_f64 now routes through the UNFUSED reference
tier instead: no libm dependence on any backend, ~7% slower than fused
on pinned builds (10.0 vs 10.7 GFLOP/s), and the backend engine is now
bit-identical to the certification reference. gemm_f64_tiled_fma stays
public for FMA-pinned consumers.

Doc fixes from the same pass: removed two stale 'gemm_f64 delegates to
matrixmultiply' claims that contradicted the swap (simd.rs re-export
comment, gemm_f64_tiled rustdoc); added # Panics to gemm_f64 (checked
preconditions vs the old wrapper's silent UB on short slices, matching
CBLAS xerbla); scoped fma determinism to per-(build,runtime) — wasm
relaxed-simd fusion is implementation-defined; corrected the AVX2
vfmadd naming (per-lane polyfill, fused semantics); fixed the
integer-corpus bound comment (k_max=128).

Gates: clippy -D warnings clean, 2185/2185 lib tests, 4 gemm doctests.

Author: claude

.claude/blackboard.md

@@ -24,11 +24,22 @@ below's foundation:
    f64→f64 store/reload never rounds ⇒ per-element op sequence
    unchanged ⇒ every bit-equality test green untouched. [MEASURED]
    ref 4.6→10.0 GF, fma 4.7→10.7 GF (2.2×).
-3. **Engine swap:** `backend::native::gemm_f64` now routes to
-   `gemm_f64_tiled_fma` — the f64 GEMM behind `BlasFloat::backend_gemm`
-   / `hpc::blas_level3::blas_gemm` / batched linalg is entirely own
-   Rust; matrixmultiply remains only in gemm_f32 and upstream
-   `Array::dot` (impl_linalg.rs, untouched at ~33 GF).
+3. **Engine swap:** `backend::native::gemm_f64` now routes to the
+   crate-native tiled kernel — the f64 GEMM behind
+   `BlasFloat::backend_gemm` / `hpc::blas_level3::blas_gemm` / batched
+   linalg is entirely own Rust; matrixmultiply remains only in gemm_f32
+   and upstream `Array::dot` (impl_linalg.rs, untouched at ~33 GF).
+   **[REVISED post-verify] Engine = the UNFUSED reference tier**, not
+   fma: the verify pass surfaced a cliff — the AVX2-polyfill/scalar
+   `mul_add` lowers to a libm `fma()` call on baseline x86-64 builds
+   (consumers do NOT inherit this repo's `.cargo` target-cpu pin; CI
+   lands exactly there). The unfused tier has no libm dependence on any
+   backend, costs only ~7% vs fused on pinned builds (10.0 vs 10.7 GF),
+   and makes the backend engine bit-identical to the certification
+   reference. `gemm_f64_tiled_fma` stays public for FMA-pinned
+   consumers. New panic contract documented on `gemm_f64` (# Panics —
+   checked preconditions vs the old wrapper's silent-UB on short
+   slices; matches CBLAS xerbla behavior).
 
 [MEASURED, 3-engine, this VM (v3 compile → AVX2 arm, PREFERRED_F64_LANES=4;
 host runtime has avx512f but committed .cargo config is v3)]:
@@ -46,6 +57,18 @@ still dead in native.rs `mod scalar` (f32 sibling completion);
 avx512f compile arm untested on CI (v3 config) — the F64x8=__m512d arm
 runs only on local v4 builds.
 
+**[VERIFY OUTCOME]** 3-angle adversarial pass on the completion diff:
+numerics PASS / swap-trace PASS / docs FAIL→fixed. Substantive P1 acted
+on: baseline-x86 libm-fma cliff → engine revised to the unfused
+reference tier (see #3 REVISED above). Doc P1s fixed: two stale
+"gemm_f64 delegates to matrixmultiply" claims (simd.rs comment +
+gemm_f64_tiled rustdoc) contradicted the swap in the same diff. P2s
+fixed: fma determinism scoped to per-(build,runtime) (wasm relaxed-simd
+fusion is implementation-defined); AVX2 vfmadd naming corrected (per-
+lane f64::mul_add polyfill, fused semantics); integer-corpus bound
+comment corrected (k_max=128, ≈2.1e4). All gates re-run green after
+fixes: clippy -D warnings, 2185/2185 lib tests, 4 gemm doctests.
+
 ## 2026-07-06 (later) — `ndarray::simd::gemm_f64_tiled` surfaced (operator directive)
 
 The crate-native tiled f64 GEMM graduated from dead code

src/backend/native.rs

@@ -238,17 +238,30 @@ pub fn gemm_f32(
 
 /// GEMM: C = alpha * A * B + beta * C (f64, row-major).
 ///
-/// Engine: the crate-native `simd_ops::gemm_f64_tiled_fma` (F64x8-vectorized
-/// fused tier — entirely in-crate Rust, no external matmul engine). The
-/// unfused sibling `ndarray::simd::gemm_f64_tiled` is the bit-exact
-/// certification reference; this entry favors throughput. Swapped from
-/// `matrixmultiply::dgemm` — results agree to BLAS tolerance (last-ulp
-/// accumulation-order/fusion differences), not bit-for-bit.
+/// Engine: the crate-native `simd_ops::gemm_f64_tiled` (F64x8-vectorized,
+/// unfused reference tier — entirely in-crate Rust, no external matmul
+/// engine for f64; `gemm_f32` below still delegates to `matrixmultiply`).
+/// The unfused tier is chosen deliberately: it is bit-identical to the
+/// certification reference on every backend AND has no dependence on the
+/// `fma` target feature (the fused tier's scalar polyfill can lower to a
+/// libm `fma()` call on baseline x86-64 builds that don't inherit this
+/// repo's `.cargo` target-cpu pin). Consumers that pin an FMA-capable
+/// target and want the fused tier call `ndarray::simd::gemm_f64_tiled_fma`
+/// directly. Swapped from `matrixmultiply::dgemm` — results agree to BLAS
+/// tolerance (last-ulp accumulation-order differences), not bit-for-bit.
+///
+/// # Panics
+///
+/// Unlike the previous `matrixmultiply` wrapper (which accepted `lda < k`
+/// silently and was UB on short slices), this entry inherits the tiled
+/// kernel's checked preconditions: panics if `lda < k`, `ldb < n`,
+/// `ldc < n`, or a slice is shorter than `(rows − 1)·ld + cols` — matching
+/// the CBLAS backends, which reject `lda < max(1, k)` via `xerbla`.
 pub fn gemm_f64(
     m: usize, n: usize, k: usize, alpha: f64, a: &[f64], lda: usize, b: &[f64], ldb: usize, beta: f64, c: &mut [f64],
     ldc: usize,
 ) {
-    crate::simd_ops::gemm_f64_tiled_fma(m, n, k, alpha, a, lda, b, ldb, beta, c, ldc);
+    crate::simd_ops::gemm_f64_tiled(m, n, k, alpha, a, lda, b, ldb, beta, c, ldc);
 }
 
 // ─── GEMV dispatch ───────────────────────────────────────────────

src/simd.rs

@@ -570,13 +570,15 @@ pub use crate::simd_ops::{array_chunks, array_chunks_checked, array_windows, arr
 // contract: unfused mul+add in ascending-k order per element → bit-identical
 // on every backend (AVX-512/AVX2/NEON/WASM/scalar) and, at α=1 β=0,
 // bit-identical to the naive triple-loop reference. This is the in-crate
-// ground-truth GEMM for probes/certification — `backend::native::gemm_f64`
-// delegates to the external `matrixmultiply` crate instead. The kernel is
-// alloc-free, but `pub mod simd`/`simd_ops` are std-gated in lib.rs, so
-// like every kernel here it is reachable only in `std` builds today.
+// ground-truth GEMM for probes/certification AND the engine behind
+// `backend::native::gemm_f64` (own Rust in the f64 BLAS path; the f32
+// sibling still delegates to the external `matrixmultiply` crate).
 // `gemm_f64_tiled_fma` is the fast fused tier (same tiling/order, one
-// rounding per step) — the engine behind `backend::native::gemm_f64`;
-// the unfused reference tier stays the certification ground truth.
+// rounding per step) for consumers on FMA-pinned targets — not the
+// backend engine, because its scalar polyfill can lower to libm `fma()`
+// on baseline builds. Both kernels are alloc-free, but `pub mod
+// simd`/`simd_ops` are std-gated in lib.rs, so they are reachable only
+// in `std` builds today.
 pub use crate::simd_ops::{gemm_f64_tiled, gemm_f64_tiled_fma};
 pub use crate::simd_soa::MultiLaneColumn;
 

src/simd_ops.rs

@@ -890,10 +890,11 @@ mod add_mul_tests {
 /// The crate-native f64 GEMM: a TILE=64-blocked loop nest with the
 /// innermost j-run vectorized on the dispatched [`F64x8`] (AVX-512 /
 /// AVX2 / NEON / WASM-SIMD128 / scalar — one source, every backend).
-/// Unlike `backend::native::gemm_f64`, which delegates to the external
-/// `matrixmultiply` crate, this kernel is entirely in-crate and carries
-/// a bit-exactness contract, so probes can certify other kernels
-/// against it as ground truth.
+/// Entirely in-crate with a bit-exactness contract, so probes can
+/// certify other kernels against it as ground truth — and since it IS
+/// the engine behind `backend::native::gemm_f64`, the native f64 BLAS
+/// path carries the same auditable numerics (the f32 sibling still
+/// delegates to the external `matrixmultiply` crate).
 ///
 /// # Bit-exactness contract
 ///
@@ -953,14 +954,21 @@ pub fn gemm_f64_tiled(
 /// Fast tier of [`gemm_f64_tiled`]: same tiling, same ascending-`p`
 /// per-element order, but each step is a **fused** multiply-add
 /// (`c = fma(alpha·A[i,p], B[p,j], c)` — one rounding per step via
-/// `F64x8::mul_add`: `vfmadd`/`vfmaq_f64` on x86/NEON, `f64::mul_add`
-/// on the scalar tail and polyfill backends).
+/// `F64x8::mul_add`; fused semantics on every backend, though the
+/// lowering varies: native `vfmadd` on the AVX-512 arm, `vfmaq_f64` on
+/// NEON, per-lane `f64::mul_add` on the AVX2 polyfill and the scalar
+/// tail — which can become a libm `fma()` call on targets built
+/// without the `fma` feature: correct, but slow. Prefer the reference
+/// tier, or an FMA-capable target pin, in that situation).
 ///
 /// Semantics relative to the reference tier:
-/// - **Deterministic per build**, but NOT bit-stable across backends:
-///   fusion quality differs (WASM without `relaxed-simd` has no fused
-///   vector `mul_add`, so its vector lanes match the unfused reference
-///   while the scalar tail is still fused).
+/// - **Deterministic per (build, runtime)**, but NOT bit-stable across
+///   backends: fusion differs (WASM without `relaxed-simd` has no
+///   fused vector `mul_add`, so its vector lanes match the unfused
+///   reference while the scalar tail is still fused; WASM WITH
+///   `relaxed-simd` uses `f64x2_relaxed_madd`, whose fusion is
+///   implementation-defined — the same wasm binary may round
+///   differently on different runtimes).
 /// - **Bit-identical to [`gemm_f64_tiled`]** whenever every
 ///   `a_val·b` product and running sum is exactly representable —
 ///   e.g. integer-valued data with products and accumulation inside
@@ -1316,8 +1324,9 @@ mod gemm_f64_tiled_tests {
     // ── FMA tier ──
 
     /// Integer-valued corpus: values in [-9, 9], so every a·b product
-    /// (≤ 81) and running sum (≤ 70·81·2) is exactly representable —
-    /// fused and unfused steps round identically → tiers bit-identical.
+    /// (≤ 81) and running sum (≤ k_max·81·|α| + 9 = 128·81·2 + 9 ≈ 2.1e4,
+    /// far below 2^53) is exactly representable — fused and unfused
+    /// steps round identically → tiers bit-identical.
     fn fill_int(len: usize, salt: u64) -> Vec<f64> {
         (0..len)
             .map(|i| ((mix(i as u64 ^ salt) % 19) as i64 - 9) as f64)