powerio-tx 0.10.0

Compiler infrastructure for power systems: parse, convert, validate, and lower grid models
Documentation
//! MATPOWER `.m` case file parser. Standard MATPOWER 7.x format.

mod locate;
mod matlab;
mod rows;
mod tokens;
mod writer;

#[cfg(test)]
mod tests;

pub use writer::write_matpower;
pub(crate) use writer::write_matpower_conversion;

use crate::network::{BalancedNetwork, BalancedNetworkTables, Generator, Hvdc, SourceFormat};
use crate::{Error, Result};

/// Owned-source entry used by the format hub: move the buffer straight into the
/// retained source (no copy) and take `name_hint` (e.g. the file stem) as the
/// network name.
pub(crate) fn parse_matpower_source(
    source: &str,
    name_hint: Option<&str>,
) -> Result<BalancedNetwork> {
    let name = name_hint
        .map(str::to_owned)
        .or_else(|| matpower_function_name(source).map(str::to_owned))
        .unwrap_or_else(|| "case".to_string());
    parse_matpower_named(source, &name)
}

fn matpower_function_name(source: &str) -> Option<&str> {
    for line in source.lines() {
        let line = line.trim_start();
        if !line.starts_with("function") {
            continue;
        }
        let Some((_, rhs)) = line.split_once('=') else {
            continue;
        };
        let rhs = rhs.trim_start();
        let end = rhs
            .find(|c: char| !(c.is_ascii_alphanumeric() || c == '_'))
            .unwrap_or(rhs.len());
        let starts_ident = rhs
            .as_bytes()
            .first()
            .is_some_and(|b| b.is_ascii_alphabetic() || *b == b'_');
        if end > 0 && starts_ident {
            return Some(&rhs[..end]);
        }
    }
    None
}

fn parse_matpower_named(source: &str, name: &str) -> Result<BalancedNetwork> {
    // Locate each assignment's text directly in `source` and build the network
    // from those borrowed slices in one pass; the typed model owns its data, so
    // the borrows end with `located` and the source Arc moves into the network.
    let net = {
        let located = locate::locate_assignments(source);
        build_case(name, |field| {
            located
                .iter()
                .find(|(f, _)| *f == field)
                .map(|(_, full)| *full)
        })?
    };
    // The other format readers validate references; the MATPOWER path must too,
    // or a duplicate or dangling bus id reaches `IndexedNetwork` as silently
    // collapsed aggregates (the dense bus-id map only debug-asserts uniqueness).
    net.check_references("MATPOWER")?;
    Ok(net)
}

/// Build a [`BalancedNetwork`] from a per-field assignment-text accessor `get`, which
/// returns the raw `mpc.<field> = …;` text for a field name. MATPOWER folds
/// demand and shunts onto the bus row; [`rows::bus_row`] splits them back out
/// into the hub's first-class [`Load`](crate::network::Load) /
/// [`Shunt`](crate::network::Shunt). The caller attaches the source afterward.
fn build_case<'a>(name: &str, get: impl Fn(&str) -> Option<&'a str>) -> Result<BalancedNetwork> {
    let base_mva = get("baseMVA")
        .and_then(|raw| matlab::scalar_from_assignment(raw, "baseMVA").transpose())
        .transpose()?
        .ok_or(Error::MissingField("baseMVA"))?;

    let bus_raw = get("bus").ok_or(Error::MissingField("bus"))?;
    let n_bus = estimate_rows(bus_raw);
    let mut buses = Vec::with_capacity(n_bus);
    let mut loads = Vec::with_capacity(n_bus);
    let mut shunts = Vec::with_capacity(n_bus);
    matlab::for_each_matrix_row(bus_raw, "bus", |row, i| {
        let (bus, load, shunt) = rows::bus_row(row, i)?;
        buses.push(bus);
        if let Some(l) = load {
            loads.push(l);
        }
        if let Some(s) = shunt {
            shunts.push(s);
        }
        Ok(())
    })?;

    let branches = parse_rows(
        get("branch").ok_or(Error::MissingField("branch"))?,
        "branch",
        rows::branch_row,
    )?;

    let generators = parse_gens(&get)?;
    let storage = parse_optional(&get, "storage", rows::storage_row)?;
    let areas = parse_optional(&get, "areas", rows::area_row)?;
    let mut hvdc = parse_optional(&get, "dcline", rows::hvdc_row)?;
    attach_dcline_costs(&get, &mut hvdc)?;

    // Bus names live in a `{...}` cell array; pull them (quotes kept) and attach
    // by position when the count matches.
    if let Some(raw) = get("bus_name") {
        let names = locate::parse_string_cell(raw);
        if names.len() == buses.len() {
            for (bus, label) in buses.iter_mut().zip(names) {
                // An empty cell is an unnamed bus, not a bus named "".
                bus.name = (!label.is_empty()).then_some(label);
            }
        }
    }

    Ok(BalancedNetwork::from_tables(BalancedNetworkTables {
        name: name.to_string(),
        base_mva,
        base_frequency: crate::network::DEFAULT_BASE_FREQUENCY,
        geo: None,
        buses: buses.into(),
        loads: loads.into(),
        shunts: shunts.into(),
        branches: branches.into(),
        switches: Vec::new().into(),
        generators: generators.into(),
        storage: storage.into(),
        hvdc: hvdc.into(),
        transformers_3w: Vec::new().into(),
        areas: areas.into(),
        solver: None,
        source_format: SourceFormat::Matpower,
    }))
}

