|
| 1 | +//! Probe 3 — per-bus inertia (`H`) ingest path. |
| 2 | +//! |
| 3 | +//! Real per-bus inertia constants `H` (seconds of stored rotational energy) are |
| 4 | +//! **not** in the PyPSA/OSM topology (`buffer.rs` flags this). They come from |
| 5 | +//! operator data: an ENTSO-E generation inventory, ESIOS dispatch, or a TSO inertia |
| 6 | +//! estimate, mapped to buses. This module is the **ingest path** — parse a `bus,H` |
| 7 | +//! table when available, align it to the grid's bus order ([`crate::PypsaImport::bus_ids`]), |
| 8 | +//! and **disclose** how many buses fell back to a proxy (mirroring `ingest.rs`'s |
| 9 | +//! `n_estimated_*` honesty). The result feeds [`crate::inertia_buffer_column`]. |
| 10 | +//! |
| 11 | +//! **No external data is bundled** (network/data policy): the caller supplies the |
| 12 | +//! file contents (the lib takes `&str`, like `from_pypsa_csv`), or uses the |
| 13 | +//! deterministic [`proxy_inertia`] for a no-data demo. A proxy is decoupled from |
| 14 | +//! wiring on purpose — the buffer axis is storage, orthogonal to topology, so a |
| 15 | +//! topology-blind proxy is the honest stand-in until measured `H` is wired. |
| 16 | +
|
| 17 | +use std::collections::HashMap; |
| 18 | + |
| 19 | +/// Provenance of an aligned per-bus inertia vector — how many buses carried a |
| 20 | +/// measured `H` vs a proxy fallback. Disclose it the way `PypsaImport` discloses |
| 21 | +/// estimated reactance/limit counts. |
| 22 | +#[derive(Debug, Clone, Copy, PartialEq, Eq)] |
| 23 | +pub struct InertiaProvenance { |
| 24 | + /// Buses whose `H` came from the measured table. |
| 25 | + pub measured: usize, |
| 26 | + /// Buses that fell back to the proxy value. |
| 27 | + pub proxy: usize, |
| 28 | +} |
| 29 | + |
| 30 | +/// Parse a `bus,H` table (CSV, `,` or `;`): a header naming a bus column |
| 31 | +/// (`bus_id` / `name` / `bus` / `bus0`) and an inertia column (`h` / `inertia` / |
| 32 | +/// `inertia_h` / `inertia_s` / `h_s`). Returns `bus_id → H`. Rows with an |
| 33 | +/// unparseable / non-positive `H` are skipped. Unknown columns are ignored, so an |
| 34 | +/// ENTSO-E/ESIOS export with extra fields parses as-is. |
| 35 | +pub fn parse_bus_inertia(text: &str) -> HashMap<String, f64> { |
| 36 | + let mut out = HashMap::new(); |
| 37 | + let mut lines = text.lines().filter(|l| !l.trim().is_empty()); |
| 38 | + let Some(header) = lines.next() else { |
| 39 | + return out; |
| 40 | + }; |
| 41 | + let delim = if header.matches(';').count() > header.matches(',').count() { |
| 42 | + ';' |
| 43 | + } else { |
| 44 | + ',' |
| 45 | + }; |
| 46 | + let cols: Vec<String> = header |
| 47 | + .split(delim) |
| 48 | + .map(|h| h.trim().to_ascii_lowercase()) |
| 49 | + .collect(); |
| 50 | + let find = |names: &[&str]| -> Option<usize> { |
| 51 | + names.iter().find_map(|w| cols.iter().position(|c| c == w)) |
| 52 | + }; |
| 53 | + let (Some(bi), Some(hi)) = ( |
| 54 | + find(&["bus_id", "name", "bus", "bus0"]), |
| 55 | + find(&["h", "inertia", "inertia_h", "inertia_s", "h_s"]), |
| 56 | + ) else { |
| 57 | + return out; |
| 58 | + }; |
| 59 | + for row in lines { |
| 60 | + let f: Vec<&str> = row.split(delim).collect(); |
| 61 | + let (Some(b), Some(h)) = (f.get(bi), f.get(hi)) else { |
| 62 | + continue; |
| 63 | + }; |
| 64 | + let bus = b.trim().to_string(); |
| 65 | + if bus.is_empty() { |
| 66 | + continue; |
| 67 | + } |
| 68 | + if let Ok(val) = h.trim().parse::<f64>() { |
| 69 | + if val > 0.0 { |
| 70 | + out.insert(bus, val); |
| 71 | + } |
| 72 | + } |
| 73 | + } |
| 74 | + out |
| 75 | +} |
| 76 | + |
| 77 | +/// Align a measured `bus → H` map to the grid's bus order, filling missing buses |
| 78 | +/// with `fallback`. Returns the per-bus `H` vector (length `bus_ids.len()`) plus its |
| 79 | +/// [`InertiaProvenance`] so the caller can disclose the proxy fraction. |
| 80 | +pub fn inertia_for_buses( |
| 81 | + bus_ids: &[String], |
| 82 | + measured: &HashMap<String, f64>, |
| 83 | + fallback: f64, |
| 84 | +) -> (Vec<f64>, InertiaProvenance) { |
| 85 | + let mut h = Vec::with_capacity(bus_ids.len()); |
| 86 | + let mut n_measured = 0usize; |
| 87 | + for id in bus_ids { |
| 88 | + match measured.get(id) { |
| 89 | + Some(&v) => { |
| 90 | + n_measured += 1; |
| 91 | + h.push(v); |
| 92 | + } |
| 93 | + None => h.push(fallback), |
| 94 | + } |
| 95 | + } |
| 96 | + ( |
| 97 | + h, |
| 98 | + InertiaProvenance { |
| 99 | + measured: n_measured, |
| 100 | + proxy: bus_ids.len() - n_measured, |
| 101 | + }, |
| 102 | + ) |
| 103 | +} |
| 104 | + |
| 105 | +/// Deterministic per-bus inertia proxy in `[base, base+span]` (SplitMix64-seeded, |
| 106 | +/// decoupled from wiring). Honest no-data stand-in: the buffer axis is storage, so a |
| 107 | +/// topology-blind proxy preserves buffer ⊥ topology. Same `(n, base, span, seed)` ⇒ |
| 108 | +/// same vector. |
| 109 | +pub fn proxy_inertia(n: usize, base: f64, span: f64, seed: u64) -> Vec<f64> { |
| 110 | + let mut state = seed; |
| 111 | + (0..n) |
| 112 | + .map(|_| { |
| 113 | + state = state.wrapping_add(0x9E37_79B9_7F4A_7C15); |
| 114 | + let mut z = state; |
| 115 | + z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9); |
| 116 | + z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB); |
| 117 | + z ^= z >> 31; |
| 118 | + base + span * (z as f64 / u64::MAX as f64) |
| 119 | + }) |
| 120 | + .collect() |
| 121 | +} |
| 122 | + |
| 123 | +#[cfg(test)] |
| 124 | +mod tests { |
| 125 | + use super::*; |
| 126 | + use crate::inertia_buffer_column; |
| 127 | + |
| 128 | + const INERTIA_CSV: &str = "\ |
| 129 | +bus_id,inertia_h,source |
| 130 | +A,5.0,nuclear |
| 131 | +B,2.5,wind |
| 132 | +C,0.0,solar |
| 133 | +D,abc,bad |
| 134 | +E,8.0,hydro |
| 135 | +"; |
| 136 | + |
| 137 | + #[test] |
| 138 | + fn parses_bus_inertia_skipping_bad_rows() { |
| 139 | + let m = parse_bus_inertia(INERTIA_CSV); |
| 140 | + assert_eq!(m.len(), 3, "C (H=0) and D (unparseable) are skipped"); |
| 141 | + assert_eq!(m["A"], 5.0); |
| 142 | + assert_eq!(m["B"], 2.5); |
| 143 | + assert_eq!(m["E"], 8.0); |
| 144 | + assert!(!m.contains_key("C") && !m.contains_key("D")); |
| 145 | + } |
| 146 | + |
| 147 | + #[test] |
| 148 | + fn aligns_to_grid_order_and_discloses_proxy_fill() { |
| 149 | + let m = parse_bus_inertia(INERTIA_CSV); |
| 150 | + // Grid has A,B,X (X has no measured H) → X uses the fallback. |
| 151 | + let bus_ids = vec!["A".to_string(), "B".to_string(), "X".to_string()]; |
| 152 | + let (h, prov) = inertia_for_buses(&bus_ids, &m, 1.0); |
| 153 | + assert_eq!(h, vec![5.0, 2.5, 1.0]); |
| 154 | + assert_eq!(prov.measured, 2); |
| 155 | + assert_eq!(prov.proxy, 1); |
| 156 | + } |
| 157 | + |
| 158 | + #[test] |
| 159 | + fn proxy_is_deterministic_and_bounded() { |
| 160 | + let a = proxy_inertia(64, 2.0, 6.0, 0xC0FFEE); |
| 161 | + let b = proxy_inertia(64, 2.0, 6.0, 0xC0FFEE); |
| 162 | + assert_eq!(a, b, "same seed ⇒ same proxy"); |
| 163 | + assert!( |
| 164 | + a.iter().all(|&h| (2.0..=8.0).contains(&h)), |
| 165 | + "in [base, base+span]" |
| 166 | + ); |
| 167 | + // A different seed gives a different vector (decoupled, not constant). |
| 168 | + assert_ne!(a, proxy_inertia(64, 2.0, 6.0, 0xBEEF)); |
| 169 | + } |
| 170 | + |
| 171 | + #[test] |
| 172 | + fn ingested_h_feeds_the_inertia_buffer_column() { |
| 173 | + let m = parse_bus_inertia(INERTIA_CSV); |
| 174 | + let bus_ids = vec!["A".to_string(), "B".to_string(), "E".to_string()]; |
| 175 | + let (h, _) = inertia_for_buses(&bus_ids, &m, 1.0); |
| 176 | + let col = inertia_buffer_column(&h, 0.2); |
| 177 | + // E (H=8) is max → 1.0, B (H=2.5) is min → 0.0; all normalized. |
| 178 | + assert_eq!(col.len(), 3); |
| 179 | + assert_eq!(col[2], 1.0); |
| 180 | + assert_eq!(col[1], 0.0); |
| 181 | + assert!(col.iter().all(|&c| (0.0..=1.0).contains(&c))); |
| 182 | + } |
| 183 | +} |
0 commit comments