Skip to main content

powerio_matrix/
pipeline.rs

1//! Orchestrates a single case → output directory.
2//!
3//! Given a parsed `BalancedNetwork`, builds the requested matrix family, writes
4//! `.mtx` files, and emits a `meta.json` sidecar describing what was
5//! produced. Used by both the `batch` CLI subcommand and the TUI's
6//! batch export screen.
7
8use std::path::{Path, PathBuf};
9
10use rand::SeedableRng;
11use rand::distr::{Distribution, StandardUniform};
12use rand_chacha::ChaCha8Rng;
13use serde::{Deserialize, Serialize};
14
15use crate::Result;
16use crate::indexed::IndexedNetwork;
17use crate::io::meta::{CaseMetadata, MatrixMetadata, write_meta_json};
18use crate::io::mtx::{write_mtx, write_vector_mtx};
19use crate::matrix::{
20    BuildOptions, MatrixStats, ZeroImpedanceRule, ZeroImpedanceSkips, build_adjacency,
21    build_bdoubleprime, build_bprime, build_lacpf, build_ybus, negate_into, sddm_check,
22    skipped_zero_impedance,
23};
24use crate::network::BalancedNetwork;
25
26#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
27#[non_exhaustive]
28pub enum MatrixKind {
29    /// MATPOWER FDPF Bp matrix.
30    BPrime,
31    /// MATPOWER FDPF Bpp matrix.
32    BDoublePrime,
33    /// `Re(Y_bus)` — full conductance matrix.
34    YbusG,
35    /// `-Im(Y_bus)` — full bus susceptance matrix (positive convention).
36    YbusB,
37    /// LACPF block: `[[G, -B], [-B, -G]]`, 2n × 2n indefinite.
38    Lacpf,
39    /// 0/1 bus adjacency matrix.
40    Adjacency,
41}
42
43impl MatrixKind {
44    pub const ALL: &'static [MatrixKind] = &[
45        Self::BPrime,
46        Self::BDoublePrime,
47        Self::YbusG,
48        Self::YbusB,
49        Self::Lacpf,
50        Self::Adjacency,
51    ];
52
53    pub fn slug(self) -> &'static str {
54        match self {
55            Self::BPrime => "bprime",
56            Self::BDoublePrime => "bdoubleprime",
57            Self::YbusG => "ybus_real",
58            Self::YbusB => "ybus_imag",
59            Self::Lacpf => "lacpf",
60            Self::Adjacency => "adjacency",
61        }
62    }
63
64    pub fn label(self) -> &'static str {
65        match self {
66            Self::BPrime => "MATPOWER Bp (FDPF)",
67            Self::BDoublePrime => "MATPOWER Bpp (FDPF)",
68            Self::YbusG => "Re(Y_bus)",
69            Self::YbusB => "-Im(Y_bus)",
70            Self::Lacpf => "LACPF block (2n×2n)",
71            Self::Adjacency => "adjacency (0/1)",
72        }
73    }
74}
75
76/// How to populate the RHS vector(s) emitted alongside each matrix.
77#[derive(Debug, Clone, Copy, PartialEq, Eq, Default, Serialize, Deserialize)]
78pub enum RhsKind {
79    #[default]
80    None,
81    /// Zero-mean Gaussian random (deterministic from `rng_seed`).
82    Random,
83    /// Power injections from the case: `b = (Pd, Qd) / baseMVA`.
84    Injection,
85}
86
87#[derive(Debug, Clone)]
88pub struct Pipeline {
89    pub matrices: Vec<MatrixKind>,
90    pub options: BuildOptions,
91    pub rhs: RhsKind,
92    pub rng_seed: u64,
93    pub source_file: Option<PathBuf>,
94}
95
96impl Default for Pipeline {
97    fn default() -> Self {
98        Self {
99            matrices: vec![MatrixKind::BPrime],
100            options: BuildOptions::default(),
101            rhs: RhsKind::None,
102            rng_seed: 0x00C0_FFEE,
103            source_file: None,
104        }
105    }
106}
107
108#[derive(Debug, Clone)]
109pub struct PipelineOutputs {
110    pub case_name: String,
111    pub files: Vec<PathBuf>,
112    pub metadata: CaseMetadata,
113}
114
115/// Longest sanitized stem before truncation. Output filenames append fixed
116/// suffixes (`_sensitivity_meta.json` is the longest at 22 bytes) and the
117/// collision hash adds 9 more, so the longest filename stays well under the
118/// common 255 byte component limit.
119const MAX_STEM_LEN: usize = 120;
120
121/// Hex digits of the SHA-256 disambiguator appended to a sanitized stem.
122/// 64 bits, so a batch export cannot be steered into an overwrite by
123/// searching for a second name that hashes into an existing stem.
124const DIGEST_LEN: usize = 16;
125
126/// Filenames Windows reserves for devices, matched against the part of a
127/// name before its first dot, case insensitively: `con` and `con.4` both
128/// resolve to the CON device.
129const WINDOWS_RESERVED: &[&str] = &[
130    "con", "prn", "aux", "nul", "com1", "com2", "com3", "com4", "com5", "com6", "com7", "com8",
131    "com9", "lpt1", "lpt2", "lpt3", "lpt4", "lpt5", "lpt6", "lpt7", "lpt8", "lpt9",
132];
133
134/// Reduces a network name to a single safe filename stem. The name comes from
135/// input content, so an unsanitized value like `../../etc/x` or `/abs/x` used
136/// in `out_dir.join(...)` would write outside `out_dir`. Keeps ASCII
137/// alphanumerics, `-`, `_`, and `.`; every other character (path separators
138/// included) becomes `_`. Leading dots are stripped so the result is never
139/// `.`, `..`, or a hidden file; trailing dots are trimmed and Windows
140/// reserved device names get a leading `_`, both invalid as Windows
141/// filenames; the length is capped at 120 bytes; and an empty result
142/// falls back to `case`.
143///
144/// A name the sanitizer had to change carries a hash of the original, so two
145/// distinct names that would sanitize identically (`a/b` and `a_b`, say)
146/// cannot silently overwrite each other's files in a multi case export. A
147/// name that already ends in the suffix shape is hashed too, so a case can
148/// never be named to impersonate another case's disambiguated stem: the
149/// suffixed and unsuffixed name spaces stay disjoint. Any other safe name
150/// passes through unchanged.
151pub fn sanitize_stem(name: &str) -> String {
152    // Skipping leading dots before mapping is equivalent to trimming them
153    // after: the map is the identity on `.` and never produces one, so the set
154    // of leading dots is the same either way.
155    let mut stem: String = name
156        .chars()
157        .skip_while(|&c| c == '.')
158        .map(|c| {
159            if c.is_ascii_alphanumeric() || matches!(c, '-' | '_' | '.') {
160                c
161            } else {
162                '_'
163            }
164        })
165        .collect();
166    while stem.ends_with('.') {
167        stem.pop();
168    }
169    let pre_dot = stem.split('.').next().unwrap_or("");
170    if WINDOWS_RESERVED
171        .iter()
172        .any(|r| pre_dot.eq_ignore_ascii_case(r))
173    {
174        stem.insert(0, '_');
175    }
176    if stem.is_empty() {
177        stem.push_str("case");
178    }
179    if stem == name && stem.len() <= MAX_STEM_LEN && !ends_with_digest(&stem) {
180        return stem;
181    }
182    // The map above only emits ASCII, so byte truncation cannot split a
183    // character, and the budget leaves room for the suffix.
184    stem.truncate(MAX_STEM_LEN - DIGEST_LEN - 1);
185    stem.push('-');
186    stem.push_str(&sha256_hex(name.as_bytes())[..DIGEST_LEN]);
187    stem
188}
189
190/// Whether a stem already carries the disambiguating suffix shape. Such a
191/// name gets a suffix of its own, so no name can be chosen to collide with
192/// another name's disambiguated stem.
193fn ends_with_digest(stem: &str) -> bool {
194    stem.len() > DIGEST_LEN
195        && stem.as_bytes()[stem.len() - DIGEST_LEN - 1] == b'-'
196        && stem[stem.len() - DIGEST_LEN..]
197            .bytes()
198            .all(|b| b.is_ascii_digit() || (b'a'..=b'f').contains(&b))
199}
200
201impl Pipeline {
202    pub fn run(&self, net: &BalancedNetwork, out_dir: impl AsRef<Path>) -> Result<PipelineOutputs> {
203        let out_dir = out_dir.as_ref();
204        std::fs::create_dir_all(out_dir)?;
205
206        let view = IndexedNetwork::new(net);
207        // The network name comes from input content, so it must not steer the
208        // output path: a name like `../../x` or `/abs/x` would otherwise write
209        // outside `out_dir`. Filenames use the sanitized stem; the metadata
210        // keeps the original name.
211        let stem = sanitize_stem(view.name());
212
213        let mut files = Vec::new();
214        let mut matrices_meta = Vec::new();
215        let mut ybus_cache = None;
216
217        for &kind in &self.matrices {
218            let matrix_path = out_dir.join(format!("{stem}_{}.mtx", kind.slug()));
219            let matrix = self.build_for_run(&view, kind, &mut ybus_cache)?;
220            write_mtx(&matrix, &matrix_path)?;
221            let stats = matrix_stats_for_kind(&matrix, &view, kind, &self.options);
222            let sddm = sddm_check(&matrix);
223            matrices_meta.push(MatrixMetadata {
224                kind: kind.slug().to_string(),
225                file: matrix_path
226                    .file_name()
227                    .and_then(|s| s.to_str())
228                    .unwrap_or("")
229                    .to_string(),
230                stats,
231                sddm,
232            });
233            files.push(matrix_path);
234
235            // RHS for matrices that take a RHS of length n (skip LACPF which is 2n).
236            if let Some(rhs) = self.build_rhs(&view, kind) {
237                let rhs_path = out_dir.join(format!("{stem}_{}_rhs.mtx", kind.slug()));
238                write_vector_mtx(&rhs, &rhs_path)?;
239                files.push(rhs_path);
240            }
241        }
242
243        // Shunt vector as a sidecar (not always meaningful, but cheap).
244        let shunt_path = out_dir.join(format!("{stem}_shunt.mtx"));
245        let base = view.per_unit_base();
246        let shunt: Vec<f64> = view.bs().iter().map(|&b| b / base).collect();
247        write_vector_mtx(&shunt, &shunt_path)?;
248        files.push(shunt_path);
249
250        let metadata = CaseMetadata {
251            case_name: view.name().to_string(),
252            source_file: self
253                .source_file
254                .as_ref()
255                .and_then(|p| p.to_str())
256                .map(str::to_string),
257            source_sha256: self
258                .source_file
259                .as_ref()
260                .and_then(|p| std::fs::read(p).ok())
261                .map(|b| sha256_hex(&b)),
262            base_mva: view.base_mva(),
263            n_buses: view.n(),
264            n_branches: view.branches().len(),
265            build_options: self.options.clone(),
266            matrices: matrices_meta,
267            powerio_version: env!("CARGO_PKG_VERSION").to_string(),
268        };
269        let meta_path = out_dir.join(format!("{stem}_meta.json"));
270        write_meta_json(&metadata, &meta_path)?;
271        files.push(meta_path);
272
273        Ok(PipelineOutputs {
274            case_name: view.name().to_string(),
275            files,
276            metadata,
277        })
278    }
279
280    fn build_for_run(
281        &self,
282        case: &IndexedNetwork,
283        kind: MatrixKind,
284        ybus_cache: &mut Option<YbusCache>,
285    ) -> Result<sprs::CsMat<f64>> {
286        match kind {
287            MatrixKind::YbusG => take_ybus_g(case, &self.options, ybus_cache),
288            MatrixKind::YbusB => take_ybus_b(case, &self.options, ybus_cache),
289            _ => build_kind(case, kind, &self.options),
290        }
291    }
292
293    fn build_rhs(&self, case: &IndexedNetwork, kind: MatrixKind) -> Option<Vec<f64>> {
294        // No meaningful RHS for the 2n LACPF block or the structural adjacency.
295        if matches!(self.rhs, RhsKind::None)
296            || matches!(kind, MatrixKind::Lacpf | MatrixKind::Adjacency)
297        {
298            return None;
299        }
300        let n = case.n();
301        let v = match self.rhs {
302            RhsKind::Random => {
303                let mut rng = ChaCha8Rng::seed_from_u64(self.rng_seed.wrapping_add(kind as u64));
304                let dist = StandardUniform;
305                let mut v: Vec<f64> = (0..n)
306                    .map(|_| {
307                        let u: f64 = dist.sample(&mut rng);
308                        u - 0.5
309                    })
310                    .collect();
311                let mean = v.iter().sum::<f64>() / n as f64;
312                for x in &mut v {
313                    *x -= mean; // zero-mean for Laplacian compatibility
314                }
315                v
316            }
317            RhsKind::Injection => {
318                let base = case.per_unit_base();
319                match kind {
320                    MatrixKind::BPrime | MatrixKind::YbusG | MatrixKind::YbusB => {
321                        case.pd().iter().map(|&p| -p / base).collect()
322                    }
323                    MatrixKind::BDoublePrime => case.qd().iter().map(|&q| -q / base).collect(),
324                    MatrixKind::Lacpf | MatrixKind::Adjacency => unreachable!(),
325                }
326            }
327            RhsKind::None => unreachable!(),
328        };
329        Some(v)
330    }
331}
332
333struct YbusCache {
334    g: Option<sprs::CsMat<f64>>,
335    b: Option<sprs::CsMat<f64>>,
336}
337
338fn fill_ybus_cache(
339    view: &IndexedNetwork,
340    opts: &BuildOptions,
341    ybus_cache: &mut Option<YbusCache>,
342) -> Result<()> {
343    let parts = build_ybus(view, opts)?;
344    *ybus_cache = Some(YbusCache {
345        g: Some(parts.g),
346        b: Some(parts.b),
347    });
348    Ok(())
349}
350
351fn take_ybus_g(
352    view: &IndexedNetwork,
353    opts: &BuildOptions,
354    ybus_cache: &mut Option<YbusCache>,
355) -> Result<sprs::CsMat<f64>> {
356    if ybus_cache.as_ref().is_none_or(|c| c.g.is_none()) {
357        fill_ybus_cache(view, opts, ybus_cache)?;
358    }
359    Ok(ybus_cache
360        .as_mut()
361        .and_then(|c| c.g.take())
362        .expect("Ybus cache was just filled, so the real part must be present"))
363}
364
365fn take_ybus_b(
366    view: &IndexedNetwork,
367    opts: &BuildOptions,
368    ybus_cache: &mut Option<YbusCache>,
369) -> Result<sprs::CsMat<f64>> {
370    if ybus_cache.as_ref().is_none_or(|c| c.b.is_none()) {
371        fill_ybus_cache(view, opts, ybus_cache)?;
372    }
373    let b = ybus_cache
374        .as_mut()
375        .and_then(|c| c.b.take())
376        .expect("Ybus cache was just filled, so the imaginary part must be present");
377    Ok(negate_into(b))
378}
379
380/// Build the square matrix for one [`MatrixKind`] from an indexed network. The
381/// single dispatch shared by the [`Pipeline`], the `verify` CLI command, and the
382/// TUI inspect screen, so the `YbusB = -Im(Y_bus)` sign lives in one place.
383pub fn build_kind(
384    view: &IndexedNetwork,
385    kind: MatrixKind,
386    opts: &BuildOptions,
387) -> Result<sprs::CsMat<f64>> {
388    match kind {
389        MatrixKind::BPrime => build_bprime(view, opts),
390        MatrixKind::BDoublePrime => build_bdoubleprime(view, opts),
391        MatrixKind::YbusG => build_ybus(view, opts).map(|p| p.g),
392        MatrixKind::YbusB => build_ybus(view, opts).map(|p| negate_into(p.b)),
393        MatrixKind::Lacpf => build_lacpf(view, opts),
394        MatrixKind::Adjacency => build_adjacency(view),
395    }
396}
397
398pub fn zero_impedance_rule_for_kind(
399    kind: MatrixKind,
400    opts: &BuildOptions,
401) -> Option<ZeroImpedanceRule> {
402    match kind {
403        MatrixKind::BPrime => Some(match opts.scheme {
404            crate::matrix::Scheme::Bx => ZeroImpedanceRule::Series,
405            crate::matrix::Scheme::Xb => ZeroImpedanceRule::Reactance,
406        }),
407        MatrixKind::BDoublePrime => Some(match opts.scheme {
408            crate::matrix::Scheme::Bx => ZeroImpedanceRule::Reactance,
409            crate::matrix::Scheme::Xb => ZeroImpedanceRule::Series,
410        }),
411        MatrixKind::YbusG | MatrixKind::YbusB | MatrixKind::Lacpf => {
412            Some(ZeroImpedanceRule::Series)
413        }
414        MatrixKind::Adjacency => None,
415    }
416}
417
418pub fn zero_impedance_skips_for_kind(
419    view: &IndexedNetwork,
420    kind: MatrixKind,
421    opts: &BuildOptions,
422) -> ZeroImpedanceSkips {
423    if !opts.skip_zero_impedance {
424        return ZeroImpedanceSkips::default();
425    }
426    zero_impedance_rule_for_kind(kind, opts).map_or_else(ZeroImpedanceSkips::default, |rule| {
427        skipped_zero_impedance(view, rule)
428    })
429}
430
431pub fn matrix_stats_for_kind(
432    matrix: &sprs::CsMat<f64>,
433    view: &IndexedNetwork,
434    kind: MatrixKind,
435    opts: &BuildOptions,
436) -> MatrixStats {
437    MatrixStats::from_csr(matrix)
438        .with_zero_impedance_skips(zero_impedance_skips_for_kind(view, kind, opts))
439}
440
441fn sha256_hex(bytes: &[u8]) -> String {
442    use sha2::{Digest, Sha256};
443    use std::fmt::Write as _;
444    let mut hasher = Sha256::new();
445    hasher.update(bytes);
446    let digest = hasher.finalize();
447    let mut out = String::with_capacity(digest.len() * 2);
448    for byte in digest {
449        let _ = write!(out, "{byte:02x}");
450    }
451    out
452}
453
454#[cfg(test)]
455mod tests {
456    use super::sanitize_stem;
457    use std::path::Path;
458
459    #[test]
460    fn sanitize_stem_confines_names_to_out_dir() {
461        // Every result must be a single component with no separators, so
462        // `out_dir.join(stem)` cannot escape `out_dir`.
463        for name in [
464            "../../etc/passwd",
465            "/abs/path",
466            "..",
467            ".",
468            "..\\..\\win",
469            "",
470            "a/b/c",
471        ] {
472            let stem = sanitize_stem(name);
473            let joined = Path::new("out").join(&stem);
474            assert_eq!(
475                joined.components().count(),
476                2,
477                "{name:?} -> {stem:?} escaped out_dir as {joined:?}"
478            );
479            assert!(!stem.is_empty());
480            assert!(stem != "." && stem != "..");
481        }
482    }
483
484    #[test]
485    fn sanitize_stem_keeps_ordinary_names() {
486        assert_eq!(sanitize_stem("case118"), "case118");
487        assert_eq!(sanitize_stem("ieee-13_bus.v2"), "ieee-13_bus.v2");
488    }
489
490    #[test]
491    fn sanitize_stem_separates_names_that_sanitize_alike() {
492        // `a/b` must not overwrite `a_b`'s files in a multi case export.
493        assert_ne!(sanitize_stem("a/b"), sanitize_stem("a_b"));
494        assert_ne!(sanitize_stem("a/b"), sanitize_stem("a\\b"));
495        assert_eq!(sanitize_stem("a_b"), "a_b");
496    }
497
498    #[test]
499    fn a_safe_name_cannot_impersonate_a_disambiguated_stem() {
500        // `.foo` sanitizes to `foo` and so carries a suffix. Naming a second
501        // case exactly that suffixed stem would overwrite it if safe names
502        // passed through unconditionally, and the suffix is derived from a
503        // published hash, so the second name takes no search to construct.
504        let unsafe_stem = sanitize_stem(".foo");
505        assert_ne!(sanitize_stem(&unsafe_stem), unsafe_stem);
506        assert_ne!(sanitize_stem(&unsafe_stem), sanitize_stem(".foo"));
507        // A name that merely contains a hyphen is untouched.
508        assert_eq!(sanitize_stem("ieee-13"), "ieee-13");
509        assert_eq!(sanitize_stem("case-deadbeef"), "case-deadbeef");
510    }
511
512    #[test]
513    fn sanitize_stem_applies_windows_filename_rules() {
514        // Trailing dots and reserved device names are invalid on Windows.
515        let trailing = sanitize_stem("case.");
516        assert!(!trailing.ends_with('.'), "{trailing:?}");
517        for name in ["con", "CON", "aux.4", "lpt9"] {
518            let stem = sanitize_stem(name);
519            let pre_dot = stem.split('.').next().unwrap();
520            assert!(
521                !super::WINDOWS_RESERVED
522                    .iter()
523                    .any(|r| pre_dot.eq_ignore_ascii_case(r)),
524                "{name:?} -> {stem:?} is still a reserved device name"
525            );
526        }
527    }
528
529    #[test]
530    fn sanitize_stem_caps_the_length() {
531        let long = "x".repeat(4096);
532        // Stem plus hash suffix stays within a filename component budget.
533        assert!(sanitize_stem(&long).len() <= super::MAX_STEM_LEN);
534    }
535}