/// A cheap upper-bound row count for an assignment (one `;` per row), used to
/// pre-size the typed vectors so parsing doesn't reallocate as it streams.
/// Capped: each `;` byte would otherwise pre-allocate a full element (~100
/// bytes), letting a small crafted file demand ~100x its size in memory up
/// front. Real cases sit far below the cap (largest vendored case: 13659
/// buses); beyond it the vectors just grow as rows actually parse.
fn estimate_rows(assignment: &str) -> usize {
    const MAX_ROW_HINT: usize = 1 << 20;
    assignment
        .bytes()
        .filter(|&b| b == b';')
        .count()
        .min(MAX_ROW_HINT)
}

/// Stream the rows of one assignment, building a typed `T` per row via `ctor`.
fn parse_rows<T>(
    assignment: &str,
    field: &str,
    ctor: impl Fn(&[f64], usize) -> Result<T>,
) -> Result<Vec<T>> {
    let mut out = Vec::with_capacity(estimate_rows(assignment));
    matlab::for_each_matrix_row(assignment, field, |row, i| {
        out.push(ctor(row, i)?);
        Ok(())
    })?;
    Ok(out)
}

/// Like [`parse_rows`] but for an optional `mpc.<field>` block (empty if absent).
fn parse_optional<'a, T>(
    get: &impl Fn(&str) -> Option<&'a str>,
    field: &str,
    ctor: impl Fn(&[f64], usize) -> Result<T>,
) -> Result<Vec<T>> {
    match get(field) {
        Some(raw) => parse_rows(raw, field, ctor),
        None => Ok(Vec::new()),
    }
}

/// Fold `mpc.dclinecost` into the dcline rows. Same row layout as `mpc.gencost`,
/// one row per dcline in order (MATPOWER's `toggle_dcline` requires full
/// coverage, so unlike `gencost` there is no reactive second block). A line
/// with no usage cost is padded with an all-zero polynomial row, and a zero
/// cost prices the line exactly as no cost term at all — so a zero row reads
/// back as no cost, and write-then-read stays stable for networks whose lines
/// carry none.
fn attach_dcline_costs<'a>(
    get: &impl Fn(&str) -> Option<&'a str>,
    hvdc: &mut [Hvdc],
) -> Result<()> {
    let Some(raw) = get("dclinecost") else {
        return Ok(());
    };
    let costs = parse_rows(raw, "dclinecost", rows::gencost_row)?;
    if costs.len() != hvdc.len() {
        return Err(Error::DcLineCostCountMismatch {
            dclines: hvdc.len(),
            dclinecost: costs.len(),
        });
    }
    for (line, cost) in hvdc.iter_mut().zip(costs) {
        let zero =
            cost.startup == 0.0 && cost.shutdown == 0.0 && cost.coeffs.iter().all(|c| *c == 0.0);
        if !zero {
            line.cost = Some(cost);
        }
    }
    Ok(())
}

/// Parse `mpc.gen` and fold in the active-power block of `mpc.gencost`.
/// Both are optional: a case with only power flow data has neither and gets no gens.
fn parse_gens<'a>(get: &impl Fn(&str) -> Option<&'a str>) -> Result<Vec<Generator>> {
    let Some(raw) = get("gen") else {
        return Ok(Vec::new());
    };
    let mut gens = parse_rows(raw, "gen", rows::gen_row)?;

    // MATPOWER lays the active-power costs first, one row per generator and in
    // the same order; reactive-power costs (if any) follow in a second block.
    if let Some(craw) = get("gencost") {
        let costs = parse_rows(craw, "gencost", rows::gencost_row)?;
        // Reject a count that is neither `n_gen` (active only) nor `2·n_gen`
        // (active + reactive). A per-row defect surfaces as `ShortRow` first.
        let n = gens.len();
        if costs.len() != n && costs.len() != 2 * n {
            return Err(Error::GenCostCountMismatch {
                gens: n,
                gencost: costs.len(),
            });
        }
        // The first `n` rows are the active-power costs in gen order; any
        // reactive-power second block is accepted but not retained.
        for (generator, cost) in gens.iter_mut().zip(costs) {
            generator.cost = Some(cost);
        }
    }

    Ok(gens)
}