pillar11: upgrade the PSD gate to Cholesky, and pin the pipeline bit-exact

Cargo.toml

@@ -91,6 +91,18 @@ name = "codec_mode_histogram"
 required-features = ["codec"]
 
 # W1a-#9 masking-primitive codegen probe imports `ndarray::simd` (std-gated).
+[[example]]
+name = "psd_bit_exactness"
+required-features = ["pillar", "std"]
+
+[[example]]
+name = "psd_jitter_sweep"
+required-features = ["pillar", "std"]
+
+[[example]]
+name = "psd_f32_vs_f64_probe"
+required-features = ["pillar", "std"]
+
 [[example]]
 name = "w1a9_codegen_probe"
 required-features = ["std"]

crates/sigker-parity/tests/w4_depth_infinity_psd.rs

@@ -27,8 +27,13 @@ use ndarray::hpc::pillar::signature::brownian_path_d2;
 use ndarray::hpc::pillar::SplitMix64;
 use sigker::signature_kernel_pde;
 
-/// Same seed as the truncated battery, so both certify the same pool.
-const PILLAR_11_SEED: u64 = 0x5EED_1111_5164_A7AB;
+// The seed is IMPORTED, not retyped. An earlier revision of this file
+// hardcoded a different value while its comment claimed "same seed as the
+// truncated battery, so both certify the same pool" — the two batteries were
+// certifying different pools, and nothing could have noticed, because both
+// sides of every comparison used the wrong constant consistently. Caught by
+// the bit-exactness gate added to the truncated battery.
+use ndarray::hpc::pillar::signature::PILLAR_11_SEED;
 /// Gram pool — O(N^2) kernel solves, so kept small.
 const N_PATHS: usize = 64;
 /// Concentration pool — O(N) solves, matched to the truncated battery's 1000

examples/psd_bit_exactness.rs

@@ -0,0 +1,70 @@
+//! Is the Pillar-11 PSD gate bit-exact, and therefore safe to gate in CI?
+//!
+//! Prints a fingerprint over the RAW BIT PATTERNS of every Gram entry, plus
+//! the Cholesky verdict. Two runs on one machine answer determinism; two runs
+//! under different `-C target-cpu` answer whether autovectorisation/FMA moves
+//! the low bits — which decides whether CI may pin exact values or only
+//! verdicts. `.cargo/config.toml` pins `target-cpu=x86-64-v4`, so this is not
+//! hypothetical.
+use ndarray::hpc::lapack::LapackOps;
+use ndarray::hpc::pillar::signature::PILLAR_11_SEED;
+use ndarray::hpc::pillar::signature::{brownian_path_d2, sigker_hl, signature_d2_deg3};
+use ndarray::hpc::pillar::SplitMix64;
+use ndarray::Array2;
+
+const N_STEPS: usize = 50;
+const SUBSET: usize = 50;
+const JITTER: f64 = 1e-4;
+
+/// FNV-1a over raw bits — any single-bit change moves it.
+fn fnv(acc: &mut u64, bits: u64) {
+    for b in bits.to_le_bytes() {
+        *acc ^= b as u64;
+        *acc = acc.wrapping_mul(0x100_0000_01b3);
+    }
+}
+
+fn main() {
+    let mut rng = SplitMix64::new(PILLAR_11_SEED);
+    let s: Vec<[f32; 15]> = (0..SUBSET)
+        .map(|_| {
+            let p = brownian_path_d2(&mut rng, N_STEPS);
+            signature_d2_deg3(&p, N_STEPS + 1)
+        })
+        .collect();
+
+    let mut sig_fp = 0xcbf2_9ce4_8422_2325u64;
+    for v in &s {
+        for x in v {
+            fnv(&mut sig_fp, x.to_bits() as u64);
+        }
+    }
+
+    let mut g = Array2::<f64>::zeros((SUBSET, SUBSET));
+    let mut k_fp = 0xcbf2_9ce4_8422_2325u64;
+    for i in 0..SUBSET {
+        for j in 0..SUBSET {
+            let v = sigker_hl(&s[i], &s[j]);
+            fnv(&mut k_fp, v.to_bits() as u64);
+            g[[i, j]] = v as f64;
+        }
+    }
+
+    let mean_diag = (0..SUBSET).map(|i| g[[i, i]]).sum::<f64>() / SUBSET as f64;
+    let eps = JITTER * mean_diag.abs().max(1.0);
+    let mut j = g.clone();
+    for i in 0..SUBSET {
+        j[[i, i]] += eps;
+    }
+    let chol = j.cholesky();
+    let mut l_fp = 0xcbf2_9ce4_8422_2325u64;
+    for v in chol.factor.iter() {
+        fnv(&mut l_fp, v.to_bits());
+    }
+
+    println!("signature bits  fnv1a = {sig_fp:#018x}");
+    println!("kernel bits     fnv1a = {k_fp:#018x}");
+    println!("cholesky L bits fnv1a = {l_fp:#018x}");
+    println!("cholesky info         = {}", chol.info);
+    println!("mean diag             = {mean_diag:.17e}");
+}

