perturbation-sim: NaN/panic hardening sweep (total_cmp + cascade NaN-abort + stats guard) before the kanban/jitson wiring

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);

crates/perturbation-sim/examples/calibrate.rs

@@ -88,7 +88,7 @@ fn quantize_rank_bits(col: &[f64], bits: u32) -> Vec<f64> {
     let n = col.len();
     let bins = (1usize << bits).min(n.max(1));
     let mut idx: Vec<usize> = (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<usize> = (0..k).map(|i| sens[(i * step).min(m - 1)].0).collect();

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<usize> = sens.iter().take(12).map(|x| x.0).collect();
     let lam2_0 = gl[1];
     let dl1: std::collections::HashMap<usize, f64> = top

crates/perturbation-sim/examples/hhtl_levels.rs

@@ -74,7 +74,7 @@ fn median(mut v: Vec<f64>) -> 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]
 }
 

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;

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();
 

crates/perturbation-sim/examples/meta_hops.rs

@@ -95,7 +95,7 @@ fn median(mut v: Vec<f64>) -> 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);

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<usize> = (0..k).map(|i| sens[(i * step).min(m - 1)].0).collect();

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,

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) ==",

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() {

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<f64>| {
-        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());

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();
 

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<usize> = (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<usize> = order.into_iter().take(m_seeds).collect();
 
     let factor_names = [

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

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}");

crates/perturbation-sim/src/cascade.rs

@@ -119,19 +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<f64>;
-    let mut flow: Vec<f64>;
+    // 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<f64> = theta_base.clone();
+    let mut flow: Vec<f64> = 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 (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
@@ -286,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<f64>, 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"
+            );
+        }
     }
 }

crates/perturbation-sim/src/stats.rs

@@ -105,10 +105,18 @@ pub fn zscore(x: &[f64]) -> Vec<f64> {
 }
 
 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::<f64>() / 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::<f64>() / x.len() as f64
 }