hdr: fix observe_batch alert gate and large-moment variance fallback

Two defects in the code recovered from #325, both found in review:

- observe_batch gated its drift alert on `old_n > 10`, while observe
  counts the new element first and tests `> 10`, i.e. >= 10 prior. At
  exactly 10 prior observations, observe_batch(&[x]) stayed silent where
  observe(x) alerted. Gate is now `old_n >= 10`.
- MomentsU32::variance's fallback (when n*sum_sq or sum^2 overflows
  u128) subtracted two ~2^64-scale f64 values and could cancel the
  variance to 0 (2^33 values split between u32::MAX and u32::MAX-1:
  0.0 instead of 0.25). It now centres on q = floor(sum/n) exactly in
  u128 first, so the float step works on variance-sized quantities.

Each fix has a regression test that failed before the change.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_019HnekoM1EidTwQLS3oFVFm

Author: claude

src/hpc/cascade.rs

@@ -220,8 +220,9 @@ impl Cascade {
     ///
     /// Drift is judged once per batch, not per element: an alert fires when
     /// the batch moves μ by more than 2σ of the pre-batch state (and the state
-    /// had > 10 observations and σ > 0), the same test `observe` applies to a
-    /// single step.
+    /// had at least 10 observations and σ > 0), the same test `observe`
+    /// applies to a single step, so `observe_batch(&[x])` alerts exactly when
+    /// `observe(x)` does.
     pub fn observe_batch(&mut self, distances: &[u32]) -> Option<ShiftAlert> {
         let b = crate::hpc::statistics::moments_u32(distances);
         if b.n == 0 {
@@ -245,7 +246,7 @@ impl Cascade {
         }
         self.observations = old_n + b.n as usize;
 
-        if old_n > 10 && old_sigma > 0.0 && (self.mu - old_mu).abs() > 2.0 * old_sigma {
+        if old_n >= 10 && old_sigma > 0.0 && (self.mu - old_mu).abs() > 2.0 * old_sigma {
             Some(ShiftAlert {
                 old_mu,
                 new_mu: self.mu,
@@ -888,6 +889,25 @@ mod tests {
         assert!(alert.new_mu > alert.old_mu + 2.0 * alert.old_sigma);
     }
 
+    /// A singleton batch is the same event as one `observe` call, so both
+    /// must agree on the alert gate at the boundary: exactly 10 prior
+    /// observations, then an outlier. `observe` counts the new element
+    /// before testing `> 10`; `observe_batch` must admit the same state.
+    #[test]
+    fn observe_batch_singleton_matches_observe_alert_gate() {
+        let mut seq = Cascade::from_threshold(8000, 2048);
+        let mut batch = Cascade::from_threshold(8000, 2048);
+        for i in 0..10u32 {
+            assert!(seq.observe(8000 + 100 * (i % 2)).is_none());
+            assert!(batch.observe(8000 + 100 * (i % 2)).is_none());
+        }
+        assert_eq!(seq.observations(), 10);
+        let a = seq.observe(20000);
+        let b = batch.observe_batch(&[20000]);
+        assert!(a.is_some(), "observe must alert on the outlier");
+        assert!(b.is_some(), "observe_batch(&[x]) must alert exactly like observe(x)");
+    }
+
     #[test]
     fn packed_database_roundtrip() {
         let vec_bytes = 256;

src/hpc/statistics.rs

@@ -410,7 +410,9 @@ impl MomentsU32 {
     /// Computed as `(n·Σx² − (Σx)²) / n²`. The numerator is formed exactly in
     /// `u128` whenever it fits (always, for Hamming-scale data: n ≤ 2³², x ≤
     /// 2¹⁷), so there is no cancellation between two large floats; only the
-    /// final division rounds. Past that range it falls back to `f64`.
+    /// final division rounds. Past that range it centres the sums on the
+    /// integer part of the mean in `u128` first, so the `f64` step still
+    /// works on small, variance-sized quantities.
     pub fn variance(&self) -> f64 {
         if self.n == 0 {
             return 0.0;
@@ -419,8 +421,18 @@ impl MomentsU32 {
         match (n.checked_mul(self.sum_sq), self.sum.checked_mul(self.sum)) {
             (Some(a), Some(b)) => (a - b) as f64 / (self.n as f64 * self.n as f64),
             _ => {
-                let mean = self.mean();
-                (self.sum_sq as f64 / self.n as f64 - mean * mean).max(0.0)
+                // Centre on the integer part of the mean, q = ⌊Σx / n⌋, with
+                // remainder r = Σx − n·q < n. Then, exactly in u128,
+                // Σ(x − q)² = Σx² − q·Σx − q·r, which is < n·2⁶⁴ and never
+                // negative at any step. The true M2 is that minus r²/n, and
+                // both terms are O(n·(σ² + 1)), so the one float subtraction
+                // cannot cancel the variance away.
+                let q = self.sum / n;
+                let r = self.sum % n;
+                let centred = self.sum_sq - q * self.sum - q * r;
+                let rf = r as f64;
+                let nf = self.n as f64;
+                ((centred as f64 - rf * rf / nf) / nf).max(0.0)
             }
         }
     }
@@ -481,6 +493,24 @@ pub fn moments_u32(values: &[u32]) -> MomentsU32 {
 mod moments_tests {
     use super::*;
 
+    /// Past the exact-`u128` range the variance must not cancel: 2³³ values
+    /// split evenly between `u32::MAX` and `u32::MAX - 1` have variance
+    /// exactly 0.25, and `n·Σx²` overflows `u128`, so this takes the
+    /// fallback path.
+    #[test]
+    fn variance_fallback_does_not_cancel() {
+        let half = 1u128 << 32;
+        let hi = u128::from(u32::MAX);
+        let lo = hi - 1;
+        let m = MomentsU32 {
+            n: 1u64 << 33,
+            sum: half * (hi + lo),
+            sum_sq: half * (hi * hi + lo * lo),
+        };
+        assert!(u128::from(m.n).checked_mul(m.sum_sq).is_none(), "fixture must take the fallback");
+        assert!((m.variance() - 0.25).abs() < 1e-9, "variance {}", m.variance());
+    }
+
     fn xorshift(n: usize, mut s: u64, mask: u32) -> Vec<u32> {
         (0..n)
             .map(|_| {