examples/psd_f32_vs_f64_probe.rs

@@ -0,0 +1,53 @@
+//! Pre-registration: does the truncated Pillar-11 Gram survive an f32
+//! Cholesky, or does it need f64 accumulation?
+//!
+//! `sigker_hl` returns f32. Cholesky on an f32 Gram can report `info > 0` on
+//! a matrix that IS positive semi-definite but ill-conditioned — a false
+//! alarm on the real battery. This measures the margin on the ACTUAL pool
+//! before any gate is pinned.
+use ndarray::hpc::lapack::LapackOps;
+use ndarray::hpc::pillar::signature::PILLAR_11_SEED;
+use ndarray::hpc::pillar::signature::{brownian_path_d2, sigker_hl, signature_d2_deg3};
+use ndarray::hpc::pillar::SplitMix64;
+use ndarray::Array2;
+
+const N_STEPS: usize = 50;
+
+fn sigs(n: usize) -> Vec<[f32; 15]> {
+    let mut rng = SplitMix64::new(PILLAR_11_SEED);
+    (0..n)
+        .map(|_| {
+            let p = brownian_path_d2(&mut rng, N_STEPS);
+            signature_d2_deg3(&p, N_STEPS + 1)
+        })
+        .collect()
+}
+
+fn main() {
+    println!(
+        "{:>7} {:>12} {:>12} {:>14} {:>14} {:>12}",
+        "subset", "f32 info", "f64 info", "min diag", "max |K_ij|", "cond-ish"
+    );
+    for &n in &[8usize, 12, 14, 15, 16, 17, 32, 50] {
+        let s = sigs(n);
+        let mut g32 = Array2::<f32>::zeros((n, n));
+        let mut g64 = Array2::<f64>::zeros((n, n));
+        for i in 0..n {
+            for j in 0..n {
+                let v = sigker_hl(&s[i], &s[j]);
+                g32[[i, j]] = v;
+                g64[[i, j]] = v as f64;
+            }
+        }
+        let i32_ = g32.cholesky().info;
+        let i64_ = g64.cholesky().info;
+        let min_diag = (0..n).map(|i| g64[[i, i]]).fold(f64::INFINITY, f64::min);
+        let max_off = (0..n)
+            .flat_map(|i| (0..n).filter(move |j| *j != i).map(move |j| (i, j)))
+            .map(|(i, j)| g64[[i, j]].abs())
+            .fold(0.0f64, f64::max);
+        let max_diag = (0..n).map(|i| g64[[i, i]]).fold(0.0f64, f64::max);
+        println!("{n:>7} {i32_:>12} {i64_:>12} {min_diag:>14.4e} {max_off:>14.4e} {:>12.2e}", max_diag / min_diag);
+    }
+    println!("\ninfo == 0 means positive definite; info > 0 is the failing leading minor.");
+}

examples/psd_jitter_sweep.rs

