From 657f0fb3bd941af447c20da327eb41bb9cb931b8 Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 19 Jun 2026 17:20:30 +0000 Subject: [PATCH 1/2] =?UTF-8?q?perturbation-sim:=20NaN/panic=20hardening?= =?UTF-8?q?=20sweep=20=E2=80=94=20total=5Fcmp=20sorts,=20cascade=20NaN-abo?= =?UTF-8?q?rt,=20stats=20empty-guard?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Foundation-hardening before the kanban↔Lance↔jitson wiring rides on it. - NaN-panic fix: 20 `partial_cmp(...).unwrap()` sites across 14 example files replaced with `f64::total_cmp` (never panics; deterministic NaN order). The lib `src/` had none in non-test code — the foundation was not panicking; this hardens the example/probe surface (the "speaking errors"). - Cascade NaN-abort: `simulate_outage` now breaks if any `theta`/`flow` value is non-finite (matching the islanding-break style) — a NaN can neither propagate into the perturbation shape nor loop to max_rounds. The "timeout for NaN". - stats.rs: empty-slice guard on `mean`/`pop_var` (0/0 = NaN). Currently unreachable (all callers validate len ≥ 2) but wired at the source so a future caller can't leak it. - NaN-source audit (src/): eigen, flow, perturbation, resilience, rolling_floor float producers all already guarded (SPECTRAL_GAP_FLOOR / rel_tol / .max(eps) / filter) — verified, no change needed. - clippy --all-targets: zero warnings (default build, no ndarray-simd/cranelift). cargo-machete: not installed; only dep is the optional feature-gated ndarray. - 95 lib tests pass. Disk-light: no cranelift/ndarray build touched. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_01CcpLeEC3XK8Eye53GKBVvi --- crates/perturbation-sim/examples/buffer.rs | 2 +- crates/perturbation-sim/examples/calibrate.rs | 4 ++-- crates/perturbation-sim/examples/explore.rs | 2 +- crates/perturbation-sim/examples/hhtl_levels.rs | 2 +- crates/perturbation-sim/examples/hhtl_resident.rs | 2 +- crates/perturbation-sim/examples/iberian.rs | 2 +- crates/perturbation-sim/examples/meta_hops.rs | 4 ++-- crates/perturbation-sim/examples/modifier_validity.rs | 2 +- crates/perturbation-sim/examples/reinforce.rs | 4 ++-- crates/perturbation-sim/examples/resilience.rs | 2 +- crates/perturbation-sim/examples/rolling_floor.rs | 2 +- crates/perturbation-sim/examples/scorecard.rs | 4 ++-- crates/perturbation-sim/examples/simulate.rs | 2 +- crates/perturbation-sim/examples/validate.rs | 2 +- crates/perturbation-sim/examples/validate_mediators.rs | 2 +- crates/perturbation-sim/examples/weakest_links.rs | 2 +- crates/perturbation-sim/src/cascade.rs | 7 +++++++ crates/perturbation-sim/src/stats.rs | 8 ++++++++ 18 files changed, 35 insertions(+), 20 deletions(-) diff --git a/crates/perturbation-sim/examples/buffer.rs b/crates/perturbation-sim/examples/buffer.rs index ebb3ee5e..66183358 100644 --- a/crates/perturbation-sim/examples/buffer.rs +++ b/crates/perturbation-sim/examples/buffer.rs @@ -214,7 +214,7 @@ fn main() { let mut weakest: Vec<(usize, f64, f64)> = (0..comps.len()) .map(|ci| (ci, mean_buf[ci], lam2s[ci])) .collect(); - weakest.sort_by(|a, b| a.1.partial_cmp(&b.1).unwrap()); + weakest.sort_by(|a, b| a.1.total_cmp(&b.1)); let mut first_yield = None; for (ci, buf, l2) in &weakest { let y = ketchup_yield(impulse, *buf); diff --git a/crates/perturbation-sim/examples/calibrate.rs b/crates/perturbation-sim/examples/calibrate.rs index a22d02e5..ae80796b 100644 --- a/crates/perturbation-sim/examples/calibrate.rs +++ b/crates/perturbation-sim/examples/calibrate.rs @@ -88,7 +88,7 @@ fn quantize_rank_bits(col: &[f64], bits: u32) -> Vec { let n = col.len(); let bins = (1usize << bits).min(n.max(1)); let mut idx: Vec = (0..n).collect(); - idx.sort_by(|&a, &b| col[a].partial_cmp(&col[b]).unwrap()); + idx.sort_by(|&a, &b| col[a].total_cmp(&col[b])); let mut out = vec![0.0; n]; for b in 0..bins { let s = b * n / bins; @@ -142,7 +142,7 @@ fn main() { (e, d * d * grid.edges[e].susceptance) }) .collect(); - sens.sort_by(|a, b| b.1.partial_cmp(&a.1).unwrap()); + sens.sort_by(|a, b| b.1.total_cmp(&a.1)); let k = 24.min(m); let step = (m / k).max(1); let cand: Vec = (0..k).map(|i| sens[(i * step).min(m - 1)].0).collect(); diff --git a/crates/perturbation-sim/examples/explore.rs b/crates/perturbation-sim/examples/explore.rs index 76392bf6..0356c941 100644 --- a/crates/perturbation-sim/examples/explore.rs +++ b/crates/perturbation-sim/examples/explore.rs @@ -338,7 +338,7 @@ fn main() { (e, d * d * grid.edges[e].susceptance) }) .collect(); - sens.sort_by(|a, b| b.1.partial_cmp(&a.1).unwrap()); + sens.sort_by(|a, b| b.1.total_cmp(&a.1)); let top: Vec = sens.iter().take(12).map(|x| x.0).collect(); let lam2_0 = gl[1]; let dl1: std::collections::HashMap = top diff --git a/crates/perturbation-sim/examples/hhtl_levels.rs b/crates/perturbation-sim/examples/hhtl_levels.rs index 4b07b4bc..cf92c7dd 100644 --- a/crates/perturbation-sim/examples/hhtl_levels.rs +++ b/crates/perturbation-sim/examples/hhtl_levels.rs @@ -74,7 +74,7 @@ fn median(mut v: Vec) -> f64 { if v.is_empty() { return 0.0; } - v.sort_by(|a, b| a.partial_cmp(b).unwrap()); + v.sort_by(|a, b| a.total_cmp(b)); v[v.len() / 2] } diff --git a/crates/perturbation-sim/examples/hhtl_resident.rs b/crates/perturbation-sim/examples/hhtl_resident.rs index cedcc520..ecb013ca 100644 --- a/crates/perturbation-sim/examples/hhtl_resident.rs +++ b/crates/perturbation-sim/examples/hhtl_resident.rs @@ -150,7 +150,7 @@ fn main() { let seed = base .iter() .enumerate() - .max_by(|x, y| x.1.abs().partial_cmp(&y.1.abs()).unwrap()) + .max_by(|x, y| x.1.abs().total_cmp(&y.1.abs())) .map(|(i, _)| i) .unwrap_or(0); let infight = simulate_outage(&sub, &p, seed, cfg).fraction_tripped; diff --git a/crates/perturbation-sim/examples/iberian.rs b/crates/perturbation-sim/examples/iberian.rs index 727dd3eb..6ba8749f 100644 --- a/crates/perturbation-sim/examples/iberian.rs +++ b/crates/perturbation-sim/examples/iberian.rs @@ -115,7 +115,7 @@ fn main() { let seed = base .iter() .enumerate() - .max_by(|a, b| a.1.abs().partial_cmp(&b.1.abs()).unwrap()) + .max_by(|a, b| a.1.abs().total_cmp(&b.1.abs())) .map(|(i, _)| i) .unwrap(); diff --git a/crates/perturbation-sim/examples/meta_hops.rs b/crates/perturbation-sim/examples/meta_hops.rs index 1fceaabc..4a5c69b3 100644 --- a/crates/perturbation-sim/examples/meta_hops.rs +++ b/crates/perturbation-sim/examples/meta_hops.rs @@ -95,7 +95,7 @@ fn median(mut v: Vec) -> f64 { if v.is_empty() { return 0.0; } - v.sort_by(|a, b| a.partial_cmp(b).unwrap()); + v.sort_by(|a, b| a.total_cmp(b)); v[v.len() / 2] } @@ -169,7 +169,7 @@ fn main() { let seed = base .iter() .enumerate() - .max_by(|x, y| x.1.abs().partial_cmp(&y.1.abs()).unwrap()) + .max_by(|x, y| x.1.abs().total_cmp(&y.1.abs())) .map(|(i, _)| i) .unwrap_or(0); infs.push(simulate_outage(&sub, &p, seed, cfg).fraction_tripped); diff --git a/crates/perturbation-sim/examples/modifier_validity.rs b/crates/perturbation-sim/examples/modifier_validity.rs index ab606a0b..99e2c323 100644 --- a/crates/perturbation-sim/examples/modifier_validity.rs +++ b/crates/perturbation-sim/examples/modifier_validity.rs @@ -131,7 +131,7 @@ fn main() { (e, d * d * grid.edges[e].susceptance) }) .collect(); - sens.sort_by(|x, y| y.1.partial_cmp(&x.1).unwrap()); + sens.sort_by(|x, y| y.1.total_cmp(&x.1)); let k = 30.min(m); let step = (m / k).max(1); let cand: Vec = (0..k).map(|i| sens[(i * step).min(m - 1)].0).collect(); diff --git a/crates/perturbation-sim/examples/reinforce.rs b/crates/perturbation-sim/examples/reinforce.rs index b2968a2c..cbc0987f 100644 --- a/crates/perturbation-sim/examples/reinforce.rs +++ b/crates/perturbation-sim/examples/reinforce.rs @@ -95,7 +95,7 @@ fn main() { } } } - cands.sort_by(|x, y| y.2.partial_cmp(&x.2).unwrap()); + cands.sort_by(|x, y| y.2.total_cmp(&x.2)); println!("\n== Reinforcement: optimal 3rd corridor across the Cheeger seam =="); println!(" base λ₂ = {lam2:.3e} (new-line susceptance b = {w_new:.3})"); @@ -141,7 +141,7 @@ fn main() { .collect(); let seed = *seam .iter() - .max_by(|&&x, &&y| base[x].abs().partial_cmp(&base[y].abs()).unwrap()) + .max_by(|&&x, &&y| base[x].abs().total_cmp(&base[y].abs())) .unwrap_or(&0); let cfg = CascadeConfig { max_rounds: 16, diff --git a/crates/perturbation-sim/examples/resilience.rs b/crates/perturbation-sim/examples/resilience.rs index d614e156..822c9556 100644 --- a/crates/perturbation-sim/examples/resilience.rs +++ b/crates/perturbation-sim/examples/resilience.rs @@ -131,7 +131,7 @@ fn main() { .map(|(i, b)| (i, cert(&induced(&grid, b)))) .collect(); // Weakest compartment first = smallest λ₂ margin. - certs.sort_by(|a, b| a.1.lambda2.partial_cmp(&b.1.lambda2).unwrap()); + certs.sort_by(|a, b| a.1.lambda2.total_cmp(&b.1.lambda2)); println!( "\n== Per-compartment certificate ({} compartments, weakest λ₂ first) ==", diff --git a/crates/perturbation-sim/examples/rolling_floor.rs b/crates/perturbation-sim/examples/rolling_floor.rs index a6c6a323..e477d1b2 100644 --- a/crates/perturbation-sim/examples/rolling_floor.rs +++ b/crates/perturbation-sim/examples/rolling_floor.rs @@ -142,7 +142,7 @@ fn main() { if v.is_empty() { return 0.0; } - v.sort_by(|a, b| a.partial_cmp(b).unwrap()); + v.sort_by(|a, b| a.total_cmp(b)); v[v.len() / 2] }; for (l, level) in levels.iter().enumerate() { diff --git a/crates/perturbation-sim/examples/scorecard.rs b/crates/perturbation-sim/examples/scorecard.rs index 8991aa87..8478448c 100644 --- a/crates/perturbation-sim/examples/scorecard.rs +++ b/crates/perturbation-sim/examples/scorecard.rs @@ -189,7 +189,7 @@ fn main() { } // Rank by exposure (most exposed first). - rows.sort_by(|a, b| b.6.partial_cmp(&a.6).unwrap()); + rows.sort_by(|a, b| b.6.total_cmp(&a.6)); println!("\n== Exposure ranking (most exposed first) =="); for (i, r) in rows.iter().enumerate() { println!(" {}. {:<10} exposure = {:.3e}", i + 1, r.0, r.6); @@ -200,7 +200,7 @@ fn main() { // the panel median: topology (mean R), buffer (1/H storage), or policy. The // binding axis dictates the intervention TYPE and the marginal exposure cut. let med = |mut v: Vec| { - v.sort_by(|a, b| a.partial_cmp(b).unwrap()); + v.sort_by(|a, b| a.total_cmp(b)); v[v.len() / 2] }; let med_r = med(rows.iter().map(|r| r.2).collect()); diff --git a/crates/perturbation-sim/examples/simulate.rs b/crates/perturbation-sim/examples/simulate.rs index b49d80a0..c9260eaa 100644 --- a/crates/perturbation-sim/examples/simulate.rs +++ b/crates/perturbation-sim/examples/simulate.rs @@ -64,7 +64,7 @@ fn main() { let seed = base .iter() .enumerate() - .max_by(|a, b| a.1.abs().partial_cmp(&b.1.abs()).unwrap()) + .max_by(|a, b| a.1.abs().total_cmp(&b.1.abs())) .map(|(i, _)| i) .unwrap(); diff --git a/crates/perturbation-sim/examples/validate.rs b/crates/perturbation-sim/examples/validate.rs index 7bd2bb3a..1cf5760a 100644 --- a/crates/perturbation-sim/examples/validate.rs +++ b/crates/perturbation-sim/examples/validate.rs @@ -95,7 +95,7 @@ fn main() { dc_flows(&grid, &alive, &eig.pseudo_apply(&p0, 1e-9)) }; let mut order: Vec = (0..base0.len()).collect(); - order.sort_by(|&a, &b| base0[b].abs().partial_cmp(&base0[a].abs()).unwrap()); + order.sort_by(|&a, &b| base0[b].abs().total_cmp(&base0[a].abs())); let seeds: Vec = order.into_iter().take(m_seeds).collect(); let factor_names = [ diff --git a/crates/perturbation-sim/examples/validate_mediators.rs b/crates/perturbation-sim/examples/validate_mediators.rs index 1cdd8754..346a2e48 100644 --- a/crates/perturbation-sim/examples/validate_mediators.rs +++ b/crates/perturbation-sim/examples/validate_mediators.rs @@ -116,7 +116,7 @@ fn main() { (e, d * d * grid.edges[e].susceptance) }) .collect(); - sens.sort_by(|x, y| y.1.partial_cmp(&x.1).unwrap()); + sens.sort_by(|x, y| y.1.total_cmp(&x.1)); // Stride-sample K lines ACROSS the full sensitivity ranking (not the top-K): // the top-K are nearly all bridges (raumgewinn ≡ 1, no variety), so a stride diff --git a/crates/perturbation-sim/examples/weakest_links.rs b/crates/perturbation-sim/examples/weakest_links.rs index 6cad96b5..5106d475 100644 --- a/crates/perturbation-sim/examples/weakest_links.rs +++ b/crates/perturbation-sim/examples/weakest_links.rs @@ -87,7 +87,7 @@ fn main() { (e, d * d * grid.edges[e].susceptance) // first-order Δλ₂ proxy }) .collect(); - struct_rank.sort_by(|x, y| y.1.partial_cmp(&x.1).unwrap()); + struct_rank.sort_by(|x, y| y.1.total_cmp(&x.1)); println!("== 1. Structural weakest links (first-order Fiedler sensitivity ∂λ₂/∂wₑ) =="); println!(" base λ₂ = {lam2:.3e}"); diff --git a/crates/perturbation-sim/src/cascade.rs b/crates/perturbation-sim/src/cascade.rs index 070a1654..79ec7ef0 100644 --- a/crates/perturbation-sim/src/cascade.rs +++ b/crates/perturbation-sim/src/cascade.rs @@ -133,6 +133,13 @@ pub fn simulate_outage( theta = eig.pseudo_apply(p, cfg.rel_tol); flow = dc_flows(grid, &alive, &theta); + // NaN guard: a non-finite angle/flow signals numerical failure; treat + // as terminal (matching the islanding-break style) so NaN cannot + // propagate into the output shape or keep the loop alive. + if theta.iter().any(|x| !x.is_finite()) || flow.iter().any(|x| !x.is_finite()) { + break; + } + if components_final > 1 { // Network fragmented: injections no longer balance per island, so // the DC solution is only the least-norm proxy. Treat as terminal diff --git a/crates/perturbation-sim/src/stats.rs b/crates/perturbation-sim/src/stats.rs index 5a1e60ec..b9d84ba5 100644 --- a/crates/perturbation-sim/src/stats.rs +++ b/crates/perturbation-sim/src/stats.rs @@ -105,10 +105,18 @@ pub fn zscore(x: &[f64]) -> Vec { } fn mean(x: &[f64]) -> f64 { + // Empty-slice guard: 0/0 = NaN. All current callers validate len ≥ 2, but + // wire the NaN at the source so a future caller can never leak it. + if x.is_empty() { + return 0.0; + } x.iter().sum::() / x.len() as f64 } fn pop_var(x: &[f64]) -> f64 { + if x.is_empty() { + return 0.0; + } let m = mean(x); x.iter().map(|&v| (v - m).powi(2)).sum::() / x.len() as f64 } From ebd24683d38b1e494ceb3569cea1104335c18b7b Mon Sep 17 00:00:00 2001 From: Claude Date: Fri, 19 Jun 2026 17:31:18 +0000 Subject: [PATCH 2/2] cascade: NaN-abort preserves last finite state, never returns NaN shape (codex #554 P2) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The previous guard broke the loop on a non-finite theta/flow, but the shape was then built from those SAME non-finite values (node_field/edge_field at the build step) — so the "abort" still leaked NaN/Inf into the perturbation shape and downstream rankings/stats. Fix: seed theta/flow with the finite base state and overwrite them each round ONLY after the new values are confirmed finite (theta_next/flow_next). A non-finite round now breaks WITHOUT adopting it, so theta/flow/components_final retain the last finite state and the shape is built from that — never from NaN/Inf. Matches the islanding-break semantics (which already left a finite proxy). New test perturbation_shape_is_always_finite locks the invariant across all terminal paths (converge / cascade / island); islanding_is_flagged also now asserts finite shape. clippy --all-targets clean; cascade suite green. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_01CcpLeEC3XK8Eye53GKBVvi --- crates/perturbation-sim/src/cascade.rs | 99 +++++++++++++++++++++++--- 1 file changed, 89 insertions(+), 10 deletions(-) diff --git a/crates/perturbation-sim/src/cascade.rs b/crates/perturbation-sim/src/cascade.rs index 79ec7ef0..8c134aea 100644 --- a/crates/perturbation-sim/src/cascade.rs +++ b/crates/perturbation-sim/src/cascade.rs @@ -119,26 +119,36 @@ pub fn simulate_outage( // Assigned on every loop iteration before any break (the loop body always // runs at least once), so no initializer is needed. - let mut theta: Vec; - let mut flow: Vec; + // Seed `theta`/`flow` with the finite base state and overwrite them each + // round ONLY after the new values are confirmed finite — so they always + // hold the last finite state. If a round goes non-finite we break without + // adopting it, and the shape below is built from that last-good state, + // never from NaN/Inf (codex #554 P2). + let mut theta: Vec = theta_base.clone(); + let mut flow: Vec = flow_base.clone(); let mut islanded = false; - let mut components_final: usize; + let mut components_final: usize = 1; let mut rounds = 0usize; loop { rounds += 1; let eig = symmetric_eigen(&grid.laplacian_of(&alive), n); - components_final = eig.nullity(cfg.rel_tol); + let nullity = eig.nullity(cfg.rel_tol); - theta = eig.pseudo_apply(p, cfg.rel_tol); - flow = dc_flows(grid, &alive, &theta); + let theta_next = eig.pseudo_apply(p, cfg.rel_tol); + let flow_next = dc_flows(grid, &alive, &theta_next); - // NaN guard: a non-finite angle/flow signals numerical failure; treat - // as terminal (matching the islanding-break style) so NaN cannot - // propagate into the output shape or keep the loop alive. - if theta.iter().any(|x| !x.is_finite()) || flow.iter().any(|x| !x.is_finite()) { + // NaN guard: a non-finite angle/flow signals numerical failure (solver + // breakdown, overflow). Break WITHOUT adopting the bad round — `theta`/ + // `flow`/`components_final` retain the last finite state, so the + // perturbation shape can never carry NaN/Inf into node_field/edge_field + // or downstream rankings/statistics (codex #554 P2). + if theta_next.iter().any(|x| !x.is_finite()) || flow_next.iter().any(|x| !x.is_finite()) { break; } + theta = theta_next; + flow = flow_next; + components_final = nullity; if components_final > 1 { // Network fragmented: injections no longer balance per island, so @@ -293,5 +303,74 @@ mod tests { let r = simulate_outage(&g, &p, 6, CascadeConfig::default()); assert!(r.islanded, "bridge trip must island the network"); assert_eq!(r.components_final, 2); + // The islanding shape (least-norm proxy) must still be finite. + assert!(r.shape.node_field.iter().all(|x| x.is_finite())); + assert!(r.shape.edge_field.iter().all(|x| x.is_finite())); + } + + #[test] + fn perturbation_shape_is_always_finite() { + // codex #554 P2: the perturbation shape must NEVER carry NaN/Inf — the + // NaN guard preserves the last finite theta/flow, so node_field / + // edge_field stay finite across every terminal path (converge, cascade, + // island, max_rounds). Locks the invariant the fix guarantees. + let cases: [(Grid, Vec, usize); 3] = [ + // generous limits: only the seed trips (converge path). + ( + Grid::new( + 4, + vec![ + Edge::new(0, 1, 1.0, 1e6), + Edge::new(1, 2, 1.0, 1e6), + Edge::new(2, 3, 1.0, 1e6), + Edge::new(3, 0, 1.0, 1e6), + ], + ), + vec![1.0, 0.0, -1.0, 0.0], + 0, + ), + // tight limits: multi-line cascade path. + ( + Grid::new( + 4, + vec![ + Edge::new(0, 1, 1.0, 0.6), + Edge::new(1, 2, 1.0, 0.6), + Edge::new(2, 3, 1.0, 0.6), + Edge::new(3, 0, 1.0, 0.6), + ], + ), + vec![1.0, 0.0, -1.0, 0.0], + 0, + ), + // two triangles + bridge: islanding path. + ( + Grid::new( + 6, + vec![ + Edge::new(0, 1, 1.0, 1e6), + Edge::new(1, 2, 1.0, 1e6), + Edge::new(2, 0, 1.0, 1e6), + Edge::new(3, 4, 1.0, 1e6), + Edge::new(4, 5, 1.0, 1e6), + Edge::new(5, 3, 1.0, 1e6), + Edge::new(2, 3, 1.0, 1e6), + ], + ), + vec![1.0, 0.0, 0.0, 0.0, 0.0, -1.0], + 6, + ), + ]; + for (g, p, seed) in cases { + let r = simulate_outage(&g, &p, seed, CascadeConfig::default()); + assert!( + r.shape.node_field.iter().all(|x| x.is_finite()), + "node_field carried a non-finite value" + ); + assert!( + r.shape.edge_field.iter().all(|x| x.is_finite()), + "edge_field carried a non-finite value" + ); + } } }