From 388d433d0bf5e7300406078a7a6119dfabd266be Mon Sep 17 00:00:00 2001 From: Claude Date: Sat, 24 Oct 2026 06:32:33 +0000 Subject: [PATCH 1/4] bgz-tensor: symmetric signed-i8 quantization probe (lab, no wire change) Isolated mathematical tests for a proposed -127..127 code space with -128 reserved as an invalid sentinel, measured against the existing FamilyGamma codec, which is not changed: - the code space: 255 values, one sentinel, one zero, odd negation; conventional i8 has no negation for -128 - saturating add over all 256x256 pairs: odd and commutative, sentinel propagates, valid inputs never produce it; conventional saturating_add breaks oddness on 16256 pairs - round-half-to-even gives symmetric ties; floor(x+0.5) does not - FamilyGamma: build() never stores -128, but encode() emits -128 for any input below the fitted range (cosine 0 in an all-positive family) and decodes it as a value; truncation shifts codes half a bin toward the mid code on both sides - four schemes over 13 synthetic distributions: symmetric keeps exact zeros and matches affine on centred data; affine beats it by about 2x RMSE on shifted and skewed data Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01G22yT6htkcdyXsihxxXdrg --- .../tests/symmetric_i8_quantization_probe.rs | 621 ++++++++++++++++++ 1 file changed, 621 insertions(+) create mode 100644 crates/bgz-tensor/tests/symmetric_i8_quantization_probe.rs diff --git a/crates/bgz-tensor/tests/symmetric_i8_quantization_probe.rs b/crates/bgz-tensor/tests/symmetric_i8_quantization_probe.rs new file mode 100644 index 000000000..c14a168fe --- /dev/null +++ b/crates/bgz-tensor/tests/symmetric_i8_quantization_probe.rs @@ -0,0 +1,621 @@ +//! Symmetric signed-i8 quantization: mathematical probe (lab only). +//! +//! Proposed code space: `-127..=127` are values, `-128` is reserved as the +//! invalid / NaN sentinel. Nothing here changes `FamilyGamma` or any wire +//! format; the existing codec is only measured. +//! +//! Hypotheses kept apart: +//! - H1 representational symmetry: zero is the exact centre of the valid +//! code range. +//! - H2 quantization symmetry: `encode(-x) == -encode(x)` under one +//! zero-centred scale and one rounding rule. +//! - H3 distributional symmetry: the data has median zero. H1 does not imply it. +//! - H4 a symmetric code removes the integer half-bin offset; it does not +//! remove calibration offsets that asymmetric data needs. +//! +//! Run with `--nocapture` to print the comparison tables. + +use bgz_tensor::fisher_z::{FamilyGamma, FisherZTable}; + +// ── the proposed codec ───────────────────────────────────────────────────── + +mod sym { + pub const SENTINEL: i8 = i8::MIN; + + /// Round-half-to-even on `x / scale`, saturating at ±127. NaN gives the + /// sentinel; ±inf saturates. + pub fn encode(x: f32, scale: f32) -> i8 { + if x.is_nan() { + return SENTINEL; + } + (x / scale).round_ties_even().clamp(-127.0, 127.0) as i8 + } + + pub fn decode(c: i8, scale: f32) -> f32 { + if c == SENTINEL { + f32::NAN + } else { + c as f32 * scale + } + } + + pub fn neg(c: i8) -> i8 { + if c == SENTINEL { + SENTINEL + } else { + -c + } + } + + /// Saturating add; either sentinel gives the sentinel. + pub fn add(a: i8, b: i8) -> i8 { + if a == SENTINEL || b == SENTINEL { + SENTINEL + } else { + (a as i16 + b as i16).clamp(-127, 127) as i8 + } + } +} + +// ── deterministic generator ──────────────────────────────────────────────── + +struct Rng(u64); +impl Rng { + fn next_u64(&mut self) -> u64 { + self.0 = self.0.wrapping_add(0x9E37_79B9_7F4A_7C15); + let mut z = self.0; + z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9); + z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB); + z ^ (z >> 31) + } + fn uniform(&mut self) -> f64 { + ((self.next_u64() >> 11) as f64 + 0.5) / (1u64 << 53) as f64 + } + fn normal(&mut self) -> f64 { + let (u, v) = (self.uniform(), self.uniform()); + (-2.0 * u.ln()).sqrt() * (std::f64::consts::TAU * v).cos() + } +} + +// ── B2: exhaustive code-space tests ──────────────────────────────────────── + +#[test] +fn proposed_code_space_has_255_values_one_sentinel_and_one_zero() { + let scale = 0.01f32; + let all: Vec = (i8::MIN..=i8::MAX).collect(); + assert_eq!(all.len(), 256); + let valid: Vec = all + .iter() + .copied() + .filter(|&c| c != sym::SENTINEL) + .collect(); + assert_eq!(valid.len(), 255); + assert!(sym::decode(sym::SENTINEL, scale).is_nan()); + assert_eq!( + valid + .iter() + .filter(|&&c| sym::decode(c, scale) == 0.0) + .count(), + 1 + ); + // H1: as many negative as positive valid codes, zero in the middle. + assert_eq!(valid.iter().filter(|&&c| c < 0).count(), 127); + assert_eq!(valid.iter().filter(|&&c| c > 0).count(), 127); + for &c in &valid { + let n = sym::neg(c); + assert_ne!( + n, + sym::SENTINEL, + "negating a valid code never hits the sentinel" + ); + assert_eq!(sym::neg(n), c); + assert_eq!( + sym::decode(n, scale), + -sym::decode(c, scale), + "decoded symmetry at {c}" + ); + assert_eq!( + sym::encode(sym::decode(c, scale), scale), + c, + "decode/encode round trip at {c}" + ); + } + assert_eq!(sym::neg(sym::SENTINEL), sym::SENTINEL); + + // Conventional i8 is not symmetric: 128 negative codes, 127 positive, + // and -128 has no negation. + assert_eq!(all.iter().filter(|&&c| c < 0).count(), 128); + assert_eq!(i8::MIN.checked_neg(), None); + assert_eq!(i8::MIN.wrapping_neg(), i8::MIN); + assert_eq!(i8::MIN.saturating_neg(), i8::MAX); + assert_eq!(i8::MIN.unsigned_abs(), 128u8, "abs(-128) does not fit i8"); +} + +#[test] +fn proposed_add_is_odd_and_commutative_over_all_256x256_pairs() { + let mut conventional_odd_failures = 0u32; + for a in i8::MIN..=i8::MAX { + for b in i8::MIN..=i8::MAX { + let s = sym::add(a, b); + assert_eq!(s, sym::add(b, a)); + if a == sym::SENTINEL || b == sym::SENTINEL { + assert_eq!(s, sym::SENTINEL, "sentinel propagates"); + continue; + } + assert_ne!(s, sym::SENTINEL, "valid inputs never produce the sentinel"); + assert_eq!( + sym::add(sym::neg(a), sym::neg(b)), + sym::neg(s), + "odd at ({a}, {b})" + ); + // Conventional saturating i8 add: -a only exists when a != -128. + if let (Some(na), Some(nb)) = (a.checked_neg(), b.checked_neg()) { + if na.saturating_add(nb).checked_neg() != Some(a.saturating_add(b)) { + conventional_odd_failures += 1; + } + } + } + } + // e.g. (-1) + (-127) saturates to -128, but 1 + 127 saturates to 127. + assert!(conventional_odd_failures > 0); + println!("conventional saturating_add oddness failures among 255x255 negatable pairs: {conventional_odd_failures}"); +} + +#[test] +fn proposed_rounding_is_symmetric_with_deterministic_ties() { + let scale = 0.01f32; + // H2: encode(-x) == -encode(x) for every tested x, including exact half codes. + let mut rng = Rng(7); + for k in -200i32..=200 { + for x in [ + (k as f32 + 0.5) * scale, + k as f32 * scale, + (rng.uniform() as f32 - 0.5) * 3.0, + ] { + assert_eq!( + sym::encode(-x, scale), + sym::neg(sym::encode(x, scale)), + "x = {x}" + ); + } + } + // Ties go to the even code, so ±0.5, ±2.5 land on 0 and ±2. + let tie = |v: f32| sym::encode(v * scale, scale); + assert_eq!( + [ + tie(0.5), + tie(-0.5), + tie(1.5), + tie(-1.5), + tie(2.5), + tie(-2.5) + ], + [0, 0, 2, -2, 2, -2] + ); + assert_eq!(sym::encode(f32::INFINITY, scale), 127); + assert_eq!(sym::encode(f32::NEG_INFINITY, scale), -127); + assert_eq!(sym::encode(f32::NAN, scale), sym::SENTINEL); + assert_eq!(sym::encode(-0.0, scale), 0); + // The common `floor(x + 0.5)` rule is not symmetric: +0.5 -> 1, -0.5 -> 0. + let floor_half = |v: f32| (v + 0.5).floor(); + assert_ne!(floor_half(-0.5), -floor_half(0.5)); +} + +// ── B4: the existing FamilyGamma contract, measured ──────────────────────── + +#[test] +fn family_gamma_emits_minus_128_only_below_its_calibrated_range() { + let g = FamilyGamma { + z_min: -0.5, + z_range: 1.0, + }; + // In range the code space is -127..=127. + for i in 0..=10_000 { + let z = g.z_min + g.z_range * (i as f32 / 10_000.0); + let c = g.encode(z.tanh()); + assert!(c >= -127, "in-range z = {z} encoded to {c}"); + } + // One full code below z_min clamps to -128, which decodes as a value + // (z_min - z_range/254), not as an invalid marker. + let below = (g.z_min - 2.0 * g.z_range / 254.0).tanh(); + assert_eq!(g.encode(below), -128); + assert!(g.decode(-128).is_finite()); + // NaN input encodes to code 0, the family's mid-z. + assert_eq!(g.encode(f32::NAN), 0); +} + +#[test] +fn fisher_z_table_build_never_stores_minus_128() { + // build() fits z_min over the table's own off-diagonal cosines and the + // diagonal is the maximum, so every stored z is >= z_min: no entry can + // clamp below -127. Property-checked over random tables. + let mut rng = Rng(1); + let mut min_seen = i8::MAX; + for trial in 0..300 { + let k = 2 + (rng.next_u64() % 40) as usize; + let d = 2 + (rng.next_u64() % 31) as usize; + let rows: Vec> = (0..k) + .map(|_| (0..d).map(|_| rng.normal() as f32).collect()) + .collect(); + let t = FisherZTable::build(&rows, k); + let m = *t.entries.iter().min().unwrap(); + min_seen = min_seen.min(m); + assert!(m >= -127, "trial {trial}: k={k} d={d} stored {m}"); + } + assert_eq!(min_seen, -127, "the lowest off-diagonal reaches -127"); +} + +#[test] +fn family_gamma_truncation_biases_codes_toward_the_mid_code() { + // encode uses `as i8`, which truncates toward zero, so every code moves + // toward code 0 (the family's mid-z) by up to one whole bin; decode adds + // no half-bin back. Measured on the exact real-valued code. + let g = FamilyGamma { + z_min: -1.0, + z_range: 2.0, + }; + let n = 200_000; + let (mut pos_sum, mut pos_n, mut neg_sum, mut neg_n) = (0.0f64, 0u32, 0.0f64, 0u32); + let (mut round_sum, mut round_n) = (0.0f64, 0u32); + for i in 0..n { + let z = g.z_min + g.z_range * ((i as f64 + 0.5) / n as f64) as f32; + let ideal = ((z - g.z_min) / g.z_range) as f64 * 254.0 - 127.0; + let got = g.encode(z.tanh()) as f64; + let err = got - ideal; + if ideal > 0.5 { + pos_sum += err; + pos_n += 1; + } else if ideal < -0.5 { + neg_sum += err; + neg_n += 1; + } + round_sum += ideal.round_ties_even() - ideal; + round_n += 1; + } + let (pos, neg, round) = ( + pos_sum / pos_n as f64, + neg_sum / neg_n as f64, + round_sum / round_n as f64, + ); + println!("FamilyGamma mean code error: ideal>0 {pos:+.4}, ideal<0 {neg:+.4}; round-to-nearest {round:+.4}"); + assert!( + (pos + 0.5).abs() < 0.02, + "positive side biased down by half a bin: {pos}" + ); + assert!( + (neg - 0.5).abs() < 0.02, + "negative side biased up by half a bin: {neg}" + ); + assert!(round.abs() < 0.01); +} + +#[test] +fn family_gamma_code_zero_is_the_mid_z_not_cosine_zero() { + // H1 does not give H3: for a family whose cosines are all positive, + // code 0 (also what lookup_i8 returns for an unknown centroid) decodes + // to a clearly positive cosine. + let g = FamilyGamma::from_cosines(&[0.1, 0.3, 0.5, 0.7, 0.9]); + let c0 = g.decode(0); + println!( + "all-positive family: code 0 decodes to cosine {c0:.4}; cosine 0 encodes to {}", + g.encode(0.0) + ); + assert!(c0 > 0.4); + // Cosine 0 lies below this family's calibrated range, so the clamp emits + // -128: an ordinary input lands on the code a sentinel scheme would reserve. + assert_eq!(g.encode(0.0), -128); +} + +// ── B3: scheme comparison over synthetic distributions ───────────────────── + +#[derive(Clone, Copy, Debug)] +enum Scheme { + /// -127..=127, scale = max|x|/127, round half to even. + Symmetric, + /// -128..=127 affine min-max, integer zero point, round half to even. + Affine, + /// 0..=127 affine min-max (7-bit unsigned). + Unsigned7, + /// FamilyGamma's mapping applied to the values: [min,max] -> -127..127 + /// by truncation, clamped to -128..127. + LegacyTrunc, +} + +fn quantize(scheme: Scheme, xs: &[f64]) -> (Vec, Vec) { + let lo = xs.iter().cloned().fold(f64::INFINITY, f64::min); + let hi = xs.iter().cloned().fold(f64::NEG_INFINITY, f64::max); + let amax = xs.iter().map(|v| v.abs()).fold(0.0, f64::max).max(1e-30); + let range = (hi - lo).max(1e-30); + let mut dec = Vec::with_capacity(xs.len()); + let mut codes = Vec::with_capacity(xs.len()); + for &x in xs { + let (c, d) = match scheme { + Scheme::Symmetric => { + let s = amax / 127.0; + let c = (x / s).round_ties_even().clamp(-127.0, 127.0); + (c as i32, c * s) + } + Scheme::Affine => { + let s = range / 255.0; + let zp = (-lo / s).round_ties_even() - 128.0; + let c = ((x / s).round_ties_even() + zp).clamp(-128.0, 127.0); + (c as i32, (c - zp) * s) + } + Scheme::Unsigned7 => { + let s = range / 127.0; + let c = ((x - lo) / s).round_ties_even().clamp(0.0, 127.0); + (c as i32, c * s + lo) + } + Scheme::LegacyTrunc => { + let v = ((x - lo) / range * 254.0 - 127.0) + .clamp(-128.0, 127.0) + .trunc(); + (v as i32, (v + 127.0) / 254.0 * range + lo) + } + }; + codes.push(c); + dec.push(d); + } + (dec, codes) +} + +fn ranks(v: &[f64]) -> Vec { + let mut idx: Vec = (0..v.len()).collect(); + idx.sort_by(|&a, &b| v[a].partial_cmp(&v[b]).unwrap()); + let mut r = vec![0.0; v.len()]; + let mut i = 0; + while i < idx.len() { + let mut j = i; + while j + 1 < idx.len() && v[idx[j + 1]] == v[idx[i]] { + j += 1; + } + let avg = (i + j) as f64 / 2.0; + for &k in &idx[i..=j] { + r[k] = avg; + } + i = j + 1; + } + r +} + +fn spearman(a: &[f64], b: &[f64]) -> f64 { + let (ra, rb) = (ranks(a), ranks(b)); + let n = a.len() as f64; + let (ma, mb) = (ra.iter().sum::() / n, rb.iter().sum::() / n); + let (mut c, mut va, mut vb) = (0.0, 0.0, 0.0); + for i in 0..a.len() { + c += (ra[i] - ma) * (rb[i] - mb); + va += (ra[i] - ma).powi(2); + vb += (rb[i] - mb).powi(2); + } + c / (va * vb).sqrt() +} + +#[derive(Debug, Clone, Copy)] +struct Metrics { + bias: f64, + median: f64, + mae: f64, + rmse: f64, + max_abs: f64, + zero_kept: f64, + sat: f64, + bins: usize, + entropy: f64, + rho: f64, + near_zero_mae: f64, +} + +fn measure(scheme: Scheme, xs: &[f64]) -> Metrics { + let (dec, codes) = quantize(scheme, xs); + let n = xs.len() as f64; + let mut errs: Vec = dec.iter().zip(xs).map(|(d, x)| d - x).collect(); + let bias = errs.iter().sum::() / n; + let mae = errs.iter().map(|e| e.abs()).sum::() / n; + let rmse = (errs.iter().map(|e| e * e).sum::() / n).sqrt(); + let max_abs = errs.iter().map(|e| e.abs()).fold(0.0, f64::max); + errs.sort_by(|a, b| a.partial_cmp(b).unwrap()); + let median = errs[errs.len() / 2]; + let zeros: Vec = (0..xs.len()).filter(|&i| xs[i] == 0.0).collect(); + let zero_kept = if zeros.is_empty() { + f64::NAN + } else { + zeros.iter().filter(|&&i| dec[i] == 0.0).count() as f64 / zeros.len() as f64 + }; + let (cmin, cmax) = match scheme { + Scheme::Symmetric => (-127, 127), + Scheme::Affine | Scheme::LegacyTrunc => (-128, 127), + Scheme::Unsigned7 => (0, 127), + }; + let sat = codes.iter().filter(|&&c| c == cmin || c == cmax).count() as f64 / n; + let mut hist = std::collections::HashMap::new(); + for &c in &codes { + *hist.entry(c).or_insert(0u32) += 1; + } + let entropy = hist + .values() + .map(|&k| { + let p = k as f64 / n; + -p * p.log2() + }) + .sum(); + let amax = xs.iter().map(|v| v.abs()).fold(0.0, f64::max); + let near: Vec = (0..xs.len()) + .filter(|&i| xs[i].abs() < 0.05 * amax) + .map(|i| (dec[i] - xs[i]).abs()) + .collect(); + let near_zero_mae = if near.is_empty() { + f64::NAN + } else { + near.iter().sum::() / near.len() as f64 + }; + Metrics { + bias, + median, + mae, + rmse, + max_abs, + zero_kept, + sat, + bins: hist.len(), + entropy, + rho: spearman(xs, &dec), + near_zero_mae, + } +} + +fn distributions() -> Vec<(&'static str, Vec)> { + let n = 20_000; + let mut r = Rng(20_261_024); + let mut gen = + |f: &mut dyn FnMut(&mut Rng) -> f64| (0..n).map(|_| f(&mut r)).collect::>(); + let sparse = |p0: f64| move |r: &mut Rng| if r.uniform() < p0 { 0.0 } else { r.normal() }; + vec![ + ("uniform symmetric", gen(&mut |r| r.uniform() * 2.0 - 1.0)), + ("gaussian centred", gen(&mut |r| r.normal())), + ("gaussian shifted +4σ", gen(&mut |r| r.normal() + 4.0)), + ( + "laplace", + gen(&mut |r| { + let u = r.uniform() - 0.5; + -u.signum() * (1.0 - 2.0 * u.abs()).ln() + }), + ), + ( + "student-t ν=3", + gen(&mut |r| { + let z = r.normal(); + let chi = (0..3).map(|_| r.normal().powi(2)).sum::(); + z / (chi / 3.0).sqrt() + }), + ), + ("exponential (skewed)", gen(&mut |r| -r.uniform().ln())), + ( + "bimodal -2 / +3", + gen(&mut |r| { + if r.uniform() < 0.5 { + -2.0 + 0.4 * r.normal() + } else { + 3.0 + 0.4 * r.normal() + } + }), + ), + ( + "multimodal 4", + gen(&mut |r| [-3.0, -1.0, 0.5, 2.5][(r.next_u64() % 4) as usize] + 0.2 * r.normal()), + ), + ("sparse 90% zero", gen(&mut sparse(0.90))), + ("sparse 99% zero", gen(&mut sparse(0.99))), + ("sparse 99.9% zero", gen(&mut sparse(0.999))), + ( + "sparse asym (+exp / -tiny)", + gen(&mut |r| { + let u = r.uniform(); + if u < 0.90 { + 0.0 + } else if u < 0.97 { + -r.uniform().ln() + } else { + -0.05 * r.uniform() + } + }), + ), + ( + "on symmetric half-codes", + gen(&mut |r| { + let k = (r.next_u64() % 254) as f64 - 127.0; + (k + 0.5).clamp(-127.0, 127.0) / 127.0 + }), + ), + ] +} + +#[test] +fn scheme_comparison_over_synthetic_distributions() { + let schemes = [ + Scheme::Symmetric, + Scheme::Affine, + Scheme::Unsigned7, + Scheme::LegacyTrunc, + ]; + let mut table = std::collections::HashMap::new(); + println!( + "{:<28} {:<11} {:>10} {:>10} {:>9} {:>9} {:>9} {:>6} {:>6} {:>4} {:>6} {:>8} {:>9}", + "distribution", + "scheme", + "bias", + "median", + "MAE", + "RMSE", + "max|e|", + "zero", + "sat", + "bins", + "H", + "rho", + "nearzero" + ); + for (name, xs) in distributions() { + for s in schemes { + let m = measure(s, &xs); + println!( + "{:<28} {:<11} {:>+10.2e} {:>+10.2e} {:>9.2e} {:>9.2e} {:>9.2e} {:>6.3} {:>6.4} {:>4} {:>6.3} {:>8.6} {:>9.2e}", + name, format!("{s:?}"), m.bias, m.median, m.mae, m.rmse, m.max_abs, m.zero_kept, m.sat, m.bins, m.entropy, m.rho, m.near_zero_mae + ); + table.insert((name, format!("{s:?}")), m); + } + } + let get = |d: &'static str, s: &str| table[&(d, s.to_string())]; + + // Symmetric keeps exact zeros; legacy truncation and 7-bit unsigned do not. + for d in ["sparse 90% zero", "sparse 99% zero", "sparse 99.9% zero"] { + assert_eq!(get(d, "Symmetric").zero_kept, 1.0, "{d}"); + assert!(get(d, "LegacyTrunc").zero_kept < 1.0, "{d}"); + } + // Zero-centred data: symmetric costs at most the one reserved code + // against the full affine range (255 vs 256 levels over a symmetric range). + let (s, a) = ( + get("gaussian centred", "Symmetric").rmse, + get("gaussian centred", "Affine").rmse, + ); + assert!( + s <= a * 1.10, + "centred gaussian: symmetric {s} vs affine {a}" + ); + // Falsifier for "symmetric is universally better": shifted data uses + // under half the symmetric code space, and affine wins clearly. + let (s, a) = ( + get("gaussian shifted +4σ", "Symmetric"), + get("gaussian shifted +4σ", "Affine"), + ); + assert!( + a.rmse * 1.5 < s.rmse, + "shifted: affine {} vs symmetric {}", + a.rmse, + s.rmse + ); + assert!( + s.bins < 128, + "shifted data occupies {} symmetric bins", + s.bins + ); + let (s, a) = ( + get("exponential (skewed)", "Symmetric"), + get("exponential (skewed)", "Affine"), + ); + assert!( + a.rmse * 1.5 < s.rmse, + "exponential: affine {} vs symmetric {}", + a.rmse, + s.rmse + ); + // Legacy truncation carries a signed bias the rounding schemes do not. + let lt = get("uniform symmetric", "LegacyTrunc"); + let sy = get("uniform symmetric", "Symmetric"); + assert!( + lt.mae > sy.mae * 1.5, + "truncation MAE {} vs symmetric {}", + lt.mae, + sy.mae + ); +} From 5bef8ef7f6201e8786339dfa7b7ea81f6d0aae8d Mon Sep 17 00:00:00 2001 From: Claude Date: Sat, 24 Oct 2026 06:34:53 +0000 Subject: [PATCH 2/4] bgz-tensor: Fisher-Z -> i8 -> tanh error decomposition example (lab) Splits FamilyGamma's round-trip error into stages that telescope exactly to the codec total (asserted per pair): input/range calibration, rounding to a code, the codec's as-i8 truncation, f32 decode, and vstanh in place of f32::tanh on the same z. On five seeded synthetic families (k=256, d=64): truncation is the largest and the only biased stage (mean |e| 2.0e-3 against 1.0e-3 for rounding in the centred family); the tanh kernel contributes about 3e-9 mean, 1.2e-7 max (<= 2 ULP from f32::tanh). A gamma fitted on 64 rows and applied to held-out rows clamps 80 pairs (27 to -128) and calibration then dominates the max error (3.6e-2). A family whose cosines all exceed the 0.9999 input clamp collapses to one z. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01G22yT6htkcdyXsihxxXdrg --- .../examples/fisher_z_error_decomposition.rs | 237 ++++++++++++++++++ 1 file changed, 237 insertions(+) create mode 100644 crates/bgz-tensor/examples/fisher_z_error_decomposition.rs diff --git a/crates/bgz-tensor/examples/fisher_z_error_decomposition.rs b/crates/bgz-tensor/examples/fisher_z_error_decomposition.rs new file mode 100644 index 000000000..5ffeb4f2a --- /dev/null +++ b/crates/bgz-tensor/examples/fisher_z_error_decomposition.rs @@ -0,0 +1,237 @@ +//! Fisher-Z -> i8 -> tanh: where the error comes from (lab probe). +//! +//! For every off-diagonal pair cosine `c` of a family, the round trip +//! `FamilyGamma::decode(FamilyGamma::encode(c))` is split into stages +//! that telescope exactly to the total error: +//! +//! | stage | what changes | error measured against | +//! |---|---|---| +//! | calibration | z clamped to the codable range `[z_min - zr/254, z_min + zr]` and the input to +-0.9999 | the true cosine | +//! | quantization | z rounded to the nearest code (round half to even), decoded in f64 | the calibration-only value | +//! | truncation | the codec's `as i8` truncation instead of rounding | the rounded value | +//! | f32 decode | the codec's f32 decode with `f32::tanh` | the same code decoded in f64 | +//! | tanh kernel | `ndarray::hpc::vml::vstanh` instead of `f32::tanh` on the same f32 z | the codec's f32 decode | +//! +//! `vstanh` is never on the codec's path today; the last row says what it +//! would cost if the decode were batched through it. +//! +//! Families are synthetic and seeded. Real GGUF families are listed in the +//! validation manifest and NOT RUN here. +//! +//! ```text +//! cargo run --release -p bgz-tensor --example fisher_z_error_decomposition +//! ``` + +use bgz_tensor::fisher_z::FamilyGamma; +use ndarray::hpc::vml::vstanh; +use ndarray::Array1; + +struct Rng(u64); +impl Rng { + fn next_u64(&mut self) -> u64 { + self.0 = self.0.wrapping_add(0x9E37_79B9_7F4A_7C15); + let mut z = self.0; + z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9); + z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB); + z ^ (z >> 31) + } + fn uniform(&mut self) -> f64 { + ((self.next_u64() >> 11) as f64 + 0.5) / (1u64 << 53) as f64 + } + fn normal(&mut self) -> f64 { + let (u, v) = (self.uniform(), self.uniform()); + (-2.0 * u.ln()).sqrt() * (std::f64::consts::TAU * v).cos() + } +} + +/// Rows `shared * s + (1 - shared) * noise`: `shared = 0` gives cosines +/// centred on zero, larger values give an all-positive family. +fn family(seed: u64, k: usize, d: usize, shared: f64) -> Vec> { + let mut r = Rng(seed); + let s: Vec = (0..d).map(|_| r.normal()).collect(); + (0..k) + .map(|_| { + (0..d) + .map(|j| (shared * s[j] + (1.0 - shared) * r.normal()) as f32) + .collect() + }) + .collect() +} + +fn cosine(a: &[f32], b: &[f32]) -> f64 { + let (mut ab, mut aa, mut bb) = (0.0f64, 0.0f64, 0.0f64); + for (&x, &y) in a.iter().zip(b) { + let (x, y) = (x as f64, y as f64); + ab += x * y; + aa += x * x; + bb += y * y; + } + ab / (aa * bb).sqrt() +} + +fn pair_cosines(rows: &[Vec]) -> Vec { + let mut out = Vec::new(); + for i in 0..rows.len() { + for j in i + 1..rows.len() { + out.push(cosine(&rows[i], &rows[j])); + } + } + out +} + +#[derive(Default)] +struct Stat { + n: usize, + sum: f64, + sum_abs: f64, + max_abs: f64, +} +impl Stat { + fn add(&mut self, e: f64) { + self.n += 1; + self.sum += e; + self.sum_abs += e.abs(); + self.max_abs = self.max_abs.max(e.abs()); + } +} + +fn ulp_distance(a: f32, b: f32) -> u32 { + let key = |x: f32| { + let i = x.to_bits() as i32; + if i < 0 { + i32::MIN - i + } else { + i + } + }; + (key(a) as i64 - key(b) as i64).unsigned_abs() as u32 +} + +/// `fit` calibrates the gamma, `eval` is the cosine set measured. +fn decompose(name: &str, fit: &[f64], eval: &[f64]) { + let fit_f32: Vec = fit.iter().map(|&c| c as f32).collect(); + let g = FamilyGamma::from_cosines(&fit_f32); + let (z_min, zr) = (g.z_min as f64, g.z_range as f64); + let z_lo = z_min - zr / 254.0; + let z_hi = z_min + zr; + + // The codec's f32 z for each pair, decoded once by f32::tanh (scalar, + // the codec) and once by vstanh (batched). + let codes: Vec = eval.iter().map(|&c| g.encode(c as f32)).collect(); + let z_f32 = Array1::from_iter(codes.iter().map(|&q| { + let normalized = (q as f32 + 127.0) / 254.0; + normalized * g.z_range + g.z_min + })); + let mut vs = Array1::::zeros(z_f32.len()); + vstanh(z_f32.view(), vs.view_mut()); + + let mut st: [Stat; 6] = Default::default(); + let mut max_tanh_ulp = 0u32; + let mut code_min = 0u8; + let mut clamped = 0usize; + for (i, &c) in eval.iter().enumerate() { + let cf = c as f32 as f64; + let z = cf.clamp(-0.9999, 0.9999).atanh(); + if !(z_lo..=z_hi).contains(&z) { + clamped += 1; + } + let s_cal = z.clamp(z_lo, z_hi).tanh(); + let ideal = (z.clamp(z_lo, z_hi) - z_min) / zr * 254.0 - 127.0; + let rounded = ideal.round_ties_even().clamp(-128.0, 127.0); + let s_round = ((rounded + 127.0) / 254.0 * zr + z_min).tanh(); + let q = codes[i]; + if q == -128 { + code_min += 1; + } + let s_trunc = ((q as f64 + 127.0) / 254.0 * zr + z_min).tanh(); + let s_codec = g.decode(q) as f64; + let s_vs = vs[i] as f64; + max_tanh_ulp = max_tanh_ulp.max(ulp_distance(vs[i], g.decode(q))); + + st[0].add(s_cal - c); + st[1].add(s_round - s_cal); + st[2].add(s_trunc - s_round); + st[3].add(s_codec - s_trunc); + st[4].add(s_vs - s_codec); + st[5].add(s_codec - c); + // Telescoping check: the stages sum to the codec's total error. + let parts = (s_cal - c) + (s_round - s_cal) + (s_trunc - s_round) + (s_codec - s_trunc); + assert!((parts - (s_codec - c)).abs() < 1e-12); + } + + let lo = eval.iter().cloned().fold(f64::INFINITY, f64::min); + let hi = eval.iter().cloned().fold(f64::NEG_INFINITY, f64::max); + println!( + "\n{name}: n={} cosine range [{lo:+.4}, {hi:+.4}] z_min={:+.4} z_range={:.4} outside codable z: {clamped} code -128: {code_min}", + eval.len(), + g.z_min, + g.z_range + ); + println!( + " {:<14} {:>12} {:>12} {:>12}", + "stage", "bias", "mean |e|", "max |e|" + ); + let names = [ + "calibration", + "quantization", + "truncation", + "f32 decode", + "tanh kernel", + "TOTAL codec", + ]; + for (n, s) in names.iter().zip(&st) { + println!( + " {:<14} {:>+12.3e} {:>12.3e} {:>12.3e}", + n, + s.sum / s.n as f64, + s.sum_abs / s.n as f64, + s.max_abs + ); + } + println!(" vstanh vs f32::tanh on the codec's z: max {max_tanh_ulp} ULP"); +} + +fn main() { + println!( + "fisher_z_error_decomposition avx512f={}", + cfg!(target_feature = "avx512f") + ); + let (k, d) = (256, 64); + + let centred = pair_cosines(&family(1, k, d, 0.0)); + decompose("centred family (shared=0), in-sample", ¢red, ¢red); + + let positive = pair_cosines(&family(2, k, d, 0.6)); + decompose( + "all-positive family (shared=0.6), in-sample", + &positive, + &positive, + ); + + // Held out: the gamma is fitted on one draw and applied to another draw + // of the same generator, so some pairs fall outside the fitted range. + let fit = pair_cosines(&family(3, 64, d, 0.6)); + let held = pair_cosines(&family(4, k, d, 0.6)); + decompose( + "all-positive, gamma fitted on 64 rows, eval on 256 held-out rows", + &fit, + &held, + ); + + // Cosines close to 1, where atanh is steep. + let near = pair_cosines(&family(5, k, d, 0.9)); + decompose( + "high-similarity family (shared=0.9), in-sample", + &near, + &near, + ); + + // Near-duplicate rows: every cosine exceeds the 0.9999 input clamp, so + // the whole family collapses onto one z and z_range hits its 1e-10 floor. + let dup = pair_cosines(&family(6, k, d, 0.995)); + decompose( + "near-duplicate family (shared=0.995), in-sample", + &dup, + &dup, + ); +} From 13922e78620e9398aa0f1209a6dd64569347ed11 Mon Sep 17 00:00:00 2001 From: Claude Date: Sat, 24 Oct 2026 06:36:37 +0000 Subject: [PATCH 3/4] board: tanh precision + symmetric i8 findings; real-data manifest NOT RUN Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01G22yT6htkcdyXsihxxXdrg --- .../2026-10-24-tanh-precision-symmetric-i8.md | 155 ++++++++++++++++++ .claude/board/entries/README.md | 3 +- 2 files changed, 157 insertions(+), 1 deletion(-) create mode 100644 .claude/board/entries/2026-10-24-tanh-precision-symmetric-i8.md diff --git a/.claude/board/entries/2026-10-24-tanh-precision-symmetric-i8.md b/.claude/board/entries/2026-10-24-tanh-precision-symmetric-i8.md new file mode 100644 index 000000000..1b4f23585 --- /dev/null +++ b/.claude/board/entries/2026-10-24-tanh-precision-symmetric-i8.md @@ -0,0 +1,155 @@ +# 2026-10-24: f32 tanh precision and symmetric signed-i8 quantization + +**Status:** MEASURED (synthetic) for A, B and C. D (real model data) is +**NOT RUN**. Nothing here changes a wire format or an existing decoder. + +## A — `vstanh` (ndarray PR #369, branch `fix/vstanh-small-arg`, head `4f6c132`) + +The old kernel computed `2·sigmoid(2x) − 1`. That cancels: every |x| below +about 1e-7 came back as 0, −0.0 came back as +0.0, and the error reached +868,220,004 ULP. Measured over all 2³² f32 inputs on AVX-512: +1,728,053,246 nonzero inputs zeroed and 975,101,106 odd-symmetry violations. + +The new kernel works on |x|: + +- |x| < 0.625: `a + a·z·P(z)`, z = a², Cephes tanhf coefficients; +- otherwise `1 − 2/(exp(2a) + 1)`; +- the sign is restored by a float select, and x == 0 returns x, so −0.0 is kept. + +`simd_exp_f32` was 39 ULP off (a degree-5 Taylor series). It is now a Cephes +expf minimax with unfused `*` and `+`, so it stays bit-equal across backends, +and measures 1 ULP over 4 M points on [−87, 88.5]. + +Exhaustive result over 2³² inputs (AVX-512): + +- max 1 ULP in every magnitude band against f64 `tanh` rounded once; +- 0 zeroed, 0 odd-symmetry violations, 0 values outside [−1, 1], 0 monotonicity violations; +- NaN preserved; 1 flip-risk input, which agrees with mpmath; +- the oracle matches mpmath correct rounding on 200,000/200,000 sampled inputs. + +Measured in the same run: `vstanh` 3.06 ns/elem, the old kernel 2.46, +scalar `f32::tanh` 11.66. + +The probe gives the same numbers on v4, v3, the x86-64 baseline, and NEON +under qemu. wasm parity is checked for exp only. The ndarray gates pass (2590 +lib tests, clippy, fmt). Disable runs D1 to D5 each turn the gate red. + +Not run: tanh under the wasm harness, the nightly-simd arm, and the +exhaustive sweep on the non-AVX-512 tiers. + +## B — symmetric signed i8 (`crates/bgz-tensor/tests/symmetric_i8_quantization_probe.rs`) + +The proposed code space has −127..127 as values and −128 as an invalid or +NaN sentinel. Values are encoded with round-half-to-even on `x/scale`, +saturating at ±127. The probe is pure math and does not touch `FamilyGamma`. + +| hypothesis | result | +|---|---| +| H1: zero is the exact centre of the valid range | holds for the proposal (127 values each side). Fails for conventional i8: 128 vs 127, and −128 has no negation (`checked_neg` is None, `saturating_neg` gives 127, `unsigned_abs` gives 128) | +| H2: `encode(−x) == −encode(x)` | holds for the proposal over ±200 codes, the half-codes and random x. Fails for `floor(x + 0.5)` and for any affine zero point | +| H3: data is median-zero | is a property of the data, not the codec. Falsified for every non-centred family below | +| H4: a symmetric code removes the half-bin offset, not calibration offset | holds. On shifted or skewed data affine beats symmetric by about 2× RMSE (gaussian +4σ: 1.06e-2 vs 2.12e-2; exponential: 1.34e-2 vs 2.68e-2), and symmetric uses fewer than 128 of its bins there | + +Over all 256×256 code pairs, the sentinel-propagating saturating add is +commutative and odd, and never produces −128 from valid inputs. +Conventional `saturating_add` breaks oddness on 16,256 of the 255×255 +negatable pairs. + +Distribution table (13 seeded distributions × 4 schemes; metrics are bias, +median, MAE, RMSE, max error, zero preservation, saturation, bins, entropy, +Spearman and near-zero error): see the test's `--nocapture` output. +Highlights: + +- **Symmetric** keeps exact zeros on every sparse set (90 / 99 / 99.9 % zero). + On centred data it is within 1 % of affine RMSE (9.32e-3 vs 9.23e-3). +- **Affine** wins on shifted, skewed and bimodal data. It also keeps zero + exactly when zero lies inside the range, because the zero point is integer. +- **Unsigned 7-bit** and **legacy truncation** both lose every exact zero. + Truncation adds a signed bias of up to +8.7e-2 on heavy tails, and doubles + MAE on uniform data. + +Disable runs, all red: round half away from zero in place of ties-to-even; +no sentinel propagation in add; add clamped to −128. + +### B4 — the FamilyGamma wire contract + +- `encode` = `(atanh(clamp(c, ±0.9999)) − z_min)/z_range·254 − 127`, + clamped to −128..127, then `as i8`, which truncates toward zero. + `decode` = `tanh((c + 127)/254·z_range + z_min)`. +- **−128 is a valid saturated value code.** `FisherZTable::build` never + stores it, because the fit covers the table's own off-diagonal and the + diagonal clamps to +127 (checked over 300 random tables). But `encode` + emits it for every input below the fitted range: cosine 0 in an + all-positive family, or a held-out pair. It decodes to the finite value + `tanh(z_min − z_range/254)`. +- Code 0 is the family's mid-z, not cosine 0. In an all-positive family it + decodes to 0.656. `lookup_i8` returns 0 for an unknown centroid, which is + therefore not neutral in such a family. +- `as i8` truncation moves codes half a bin toward code 0 on both sides + (measured −0.501 above zero and +0.501 below). +- The 8-byte gamma plus k² raw i8 carries no version byte. It is written + durably: thinking-engine `tts_17b_hhtld_encode.rs:676-680` stores + `{role}.fisher_z` and `{role}.fisher_z_gamma` in safetensors, and + `tts_17b_hhtld_decode.rs:298-300` reads them back. + `palette256_reference_contract.rs` pins `to_bytes` byte for byte. + +**Conclusion:** a −128 sentinel cannot be introduced in place. Doing so would +change the meaning of a byte that existing artifacts may contain (codes +encoded below the range). Any symmetric codec needs a new, versioned format +and leaves `FamilyGamma` and its decoder untouched. Not proposed for merge +here; this is a decision for the format's owner. + +Other i8 conventions seen in this workspace (not changed): + +- thinking-engine `signed_engine` and `tensor_bridge` round and clamp to −128..127; +- `l4.rs` maps `i8::MIN` to magnitude 127; +- bgz-tensor `zipper`, `slot_l`, `matryoshka` and `xor_adaptive` clamp to −127..127. + +The workspace has no single i8 convention. + +## C — error decomposition (`crates/bgz-tensor/examples/fisher_z_error_decomposition.rs`) + +Each pair's round-trip error is split into stages that sum exactly to the +codec total; this is asserted per pair. Mean |e|, k = 256, d = 64: + +| family | calibration | rounding | truncation (extra) | f32 decode | tanh kernel (vstanh − f32::tanh) | total | +|---|---|---|---|---|---|---| +| centred | 2.2e-9 | 1.0e-3 | 2.0e-3 | 1.5e-8 | 3.2e-9 (≤2 ULP) | 2.0e-3 | +| all-positive | 1.5e-8 | 4.7e-4 | 9.3e-4 | 3.4e-8 | 3.1e-8 (≤2 ULP) | 9.3e-4 | +| held-out gamma | 2.6e-5 (max **3.6e-2**) | 3.2e-4 | 6.5e-4 | 3.1e-8 | 3.0e-8 | 6.7e-4 (max 3.6e-2) | +| high similarity | 1.5e-8 | 1.7e-5 | 3.4e-5 | 1.5e-8 | 0 | 3.4e-5 | +| near-duplicate | 7.3e-5 | 0 | 0 | 0 | 0 | 7.3e-5 | + +The ranking is: calibration on out-of-range inputs, then truncation, then +quantization, then f32 arithmetic and the tanh kernel. The tanh kernel sits +five orders of magnitude below quantization. `vstanh` is not on the codec's +path today, and batching the decode through it would cost nothing measurable. + +The held-out case is the realistic risk: 80 of 32,640 pairs fall outside the +fitted range, 27 of them land on −128, and none is flagged. When every +cosine of a family exceeds 0.9999, the family collapses to one z. + +## D — real-data validation manifest: **NOT RUN** + +| # | input | command | measures | status | +|---|---|---|---|---| +| D1 | jina-v3 Q8_0 GGUF (`gaianet/jina-embeddings-v3-GGUF`, snapshot `d7d998aab1a7f2ea0aa256d8b3f035cbd0af682a`, the default path in `bgz-tensor/examples/gguf_families.rs`) | build `FisherZTable` per role group from real representatives, then the C decomposition | stage errors on real families | NOT RUN — file not on disk; record its sha256 when fetched | +| D2 | the same GGUF | held-out split per group (fit on half the rows, evaluate on the rest) | how often real data reaches −128 or the clamp | NOT RUN | +| D3 | existing thinking-engine `*.fisher_z` safetensors artifacts | histogram of stored codes | whether any durable artifact holds −128 | NOT RUN — no artifact on disk | +| D4 | Jina v5 safetensors (the workspace ground truth) | same as D1 | the same on a second model | NOT RUN | +| D5 | ONNX Runtime, offline | `Tanh` and `Exp` operator tests against `vstanh` and `vsexp` | operator conformance | NOT RUN — blocked on ndarray #369 merging | + +Each run must record: the file sha256, the git SHA of lance-graph and +ndarray, the target-cpu the arm printed, and the exact command. + +## Proposed handoff (not posted to EXO-50) + +1. Merge ndarray #369 (Workstream A). +2. Then the lance-graph-onnx `Tanh` → `vstanh` and `Exp` → `vsexp` alias, + removing the lib.rs:35-36 sentence that calls them approximate, and + adding operator tests (Stage 4). +3. The decision to make: whether a versioned symmetric Fisher-Z format is + wanted. If yes, it gets a new format tag and its own reader; + `FamilyGamma` and its decoder stay as they are. +4. Separately: report out-of-range codes instead of clamping silently, and + run D1–D3. diff --git a/.claude/board/entries/README.md b/.claude/board/entries/README.md index 9a3a2b2de..b9a2416de 100644 --- a/.claude/board/entries/README.md +++ b/.claude/board/entries/README.md @@ -25,10 +25,11 @@ index row, (3) no duplicate entry id. Checks 1 and 2 are deliberately opposite directions; the stranding this convention prevents shows up in exactly one of them, never both. -462 entries, 2026-08-06 .. 2026-10-23. +463 entries, 2026-08-06 .. 2026-10-24. | date | entry id | finding | file | |---|---|---|---| +| 2026-10-24 | `tanh-precision-symmetric-i8` | | [2026-10-24-tanh-precision-symmetric-i8.md](2026-10-24-tanh-precision-symmetric-i8.md) | | 2026-10-23 | `ps2-compose-chain-preregistration` | | [2026-10-23-ps2-compose-chain-preregistration.md](2026-10-23-ps2-compose-chain-preregistration.md) | | 2026-10-23 | `ledger-revision-1635-1726` | | [2026-10-23-ledger-revision-1635-1726.md](2026-10-23-ledger-revision-1635-1726.md) | | 2026-10-23 | `EXP-GF2-ELIM-PREREG` | Gaussian elimination over GF(2) with a replayable row-operation witness — pre-registration | [2026-10-23-gf2-elimination-preregistration.md](2026-10-23-gf2-elimination-preregistration.md) | From dda286132131b0d8d56c08fb27d12ed1e2a71d48 Mon Sep 17 00:00:00 2001 From: Claude Date: Sat, 24 Oct 2026 06:53:57 +0000 Subject: [PATCH 4/4] lance-graph-onnx: replace the vsneg tripwire now that ndarray #368 fixed it vml_vsneg_is_not_onnx_neg_on_positive_zero pinned vsneg(+0) == +0 so it would fail once ndarray fixed the sign of zero. ndarray #368 did, and the test fails on main (Rust Tests run 40047762268). The replacement pins that vsneg and Neg both give -0.0 for +0.0 in the SIMD body and the tail. Neg stays inline: vsneg is f32-only and writes into a view. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01G22yT6htkcdyXsihxxXdrg --- .claude/plans/d-onnx-native-v1.md | 1 + crates/lance-graph-onnx/src/lib.rs | 6 ++--- crates/lance-graph-onnx/tests/conformance.rs | 23 ++++++++++---------- 3 files changed, 16 insertions(+), 14 deletions(-) diff --git a/.claude/plans/d-onnx-native-v1.md b/.claude/plans/d-onnx-native-v1.md index 18f582320..301b6f350 100644 --- a/.claude/plans/d-onnx-native-v1.md +++ b/.claude/plans/d-onnx-native-v1.md @@ -117,6 +117,7 @@ The old 50 ns target does not apply to aliases. An alias's cost is the native fu ## 6. Proposed (not implemented) - An ndarray owner fix for `vml::vsneg` (+0) and `vml::vsexp` (±inf, overflow halving): AdaWorldAPI/ndarray#368, open. After it lands, `Neg` could become an ALIAS. + - ⊘ 2026-10-24: #368 merged. The tripwire `vml_vsneg_is_not_onnx_neg_on_positive_zero` failed as designed and was replaced by `vml_vsneg_matches_onnx_neg_on_positive_zero` (lance-graph#1757). `Neg` stays INLINE for now: `vsneg` is f32-only and writes into a view, while `neg` is generic over f32/f64. - A native unmasked `i32 → i64` sum in `ndarray::simd`, so that `ReduceSum::sum_i32` becomes a direct call. - The NATIVE rows in §2, each with its open decision named. - `#1752` `onnx_compat` stays as the independent probe. This crate is the facade, and #1752's measurements are kept as the diagnostic baseline. diff --git a/crates/lance-graph-onnx/src/lib.rs b/crates/lance-graph-onnx/src/lib.rs index b9a313f87..bf99774e9 100644 --- a/crates/lance-graph-onnx/src/lib.rs +++ b/crates/lance-graph-onnx/src/lib.rs @@ -194,9 +194,9 @@ binary!( /// `Neg`. INLINE: `ndarray`'s element map with `-v`, the same body as its /// `Neg for &ArrayRef` impl. /// -/// Not an alias of `ndarray::hpc::vml::vsneg`: its 16-lane body computes -/// `0.0 - v`, so `Neg(+0.0)` gives `+0.0` instead of `-0.0` whenever the input -/// has at least 16 elements (its scalar tail computes `-v` and is correct). +/// Not an alias of `ndarray::hpc::vml::vsneg`, which is `f32`-only and writes +/// into an output view. Both give `-0.0` for `+0.0` (vsneg since ndarray +/// #368; pinned by the conformance test). #[inline] pub fn neg, D: Dimension>(a: &ArrayBase) -> Array { a.map(|&v| -v) diff --git a/crates/lance-graph-onnx/tests/conformance.rs b/crates/lance-graph-onnx/tests/conformance.rs index e4a004080..3fd9b6e57 100644 --- a/crates/lance-graph-onnx/tests/conformance.rs +++ b/crates/lance-graph-onnx/tests/conformance.rs @@ -71,19 +71,20 @@ fn unary_aliases_are_exact_on_special_values_and_every_tail_length() { } /// The aliases are the native items, not wrappers around them. -/// Why `Neg` is not an alias: the vml kernel's SIMD body loses the sign of -/// zero. If ndarray fixes it, this test fails and `Neg` can become an alias. +/// `vsneg` used to lose the sign of zero in its SIMD body, which is why `Neg` +/// is not an alias. ndarray #368 fixed it; this test now pins that `vsneg` +/// and `Neg` agree on +0 in the SIMD body and the tail. Making `Neg` an +/// alias is a separate change: `vsneg` is `f32`-only and writes into a view. #[test] -fn vml_vsneg_is_not_onnx_neg_on_positive_zero() { - let x = Array1::::zeros(16); - let mut out = Array1::::zeros(16); +fn vml_vsneg_matches_onnx_neg_on_positive_zero() { + let x = Array1::::zeros(17); + let mut out = Array1::::zeros(17); ndarray::hpc::vml::vsneg(x.view(), out.view_mut()); - assert_eq!( - out[0].to_bits(), - 0.0f32.to_bits(), - "vsneg(+0) is +0 on the SIMD path" - ); - assert_eq!(onnx::neg(&x)[0].to_bits(), (-0.0f32).to_bits()); + let y = onnx::neg(&x); + for i in [0, 16] { + assert_eq!(out[i].to_bits(), (-0.0f32).to_bits(), "vsneg(+0) at {i}"); + assert_eq!(y[i].to_bits(), (-0.0f32).to_bits(), "Neg(+0) at {i}"); + } } #[test]