@@ -0,0 +1,101 @@
+//! Pick the PSD jitter by measurement, and prove it still discriminates.
+//!
+//! The truncated signature is a 15-dimensional feature map (1+2+4+8), so the
+//! Gram of N > 15 paths is rank-deficient BY CONSTRUCTION and plain Cholesky
+//! — which tests positive DEFINITE — must fail at leading minor 16. The
+//! property the battery actually wants is positive SEMI-definite, i.e.
+//! Cholesky with a relative jitter on the diagonal.
+//!
+//! A jitter large enough to admit the real (singular, PSD) Gram must still
+//! REJECT a genuinely indefinite one. This sweeps both together; a jitter
+//! that admits both is worthless.
+use ndarray::hpc::lapack::LapackOps;
+use ndarray::hpc::pillar::signature::PILLAR_11_SEED;
+use ndarray::hpc::pillar::signature::{brownian_path_d2, sigker_hl, signature_d2_deg3};
+use ndarray::hpc::pillar::SplitMix64;
+use ndarray::Array2;
+
+const N_STEPS: usize = 50;
+const SUBSET: usize = 50;
+
+fn real_gram(n: usize) -> Array2<f64> {
+    let mut rng = SplitMix64::new(PILLAR_11_SEED);
+    let s: Vec<[f32; 15]> = (0..n)
+        .map(|_| {
+            let p = brownian_path_d2(&mut rng, N_STEPS);
+            signature_d2_deg3(&p, N_STEPS + 1)
+        })
+        .collect();
+    let mut g = Array2::<f64>::zeros((n, n));
+    for i in 0..n {
+        for j in 0..n {
+            g[[i, j]] = sigker_hl(&s[i], &s[j]) as f64;
+        }
+    }
+    g
+}
+
+/// A genuinely indefinite matrix that PASSES both weak criteria: every
+/// diagonal positive, Cauchy-Schwarz satisfied for every pair.
+fn indefinite_gram(n: usize) -> Array2<f64> {
+    let mut g = real_gram(n);
+    let s01 = (g[[0, 0]] * g[[1, 1]]).sqrt();
+    let s02 = (g[[0, 0]] * g[[2, 2]]).sqrt();
+    let s12 = (g[[1, 1]] * g[[2, 2]]).sqrt();
+    g[[0, 1]] = -0.999 * s01;
+    g[[1, 0]] = g[[0, 1]];
+    g[[0, 2]] = 0.999 * s02;
+    g[[2, 0]] = g[[0, 2]];
+    g[[1, 2]] = 0.999 * s12;
+    g[[2, 1]] = g[[1, 2]];
+    g
+}
+
+fn chol_info(g: &Array2<f64>, jitter: f64) -> i32 {
+    let n = g.nrows();
+    let mean_diag = (0..n).map(|i| g[[i, i]]).sum::<f64>() / n as f64;
+    let eps = jitter * mean_diag.abs().max(1.0);
+    let mut j = g.clone();
+    for i in 0..n {
+        j[[i, i]] += eps;
+    }
+    j.cholesky().info
+}
+
+fn weak_criteria_hold(g: &Array2<f64>) -> (bool, bool) {
+    let n = g.nrows();
+    let diag_ok = (0..n).all(|i| g[[i, i]] > 0.0);
+    let mut cs_ok = true;
+    for i in 0..n {
+        for j in i + 1..n {
+            if g[[i, j]] * g[[i, j]] > g[[i, i]] * g[[j, j]] * 1.001 {
+                cs_ok = false;
+            }
+        }
+    }
+    (diag_ok, cs_ok)
+}
+
+fn main() {
+    let real = real_gram(SUBSET);
+    let bad = indefinite_gram(SUBSET);
+
+    let (rd, rc) = weak_criteria_hold(&real);
+    let (bd, bc) = weak_criteria_hold(&bad);
+    println!("weak criteria (diag>0, Cauchy-Schwarz):");
+    println!("  real Gram       diag={rd} cs={rc}");
+    println!("  indefinite Gram diag={bd} cs={bc}   <- both PASS, which is the point");
+
+    println!("\n{:>12} {:>14} {:>18} {:>10}", "jitter", "real info", "indefinite info", "verdict");
+    for &j in &[0.0f64, 1e-12, 1e-9, 1e-7, 1e-6, 1e-5, 1e-4, 1e-3, 1e-2, 1e-1] {
+        let r = chol_info(&real, j);
+        let b = chol_info(&bad, j);
+        let verdict = match (r == 0, b == 0) {
+            (true, false) => "USABLE",
+            (false, false) => "too tight",
+            (true, true) => "TOO LOOSE",
+            (false, true) => "impossible?",
+        };
+        println!("{j:>12.0e} {r:>14} {b:>18} {verdict:>10}");
+    }
+}