powerio_matrix/
pipeline.rs

1//! Orchestrates a single case → output directory.
2//!
3//! Given a parsed `Network`, 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, build_adjacency, build_bdoubleprime, build_bprime, build_lacpf,
21    build_ybus, negate_into, sddm_check,
22};
23use crate::network::Network;
24
25#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
26#[non_exhaustive]
27pub enum MatrixKind {
28    /// FDPF B' (shuntless, taps=1, shifts=0).
29    BPrime,
30    /// FDPF B'' (with shunts, taps; r=0 if BX scheme).
31    BDoublePrime,
32    /// `Re(Y_bus)` — full conductance matrix.
33    YbusG,
34    /// `-Im(Y_bus)` — full susceptance Laplacian (positive convention).
35    YbusB,
36    /// LACPF block: `[[G, -B], [-B, -G]]`, 2n × 2n indefinite.
37    Lacpf,
38    /// 0/1 bus adjacency matrix.
39    Adjacency,
40}
41
42impl MatrixKind {
43    pub const ALL: &'static [MatrixKind] = &[
44        Self::BPrime,
45        Self::BDoublePrime,
46        Self::YbusG,
47        Self::YbusB,
48        Self::Lacpf,
49        Self::Adjacency,
50    ];
51
52    pub fn slug(self) -> &'static str {
53        match self {
54            Self::BPrime => "bprime",
55            Self::BDoublePrime => "bdoubleprime",
56            Self::YbusG => "ybus_real",
57            Self::YbusB => "ybus_imag",
58            Self::Lacpf => "lacpf",
59            Self::Adjacency => "adjacency",
60        }
61    }
62
63    pub fn label(self) -> &'static str {
64        match self {
65            Self::BPrime => "B' (FDPF, shuntless)",
66            Self::BDoublePrime => "B'' (FDPF, with shunts)",
67            Self::YbusG => "Re(Y_bus)",
68            Self::YbusB => "-Im(Y_bus)",
69            Self::Lacpf => "LACPF block (2n×2n)",
70            Self::Adjacency => "adjacency (0/1)",
71        }
72    }
73}
74
75/// How to populate the RHS vector(s) emitted alongside each matrix.
76#[derive(Debug, Clone, Copy, PartialEq, Eq, Default, Serialize, Deserialize)]
77pub enum RhsKind {
78    #[default]
79    None,
80    /// Zero-mean Gaussian random (deterministic from `rng_seed`).
81    Random,
82    /// Power injections from the case: `b = (Pd, Qd) / baseMVA`.
83    Injection,
84}
85
86#[derive(Debug, Clone)]
87pub struct Pipeline {
88    pub matrices: Vec<MatrixKind>,
89    pub options: BuildOptions,
90    pub rhs: RhsKind,
91    pub rng_seed: u64,
92    pub source_file: Option<PathBuf>,
93}
94
95impl Default for Pipeline {
96    fn default() -> Self {
97        Self {
98            matrices: vec![MatrixKind::BPrime],
99            options: BuildOptions::default(),
100            rhs: RhsKind::None,
101            rng_seed: 0x00C0_FFEE,
102            source_file: None,
103        }
104    }
105}
106
107#[derive(Debug, Clone)]
108pub struct PipelineOutputs {
109    pub case_name: String,
110    pub files: Vec<PathBuf>,
111    pub metadata: CaseMetadata,
112}
113
114impl Pipeline {
115    pub fn run(&self, net: &Network, out_dir: impl AsRef<Path>) -> Result<PipelineOutputs> {
116        let out_dir = out_dir.as_ref();
117        std::fs::create_dir_all(out_dir)?;
118
119        let view = IndexedNetwork::new(net);
120
121        let mut files = Vec::new();
122        let mut matrices_meta = Vec::new();
123
124        for &kind in &self.matrices {
125            let matrix_path = out_dir.join(format!("{}_{}.mtx", view.name(), kind.slug()));
126            let matrix = self.build(&view, kind)?;
127            write_mtx(&matrix, &matrix_path)?;
128            let stats = MatrixStats::from_csr(&matrix);
129            let sddm = sddm_check(&matrix);
130            matrices_meta.push(MatrixMetadata {
131                kind: kind.slug().to_string(),
132                file: matrix_path
133                    .file_name()
134                    .and_then(|s| s.to_str())
135                    .unwrap_or("")
136                    .to_string(),
137                stats,
138                sddm,
139            });
140            files.push(matrix_path);
141
142            // RHS for matrices that take a RHS of length n (skip LACPF which is 2n).
143            if let Some(rhs) = self.build_rhs(&view, kind) {
144                let rhs_path = out_dir.join(format!("{}_{}_rhs.mtx", view.name(), kind.slug()));
145                write_vector_mtx(&rhs, &rhs_path)?;
146                files.push(rhs_path);
147            }
148        }
149
150        // Shunt vector as a sidecar (not always meaningful, but cheap).
151        let shunt_path = out_dir.join(format!("{}_shunt.mtx", view.name()));
152        let base = view.per_unit_base();
153        let shunt: Vec<f64> = view.bs().iter().map(|&b| b / base).collect();
154        write_vector_mtx(&shunt, &shunt_path)?;
155        files.push(shunt_path);
156
157        let metadata = CaseMetadata {
158            case_name: view.name().to_string(),
159            source_file: self
160                .source_file
161                .as_ref()
162                .and_then(|p| p.to_str())
163                .map(str::to_string),
164            source_sha256: self
165                .source_file
166                .as_ref()
167                .and_then(|p| std::fs::read(p).ok())
168                .map(|b| sha256_hex(&b)),
169            base_mva: view.base_mva(),
170            n_buses: view.n(),
171            n_branches: view.branches().len(),
172            build_options: self.options.clone(),
173            matrices: matrices_meta,
174            powerio_version: env!("CARGO_PKG_VERSION").to_string(),
175        };
176        let meta_path = out_dir.join(format!("{}_meta.json", view.name()));
177        write_meta_json(&metadata, &meta_path)?;
178        files.push(meta_path);
179
180        Ok(PipelineOutputs {
181            case_name: view.name().to_string(),
182            files,
183            metadata,
184        })
185    }
186
187    fn build(&self, case: &IndexedNetwork, kind: MatrixKind) -> Result<sprs::CsMat<f64>> {
188        build_kind(case, kind, &self.options)
189    }
190
191    fn build_rhs(&self, case: &IndexedNetwork, kind: MatrixKind) -> Option<Vec<f64>> {
192        // No meaningful RHS for the 2n LACPF block or the structural adjacency.
193        if matches!(self.rhs, RhsKind::None)
194            || matches!(kind, MatrixKind::Lacpf | MatrixKind::Adjacency)
195        {
196            return None;
197        }
198        let n = case.n();
199        let v = match self.rhs {
200            RhsKind::Random => {
201                let mut rng = ChaCha8Rng::seed_from_u64(self.rng_seed.wrapping_add(kind as u64));
202                let dist = StandardUniform;
203                let mut v: Vec<f64> = (0..n)
204                    .map(|_| {
205                        let u: f64 = dist.sample(&mut rng);
206                        u - 0.5
207                    })
208                    .collect();
209                let mean = v.iter().sum::<f64>() / n as f64;
210                for x in &mut v {
211                    *x -= mean; // zero-mean for Laplacian compatibility
212                }
213                v
214            }
215            RhsKind::Injection => {
216                let base = case.per_unit_base();
217                match kind {
218                    MatrixKind::BPrime | MatrixKind::YbusG | MatrixKind::YbusB => {
219                        case.pd().iter().map(|&p| -p / base).collect()
220                    }
221                    MatrixKind::BDoublePrime => case.qd().iter().map(|&q| -q / base).collect(),
222                    MatrixKind::Lacpf | MatrixKind::Adjacency => unreachable!(),
223                }
224            }
225            RhsKind::None => unreachable!(),
226        };
227        Some(v)
228    }
229}
230
231/// Build the square matrix for one [`MatrixKind`] from an indexed network. The
232/// single dispatch shared by the [`Pipeline`], the `verify` CLI command, and the
233/// TUI inspect screen, so the `YbusB = -Im(Y_bus)` sign lives in one place.
234pub fn build_kind(
235    view: &IndexedNetwork,
236    kind: MatrixKind,
237    opts: &BuildOptions,
238) -> Result<sprs::CsMat<f64>> {
239    match kind {
240        MatrixKind::BPrime => build_bprime(view, opts),
241        MatrixKind::BDoublePrime => build_bdoubleprime(view, opts),
242        MatrixKind::YbusG => build_ybus(view, opts).map(|p| p.g),
243        MatrixKind::YbusB => build_ybus(view, opts).map(|p| negate_into(p.b)),
244        MatrixKind::Lacpf => build_lacpf(view, opts),
245        MatrixKind::Adjacency => build_adjacency(view),
246    }
247}
248
249fn sha256_hex(bytes: &[u8]) -> String {
250    use sha2::{Digest, Sha256};
251    use std::fmt::Write as _;
252    let mut hasher = Sha256::new();
253    hasher.update(bytes);
254    let digest = hasher.finalize();
255    let mut out = String::with_capacity(digest.len() * 2);
256    for byte in digest {
257        let _ = write!(out, "{byte:02x}");
258    }
259    out
260}
powerio_matrix/pipeline.rs

powerio_matrix/
pipeline.rs
