diff --git a/uq-usecase/notebooks/m_viz.ipynb b/uq-usecase/notebooks/m_viz.ipynb new file mode 100644 index 000000000..07264453c --- /dev/null +++ b/uq-usecase/notebooks/m_viz.ipynb @@ -0,0 +1,850 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "4ac1d3fd", + "metadata": {}, + "source": [ + "# m_viz: MATPOWER case geovisualization\n", + "* conda: h-py312-basic\n", + "\n", + "Parses MATPOWER `.m` case files (ACOPF or PF solutions) into normalized pandas DataFrames\n", + "(bus, gen, branch) and produces interactive Plotly geo maps. Geographic coordinates come from a\n", + "`.gic` or `.AUX` file. Generator fuel labels and multi-circuit branch data are read from the `.m`\n", + "file.\n", + "\n", + "The map visualizes the solution stored in the `.m` file: bus loads, dispatched generation (marker\n", + "size proportional to MW), and branch loading percentage (viridis color scale).\n", + "\n", + "Optionally, point `GRIDKIT_REPO` to a local GridKit repository clone to augment hover labels with\n", + "GridKit-assigned bus and branch IDs.\n", + "\n", + "Key steps:\n", + "- **Colocated bus splitting**: buses sharing identical lat/lon are spread on a small circle so each is individually hoverable.\n", + "- **Generator fan-out**: when a bus has more than one generator, all generators at that bus are spread radially so each marker is individually hoverable; thin connector lines link them back to the bus. Optionally, single-generator buses can also be offset (`FAN_SINGLE_GEN = True`) to keep the bus marker visible underneath.\n", + "- **Fault bus overlay**: optionally highlights a faulted bus with a red marker.\n", + "- **Branch loading coloring**: branches colored by loading percentage when flow data is available.\n", + "\n", + "Currently supported cases: **Hawaii40**, **Illinois (ACTIVSg200)**, **Texas (ACTIVSg2000)**, **WECC (ACTIVSg10k)**.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "49ea1a74", + "metadata": {}, + "outputs": [], + "source": [ + "from pathlib import Path\n", + "import importlib\n", + "import os\n", + "import sys\n", + "\n", + "# py-utils is one level up from the notebooks/ directory\n", + "_py_utils = Path.cwd().parent / \"py-utils\"\n", + "if str(_py_utils) not in sys.path:\n", + " sys.path.insert(0, str(_py_utils))\n", + "\n", + "import pandas as pd\n", + "\n", + "import m_viz_utils\n", + "from m_viz_utils import read_matpower_case, summarize_case, validation_report_df" + ] + }, + { + "cell_type": "markdown", + "id": "0834cfda", + "metadata": {}, + "source": [ + "## environment and display setup\n", + "\n", + "Keep the notebook lightweight for now. We only need enough setup to inspect the `.m` case tables cleanly." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e3e9c781", + "metadata": {}, + "outputs": [], + "source": [ + "from IPython.core.interactiveshell import InteractiveShell\n", + "\n", + "InteractiveShell.ast_node_interactivity = \"all\"\n", + "\n", + "pd.set_option(\"display.max_rows\", 10)\n", + "pd.set_option(\"display.max_columns\", 100)\n", + "pd.set_option(\"display.max_colwidth\", 80)\n", + "\n", + "onkestrel = \"NREL_CLUSTER\" in os.environ and os.environ[\"NREL_CLUSTER\"] == \"kestrel\"\n", + "print(f\"onkestrel: {onkestrel}\")" + ] + }, + { + "cell_type": "markdown", + "id": "adefd4f0", + "metadata": {}, + "source": [ + "# configure input paths: select case\n", + "\n", + "Download TAMU synthetic grid cases from\n", + "[Texas A&M Electric Grid Test Cases](https://electricgrids.engr.tamu.edu/electric-grid-test-cases/).\n", + "After unzipping, point `CASE_DATA_DIR` to the folder containing the extracted files.\n", + "\n", + "Files used (Hawaii40 example):\n", + "\n", + "| File | Used for |\n", + "|---|---|\n", + "| `Hawaii40_20231026.m` | MATPOWER case (buses, generators, branches) |\n", + "| `Hawaii40_GIC_data.gic` | Geographic coordinates (preferred) |\n", + "| `Hawaii40_20231026.AUX` | Geographic coordinates (fallback if no `.gic`) |\n", + "\n", + "Files used (ACTIVSg200 / Illinois example):\n", + "\n", + "| File | Used for |\n", + "|---|---|\n", + "| `case_ACTIVSg200.m` | MATPOWER case (buses, generators, branches) |\n", + "| `ACTIVSg200_GIC_data.gic` | Geographic coordinates (preferred) |\n", + "| `ACTIVSg200.AUX` | Geographic coordinates (fallback if no `.gic`) |\n", + "\n", + "**To use your own data:** set `CASE_DATA_DIR` to your extracted folder and `CASE_NAME` to match. Everything else is auto-detected.\n", + "\n", + "Optionally, set `GRIDKIT_REPO` to the root of a local GridKit repository clone to augment hover labels with GridKit-assigned bus and branch IDs. JSON case files are read from `GRIDKIT_REPO/examples/PhasorDynamics/`. Set `GRIDKIT_REPO = None` to skip this augmentation.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "837e2545", + "metadata": {}, + "outputs": [], + "source": [ + "# ── USER CONFIGURATION ──────────────────────────────────────────────────────\n", + "# Set CASE_NAME, CASE_DATA_DIR, and CASE_JSON_PATH. Everything else is auto-detected.\n", + "\n", + "# CASE_NAME = \"Hawaii40\"\n", + "# CASE_NAME = \"ACTIVSg200\"\n", + "CASE_NAME = \"ACTIVSg2000\"\n", + "# CASE_NAME = \"ACTIVSg10k\"\n", + "\n", + "# CASE_DATA_DIR: folder containing your extracted .m, .gic, and .AUX files.\n", + "# Uncomment and set your own path:\n", + "# CASE_DATA_DIR = Path(\"/path/to/my/ACTIVSg200\")\n", + "\n", + "# My paths (Kestrel / Mac) — comment out if setting CASE_DATA_DIR above:\n", + "if onkestrel:\n", + " CASE_DATA_DIR = (\n", + " Path(\"/kfs2/projects/scidac/scidac-data\") / CASE_NAME / \"raw-tamu-data\"\n", + " )\n", + "else:\n", + " CASE_DATA_DIR = (\n", + " Path(\"/Users/isatkaus/projects/scidac/scidac-data\")\n", + " / CASE_NAME\n", + " / \"raw-tamu-data\"\n", + " )\n", + "\n", + "# GRIDKIT_REPO: root of your local GridKit repository clone.\n", + "# Optional — set to None to skip GridKit JSON ID augmentation on hover labels.\n", + "# JSON case files live inside the GridKit repo under examples/PhasorDynamics/.\n", + "GRIDKIT_REPO = Path(\"/home/isatkaus/gridkit\")\n", + "\n", + "_case_json_paths = (\n", + " {\n", + " \"Hawaii40\": GRIDKIT_REPO / \"examples/PhasorDynamics/Medium/Hawaii/hawaii.json\",\n", + " \"ACTIVSg200\": GRIDKIT_REPO\n", + " / \"examples/PhasorDynamics/Large/Illinois/illinois.json\",\n", + " \"ACTIVSg2000\": GRIDKIT_REPO / \"examples/PhasorDynamics/Large/Texas/texas.json\",\n", + " \"ACTIVSg10k\": GRIDKIT_REPO / \"examples/PhasorDynamics/Large/WECC/wecc.json\",\n", + " }\n", + " if GRIDKIT_REPO is not None\n", + " else {}\n", + ")\n", + "CASE_JSON_PATH = _case_json_paths.get(CASE_NAME, None)\n", + "# ────────────────────────────────────────────────────────────────────────────\n", + "\n", + "case_name = CASE_NAME\n", + "case_data_dir = CASE_DATA_DIR\n", + "\n", + "# File detection: set to None for auto-detect, or specify explicitly.\n", + "# Auto-detect looks for case_.m first, then any single .m in case_data_dir.\n", + "m_file_name = None\n", + "# Hawaii explicit example: m_file_name = \"Hawaii40_20231026.m\"\n", + "gic_file_name = f\"{case_name}_GIC_data.gic\"\n", + "aux_file_name = None # auto: .AUX/.aux or .AUX/.aux\n", + "\n", + "if m_file_name is None:\n", + " default_m = case_data_dir / f\"case_{case_name}.m\"\n", + " if default_m.exists():\n", + " m_file_path = default_m\n", + " else:\n", + " m_candidates = sorted(case_data_dir.glob(\"*.m\"))\n", + " if len(m_candidates) == 1:\n", + " m_file_path = m_candidates[0]\n", + " else:\n", + " raise FileNotFoundError(\n", + " f\"Could not uniquely determine .m file in {case_data_dir}. \"\n", + " f\"Candidates: {[p.name for p in m_candidates]}\"\n", + " )\n", + "else:\n", + " m_file_path = case_data_dir / m_file_name\n", + "\n", + "gic_file_path = case_data_dir / gic_file_name\n", + "\n", + "if aux_file_name is None:\n", + " # case-insensitive: try .AUX then .aux\n", + " def _find_aux(stem, directory):\n", + " for ext in (\".AUX\", \".aux\"):\n", + " p = directory / f\"{stem}{ext}\"\n", + " if p.exists():\n", + " return p\n", + " return directory / f\"{stem}.AUX\" # canonical non-existent path for warnings\n", + "\n", + " default_aux = _find_aux(m_file_path.stem, case_data_dir)\n", + " if default_aux.exists():\n", + " aux_file_path = default_aux\n", + " else:\n", + " aux_file_path = _find_aux(case_name, case_data_dir)\n", + "else:\n", + " aux_file_path = case_data_dir / aux_file_name\n", + "\n", + "gic_exists = gic_file_path.exists()\n", + "aux_exists = aux_file_path.exists()\n", + "\n", + "geo_file_path = gic_file_path if gic_exists else aux_file_path\n", + "\n", + "if not m_file_path.exists():\n", + " raise FileNotFoundError(f\"MATPOWER .m file not found: {m_file_path}\")\n", + "if not gic_exists:\n", + " if aux_exists:\n", + " print(f\"INFO: GIC not found; using AUX for geo: {aux_file_path}\")\n", + " else:\n", + " print(f\"WARNING: GIC file not found at {gic_file_path}\")\n", + " print(f\"WARNING: AUX file not found at {aux_file_path}\")\n", + " print(\"WARNING: No geo source file found (.gic or .AUX)\")\n", + "\n", + "print(f\"case_data_dir: {case_data_dir}\")\n", + "print(f\"m_file_path: {m_file_path}\")\n", + "print(f\"gic_file_path: {gic_file_path} (exists: {gic_exists})\")\n", + "print(f\"aux_file_path: {aux_file_path} (exists: {aux_exists})\")\n", + "print(f\"geo_file_path: {geo_file_path}\")\n", + "print(\n", + " f\"case_json_path: {CASE_JSON_PATH} (exists: {CASE_JSON_PATH is not None and CASE_JSON_PATH.exists()})\"\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "1c0ab4b7", + "metadata": {}, + "source": [ + "## load MATPOWER case\n", + "\n", + "Read the `.m` file into normalized `bus` / `gen` / `branch` tables via `m_viz_utils.py`." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "aeed51df", + "metadata": {}, + "outputs": [], + "source": [ + "importlib.reload(m_viz_utils)\n", + "from m_viz_utils import read_matpower_case, summarize_case, validation_report_df\n", + "\n", + "case_data = read_matpower_case(m_file_path)\n", + "summary = summarize_case(case_data)\n", + "validation_df = validation_report_df(case_data)\n", + "\n", + "summary" + ] + }, + { + "cell_type": "markdown", + "id": "aa687e44", + "metadata": {}, + "source": [ + "## validation\n", + "\n", + "These checks confirm that generator and branch bus references map back to valid `BUS_I` values." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9eaa4859", + "metadata": {}, + "outputs": [], + "source": [ + "validation_df" + ] + }, + { + "cell_type": "markdown", + "id": "0f021800", + "metadata": {}, + "source": [ + "## quick previews\n", + "\n", + "Preview the normalized tables before adding geo joins or plotting logic." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "591a4f1b", + "metadata": {}, + "outputs": [], + "source": [ + "pd.set_option(\"display.max_rows\", 6)\n", + "print(\"bus table\")\n", + "case_data.bus\n", + "\n", + "print(\"gen table\")\n", + "case_data.gen\n", + "\n", + "print(\"branch table\")\n", + "case_data.branch" + ] + }, + { + "cell_type": "markdown", + "id": "225d6abb", + "metadata": {}, + "source": [ + "## optional auxiliary tables\n", + "\n", + "If the MATPOWER file exposes generator fuel labels or cost tables, inspect them here." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b0852e8e", + "metadata": {}, + "outputs": [], + "source": [ + "if case_data.genfuel is not None:\n", + " print(\"genfuel index\")\n", + " case_data.genfuel\n", + "\n", + "if case_data.gencost is not None:\n", + " print(\"gencost table\")\n", + " case_data.gencost.head()" + ] + }, + { + "cell_type": "markdown", + "id": "5b92bd19", + "metadata": {}, + "source": [ + "## geo merge and coordinate loading\n", + "\n", + "Loads bus coordinates from `.gic` or `.AUX`, merges into the MATPOWER tables, splits colocated buses, and fans out generators that share a bus.\n", + "\n", + "The geo merge result is passed to the geo plot cell below.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "88047800", + "metadata": {}, + "outputs": [], + "source": [ + "importlib.reload(m_viz_utils)\n", + "import plotly.graph_objects as go\n", + "from m_viz_utils import attach_geo_to_case\n", + "\n", + "# ── USER CONFIGURATION ──────────────────────────────────────────────────────\n", + "# Circle radius (degrees lat/lon, a distance) for spreading colocated buses.\n", + "# Per-case defaults; override by setting SPLIT_RADIUS to a fixed value instead.\n", + "_split_radius_defaults = {\n", + " \"Hawaii40\": 0.006,\n", + " \"ACTIVSg200\": 0.006,\n", + " \"ACTIVSg2000\": 0.02,\n", + " \"ACTIVSg10k\": 0.04,\n", + "}\n", + "SPLIT_RADIUS = _split_radius_defaults.get(case_name, 0.006)\n", + "\n", + "# Circle radius (degrees lat/lon, a distance) for spreading colocated buses.\n", + "# Per-case defaults; override by setting GEN_FANOUT_RADIUS to a fixed value instead.\n", + "_fanout_radius_defaults = {\n", + " \"Hawaii40\": 0.0025,\n", + " \"ACTIVSg200\": 0.004,\n", + " \"ACTIVSg2000\": 0.008,\n", + " \"ACTIVSg10k\": 0.015,\n", + "}\n", + "GEN_FANOUT_RADIUS = _fanout_radius_defaults.get(case_name, 0.004)\n", + "\n", + "# FAN_SINGLE_GEN: offset single-generator buses by GEN_FANOUT_RADIUS (eastward)\n", + "# so the bus marker underneath remains individually hoverable.\n", + "# False — single-gen buses stay at bus center (default; consistent across cases)\n", + "# True — single-gen buses also get a small eastward offset\n", + "FAN_SINGLE_GEN = False\n", + "# ────────────────────────────────────────────────────────────────────────────\n", + "\n", + "# Fanout is enabled automatically when any bus has more than one generator.\n", + "n_multi_gen_buses = int((case_data.gen.groupby(\"GEN_BUS\").size() > 1).sum())\n", + "enable_gen_fanout = n_multi_gen_buses > 0 or FAN_SINGLE_GEN\n", + "split_radius = SPLIT_RADIUS\n", + "gen_fanout_radius = GEN_FANOUT_RADIUS\n", + "fan_single_gen = FAN_SINGLE_GEN\n", + "\n", + "if geo_file_path.exists():\n", + " geo_result = attach_geo_to_case(\n", + " case_data,\n", + " geo_file_path,\n", + " split_colocated=True,\n", + " split_radius=split_radius,\n", + " gen_fanout=enable_gen_fanout,\n", + " gen_fanout_radius=gen_fanout_radius,\n", + " gen_fanout_singles=fan_single_gen,\n", + " )\n", + " geo_case_data = geo_result.case_data\n", + "\n", + " # --- optional: attach GridKit JSON ids to hover text ---\n", + " import gridkit_utils\n", + "\n", + " importlib.reload(gridkit_utils)\n", + " from gridkit_utils import attach_json_ids\n", + "\n", + " case_json_path = CASE_JSON_PATH\n", + " if case_json_path is not None and case_json_path.exists():\n", + " attach_json_ids(geo_case_data, case_json_path)\n", + " print(f\"JSON ids attached from {case_json_path.name}\")\n", + " else:\n", + " print(\"JSON ids not attached (CASE_JSON_PATH is None or file not found)\")\n", + "\n", + " pd.Series(\n", + " {\n", + " \"geo_source\": str(geo_file_path),\n", + " \"split_applied\": geo_result.split_applied,\n", + " \"split_radius\": split_radius,\n", + " \"gen_fanout_enabled\": enable_gen_fanout,\n", + " \"gen_fanout_radius\": gen_fanout_radius,\n", + " \"fan_single_gen\": fan_single_gen,\n", + " \"n_buses_with_multiple_gens\": n_multi_gen_buses,\n", + " \"unique_bus_locations_before\": geo_result.n_unique_bus_locations_before,\n", + " \"unique_bus_locations_after\": geo_result.n_unique_bus_locations_after,\n", + " }\n", + " )\n", + "else:\n", + " geo_result = None\n", + " geo_case_data = None\n", + " print(\n", + " \"Geo prep skipped: no geo source file found. \"\n", + " \"Set geo_file_path to a valid .gic or .AUX file, then rerun this cell.\"\n", + " )" + ] + }, + { + "cell_type": "markdown", + "id": "2427473c", + "metadata": {}, + "source": [ + "## geo plot\n", + "\n", + "Map with branches, buses (sized by PD), and generators (fuel-colored when available).\n", + "\n", + "Branch **line loading** (%) = max(|S_F|, |S_T|) / RATE_A × 100, where S_F = √(PF²+QF²) and S_T = √(PT²+QT²). Requires PF, QF, PT, QT, and RATE_A columns in the `.m` branch table. Branches with RATE_A = 0 are shown in gray.\n", + "\n", + "### Fault bus selection (optional)\n", + "\n", + "Set `fault_bus` to highlight a bus with a red marker. Bus faults only for now (branch faults not yet supported). Set to None to skip the fault overlay entirely.\n", + "\n", + "`fault_bus` accepts any of four identifiers:\n", + "\n", + "| Identifier | Source | Example (Hawaii40) | Example (ACTIVSg200) |\n", + "|---|---|---|---|\n", + "| `BUS_I` (int) | `.m` file `mpc.bus` col 1 | `1` | `49` |\n", + "| `bus_name` (str) | `.m` file `mpc.bus_name` | `\"ALOHA138\"` | `\"RANTOUL 2 1\"` |\n", + "| JSON bus number (int) | `\"number\"` field in JSON `Bus` device | `1` (same as `BUS_I`) | `49` (same as `BUS_I`) |\n", + "| JSON bus name (str) | `\"name\"` field in JSON `Bus` device | `\"ALOHA138\"` (same as `bus_name`) | `\"RANTOUL 2 1\"` (same) |\n", + "\n", + "For both Hawaii40 and ACTIVSg200, the JSON `\"number\"` and `\"name\"` fields match the `.m` file's `BUS_I` and `bus_name` exactly, so all four forms are equivalent.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8888c8eb", + "metadata": {}, + "outputs": [], + "source": [ + "importlib.reload(m_viz_utils)\n", + "from m_viz_utils import plot_grid, lookup_fault_bus\n", + "\n", + "# ── USER CONFIGURATION ──────────────────────────────────────────────────────\n", + "# SHOW_LOADING = False # False: uniform gray branches\n", + "SHOW_LOADING = True # True: viridis-colored branches by loading_pct\n", + "\n", + "# SAVE_HTML = False\n", + "SAVE_HTML = True\n", + "\n", + "# SAVE_FIGS_DIR: directory where HTML figures are saved.\n", + "# Uncomment and set your own path:\n", + "# SAVE_FIGS_DIR = Path(\"/path/to/my/figs\")\n", + "\n", + "# My paths (Kestrel / Mac):\n", + "if onkestrel:\n", + " SAVE_FIGS_DIR = Path(\"/home/isatkaus/gridkit/uq-usecase/figs\")\n", + "else:\n", + " SAVE_FIGS_DIR = Path(\"/Users/isatkaus/projects/scidac/scidac-data/figs\")\n", + "\n", + "# FIG_WIDTH / FIG_HEIGHT: figure size in pixels. None = plot_grid defaults (width=780, height=650).\n", + "FIG_WIDTH = 920\n", + "FIG_HEIGHT = None\n", + "\n", + "# MARKER_SCALE: global multiplier for all circle marker sizes (buses and generators).\n", + "# 1.0 — default sizes\n", + "# 0.5 — half size (useful for large cases with many overlapping markers)\n", + "# 2.0 — double size\n", + "MARKER_SCALE = 1.0\n", + "\n", + "# FAULT_BUS: bus to highlight with a red fault marker. Set to None to skip.\n", + "# Bus faults only for now (branch faults not yet supported).\n", + "# Accepts any of four identifiers (all equivalent for Hawaii40 and ACTIVSg200):\n", + "# BUS_I (int) from .m file mpc.bus col 1 e.g. 1 or 49\n", + "# bus_name (str) from .m file mpc.bus_name e.g. \"ALOHA138\" or \"RANTOUL 2 1\"\n", + "# JSON \"number\" Bus device \"number\" field e.g. 1 or 49 (= BUS_I)\n", + "# JSON \"name\" Bus device \"name\" field e.g. \"ALOHA138\" or \"RANTOUL 2 1\" (= bus_name)\n", + "_fault_defaults = {\n", + " \"Hawaii40\": \"ALOHA138\",\n", + " \"ACTIVSg200\": 147, # parametric — all 200 buses are wired with BusFault devices\n", + " \"ACTIVSg2000\": 6349, ### 6349, 7428, 7422\n", + " \"ACTIVSg10k\": 100,\n", + "}\n", + "FAULT_BUS = _fault_defaults.get(case_name, None)\n", + "\n", + "# INCLUDE_PLOTLYJS: controls how plotly.js is included in the saved HTML.\n", + "# True — bundles plotly.js (~3 MB); required for MkDocs iframes and offline viewing (default)\n", + "# \"cdn\" — loads plotly.js from CDN at open time; ~50 KB file, but requires internet access\n", + "INCLUDE_PLOTLYJS = True\n", + "# INCLUDE_PLOTLYJS = \"cdn\"\n", + "# ────────────────────────────────────────────────────────────────────────────\n", + "\n", + "show_loading = SHOW_LOADING\n", + "save_html = SAVE_HTML\n", + "save_figs_dir = SAVE_FIGS_DIR\n", + "fig_width = FIG_WIDTH\n", + "fig_height = FIG_HEIGHT\n", + "marker_scale = MARKER_SCALE\n", + "fault_bus = FAULT_BUS\n", + "\n", + "# Show gen connectors whenever gen fanout was enabled in the geo merge cell.\n", + "show_gen_connectors = enable_gen_fanout\n", + "\n", + "if geo_case_data is not None:\n", + " fig_geo = plot_grid(\n", + " geo_case_data,\n", + " zoom=4,\n", + " show_loading=show_loading,\n", + " show_gen_connectors=show_gen_connectors,\n", + " marker_scale=marker_scale,\n", + " )\n", + " if fig_width is not None or fig_height is not None:\n", + " _ = fig_geo.update_layout(\n", + " width=fig_width if fig_width is not None else fig_geo.layout.width,\n", + " height=fig_height if fig_height is not None else fig_geo.layout.height,\n", + " )\n", + "\n", + " # --- fault bus overlay ---\n", + " if fault_bus is not None:\n", + " bus_row = lookup_fault_bus(geo_case_data.bus, fault_bus)\n", + " if bus_row.empty:\n", + " bus_df = geo_case_data.bus\n", + " print(\n", + " f\"WARNING: fault_bus={fault_bus!r} not found.\\n\"\n", + " f\"Available bus_name values: {bus_df['bus_name'].tolist() if 'bus_name' in bus_df.columns else 'N/A'}\\n\"\n", + " f\"BUS_I range: {bus_df['BUS_I'].min()} – {bus_df['BUS_I'].max()}\"\n", + " )\n", + " else:\n", + " fault_lat = float(bus_row[\"lat\"].iloc[0])\n", + " fault_lon = float(bus_row[\"lon\"].iloc[0])\n", + " fault_bus_i = int(bus_row[\"BUS_I\"].iloc[0])\n", + " fault_label = (\n", + " bus_row[\"bus_name\"].iloc[0]\n", + " if \"bus_name\" in bus_row.columns\n", + " else str(fault_bus_i)\n", + " )\n", + " fault_group = f\"fault_{fault_bus_i}\"\n", + "\n", + " _ = fig_geo.add_trace(\n", + " go.Scattermap(\n", + " mode=\"markers\",\n", + " lat=[fault_lat],\n", + " lon=[fault_lon],\n", + " hoverinfo=\"skip\",\n", + " marker=dict(size=36, color=\"black\", opacity=1.0),\n", + " legendgroup=fault_group,\n", + " showlegend=False,\n", + " )\n", + " )\n", + " _ = fig_geo.add_trace(\n", + " go.Scattermap(\n", + " name=f\"⚡ fault: {fault_label}\",\n", + " mode=\"markers\",\n", + " lat=[fault_lat],\n", + " lon=[fault_lon],\n", + " hovertext=[\n", + " f\"⚡ FAULT
BUS_I: {fault_bus_i}
Name: {fault_label}\"\n", + " ],\n", + " hoverinfo=\"text\",\n", + " marker=dict(size=26, color=\"red\", opacity=1.0),\n", + " legendgroup=fault_group,\n", + " showlegend=True,\n", + " )\n", + " )\n", + "\n", + " if save_html:\n", + " save_figs_dir.mkdir(parents=True, exist_ok=True)\n", + " suffix = \"_loading\" if show_loading else \"\"\n", + " fault_suffix = f\"_fault{fault_bus}\" if fault_bus is not None else \"\"\n", + " html_stem = f\"{case_name}_geo{suffix}{fault_suffix}\"\n", + " html_path = save_figs_dir / f\"{html_stem}.html\"\n", + " # height=None + viewport CSS fills browser window in new tab\n", + " fig_html = go.Figure(fig_geo)\n", + " _ = fig_html.update_layout(height=None, width=None, autosize=True)\n", + " _viewport_css = (\n", + " \"var s=document.createElement('style');\"\n", + " \"s.textContent='html,body{height:100vh;margin:0;padding:0;overflow:hidden;}';\"\n", + " \"document.head.appendChild(s);\"\n", + " )\n", + " fig_html.write_html(\n", + " str(html_path),\n", + " full_html=True,\n", + " include_plotlyjs=INCLUDE_PLOTLYJS,\n", + " post_script=_viewport_css,\n", + " )\n", + " print(f\"figure saved to {html_path}\")\n", + "\n", + " fig_geo\n", + "else:\n", + " print(\n", + " \"Geo plot skipped because geo_case_data is not available. \"\n", + " \"Provide a valid geo_file_path and rerun geo prep cell first.\"\n", + " )" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3a839cd2", + "metadata": {}, + "outputs": [], + "source": [ + "#### check file sizes for notebook and saved HTML\n", + "def _fmt_mb(path):\n", + " mb = Path(path).stat().st_size / 1024**2\n", + " return f\"{mb:.2f} MB ({path})\"\n", + "\n", + "\n", + "nb_path = Path.cwd() / \"m_viz.ipynb\"\n", + "print(_fmt_mb(nb_path))\n", + "\n", + "if save_html and html_path.exists():\n", + " print(_fmt_mb(html_path))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "95258def", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "356d7b66", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "08b5b406", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b45d916d", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "30f8ede3", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a5b9c84a", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "12efbc6d", + "metadata": {}, + "outputs": [], + "source": [ + "# importlib.reload(m_viz_utils)\n", + "# from m_viz_utils import plot_grid, lookup_fault_bus\n", + "\n", + "# # --- plot options ---\n", + "# show_loading = (\n", + "# True # True: viridis-colored branches by loading_pct (requires PF/QF/PT/QT/RATE_A)\n", + "# )\n", + "# # show_loading = False # False: uniform gray branches\n", + "\n", + "# # save_html = False # True: write figure to save_figs_dir as .html\n", + "# save_html = True\n", + "\n", + "# save_figs_dir = scidac_data_dir / \"figs\" # adjust as needed\n", + "# # save_figs_dir = Path(\"/home/isatkaus/projects/scidac/isatkaus/scidac-notebooks/figs\")\n", + "\n", + "# show_gen_connectors = \"hawaii\" in str(case_name).lower()\n", + "\n", + "# # --- fault bus (optional) ---\n", + "# # Identify a faulted bus by BUS_I (int) or bus_name (str, from mpc.bus_name).\n", + "# # Set to None to skip.\n", + "# fault_bus = None\n", + "# # fault_bus = 1 # ACTIVSg200: by BUS_I integer\n", + "# fault_bus = \"ALOHA138\" # Hawaii40: by bus_name string\n", + "\n", + "# if geo_case_data is not None:\n", + "# fig_geo = plot_grid(\n", + "# geo_case_data,\n", + "# zoom=7,\n", + "# show_loading=show_loading,\n", + "# show_gen_connectors=show_gen_connectors,\n", + "# )\n", + "\n", + "# # --- fault bus overlay ---\n", + "# if fault_bus is not None:\n", + "# bus_row = lookup_fault_bus(geo_case_data.bus, fault_bus)\n", + "# if bus_row.empty:\n", + "# bus_df = geo_case_data.bus\n", + "# print(\n", + "# f\"WARNING: fault_bus={fault_bus!r} not found.\\n\"\n", + "# f\"Available bus_name values: {bus_df['bus_name'].tolist() if 'bus_name' in bus_df.columns else 'N/A'}\\n\"\n", + "# f\"BUS_I range: {bus_df['BUS_I'].min()} – {bus_df['BUS_I'].max()}\"\n", + "# )\n", + "# else:\n", + "# fault_lat = float(bus_row[\"lat\"].iloc[0])\n", + "# fault_lon = float(bus_row[\"lon\"].iloc[0])\n", + "# fault_bus_i = int(bus_row[\"BUS_I\"].iloc[0])\n", + "# fault_label = (\n", + "# bus_row[\"bus_name\"].iloc[0]\n", + "# if \"bus_name\" in bus_row.columns\n", + "# else str(fault_bus_i)\n", + "# )\n", + "# fault_group = f\"fault_{fault_bus_i}\"\n", + "\n", + "# _ = fig_geo.add_trace(\n", + "# go.Scattermap(\n", + "# mode=\"markers\",\n", + "# lat=[fault_lat],\n", + "# lon=[fault_lon],\n", + "# hoverinfo=\"skip\",\n", + "# marker=dict(size=36, color=\"black\", opacity=1.0),\n", + "# legendgroup=fault_group,\n", + "# showlegend=False,\n", + "# )\n", + "# )\n", + "# _ = fig_geo.add_trace(\n", + "# go.Scattermap(\n", + "# name=f\"⚡ fault: {fault_label}\",\n", + "# mode=\"markers\",\n", + "# lat=[fault_lat],\n", + "# lon=[fault_lon],\n", + "# hovertext=[\n", + "# f\"⚡ FAULT
BUS_I: {fault_bus_i}
Name: {fault_label}\"\n", + "# ],\n", + "# hoverinfo=\"text\",\n", + "# marker=dict(size=26, color=\"red\", opacity=1.0),\n", + "# legendgroup=fault_group,\n", + "# showlegend=True,\n", + "# )\n", + "# )\n", + "\n", + "# if save_html:\n", + "# save_figs_dir.mkdir(parents=True, exist_ok=True)\n", + "# suffix = \"_loading\" if show_loading else \"\"\n", + "# fault_suffix = f\"_fault{fault_bus}\" if fault_bus is not None else \"\"\n", + "# html_stem = f\"{case_name}_geo{suffix}{fault_suffix}\"\n", + "# html_path = save_figs_dir / f\"{html_stem}.html\"\n", + "# # height=None → 100% of container; fills browser window in new tab,\n", + "# # fits the iframe height in the MkDocs page without internal scrollbars\n", + "# fig_html = go.Figure(fig_geo)\n", + "# _ = fig_html.update_layout(height=None)\n", + "# fig_html.write_html(str(html_path), full_html=True, include_plotlyjs=True)\n", + "# print(f\"figure saved to {html_path}\")\n", + "\n", + "# fig_geo\n", + "# else:\n", + "# print(\n", + "# \"Geo plot skipped because geo_case_data is not available. \"\n", + "# \"Provide a valid geo_file_path and rerun geo prep cell first.\"\n", + "# )" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b565afb9", + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3666fd26", + "metadata": {}, + "outputs": [], + "source": [ + "show_gen_connectors = \"hawaii\" in str(case_name).lower()\n", + "\n", + "# --- figure size (pixels; None = use plot_grid defaults: width=780, height=650) ---\n", + "fig_width = None\n", + "fig_height = None" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "h-py312-basic", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.13" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/uq-usecase/notebooks/setup_env.md b/uq-usecase/notebooks/setup_env.md new file mode 100644 index 000000000..4ad75fda3 --- /dev/null +++ b/uq-usecase/notebooks/setup_env.md @@ -0,0 +1,79 @@ +# Conda environment setup + +The notebooks were developed and tested with Python 3.12. Other versions (3.11+) should also work. + +More notebooks will be added in future PRs. The environment below covers all currently included +notebooks. + +## Packages + +| Channel | Packages | +|---|---| +| conda | `python=3.12 numpy pandas plotly scipy jupyter` | +| pip | `matpowercaseframes==2.1.0` | + +```bash +conda create --prefix ./my-env python=3.12 numpy pandas plotly scipy jupyter --yes +conda activate ./my-env +pip install matpowercaseframes==2.1.0 +``` + +## Notebooks + +### m_viz + +`m_viz.ipynb` parses MATPOWER `.m` case files (ACOPF or PF solutions) into normalized pandas +DataFrames and produces interactive Plotly geo maps. Geographic coordinates come from a `.gic` +or `.AUX` file. The map shows bus loads (PD), dispatched generation (fuel-colored, sized by MW), +and branch loading percentage (viridis color scale). + +Optionally, set `GRIDKIT_REPO` to the root of a local GridKit repository clone to augment hover +labels with GridKit-assigned bus and branch IDs (`gridkit_utils.py`, requires `scipy`). Set +`GRIDKIT_REPO = None` to skip this. + +Supported cases: **Hawaii40**, **Illinois (ACTIVSg200)**, **Texas (ACTIVSg2000)**, **WECC (ACTIVSg10k)**. + +#### Input data + +Download TAMU synthetic grid cases from +[Texas A&M Electric Grid Test Cases](https://electricgrids.engr.tamu.edu/electric-grid-test-cases/). +After unzipping, set `CASE_DATA_DIR` to the folder containing the extracted files and `CASE_NAME` +to match. Everything else is auto-detected. + +##### Hawaii40 + +| File | Used for | +|---|---| +| `Hawaii40_20231026.m` | MATPOWER case (buses, generators, branches) | +| `Hawaii40_GIC_data.gic` | Geographic coordinates (preferred) | +| `Hawaii40_20231026.AUX` | Geographic coordinates (fallback if no `.gic`) | + +##### Illinois (ACTIVSg200) + +| File | Used for | +|---|---| +| `case_ACTIVSg200.m` | MATPOWER case (buses, generators, branches) | +| `ACTIVSg200_GIC_data.gic` | Geographic coordinates (preferred) | +| `ACTIVSg200.AUX` | Geographic coordinates (fallback if no `.gic`) | + +##### Texas (ACTIVSg2000) + +| File | Used for | +|---|---| +| `case_ACTIVSg2000.m` | MATPOWER case (buses, generators, branches) | +| `ACTIVSg2000_GIC_data.gic` | Geographic coordinates (preferred) | +| `ACTIVSg2000_dynamics.AUX` | Geographic coordinates (fallback if no `.gic`) | + +##### WECC (ACTIVSg10k) + +| File | Used for | +|---|---| +| `case_ACTIVSg10k.m` | MATPOWER case (buses, generators, branches) | +| `ACTIVSg10k_GIC_data.gic` | Geographic coordinates | + +#### Supporting files + +| File | Purpose | +|---|---| +| `py-utils/m_viz_utils.py` | MATPOWER `.m` parsing, geo merge, colocated bus splitting, generator fan-out, plotting | +| `py-utils/gridkit_utils.py` | GridKit JSON ID augmentation for hover labels (optional) | diff --git a/uq-usecase/py-utils/gridkit_utils.py b/uq-usecase/py-utils/gridkit_utils.py new file mode 100644 index 000000000..e57545089 --- /dev/null +++ b/uq-usecase/py-utils/gridkit_utils.py @@ -0,0 +1,641 @@ +""" +GridKit helper utilities for use in gridkit_helper.ipynb and related notebooks. +""" + +import copy +import json +import os +import shutil +import subprocess +from numbers import Number + +import numpy as np +import pandas as pd +from scipy.stats import qmc, norm as sp_norm + +# --------------------------------------------------------------------------- +# Monitorable variables by element class +# Source: GridKit docs INPUT_FORMAT.md +# - Bus classes table (init vars + other variables available to monitor) +# - Device classes table (Variables available to monitor column) +# --------------------------------------------------------------------------- +MONITORABLE_VARS_BY_ELEMENT = { + # --- Bus classes --- + "Bus": ["Vr", "Vi", "Vm", "Va"], + # --- Device classes --- + "Branch": ["ir1", "ii1", "im1", "p1", "q1", "ir2", "ii2", "im2", "p2", "q2"], + "Load": ["p", "q"], + "Genrou": ["ir", "ii", "p", "q", "delta", "omega", "speed"], + "Gensal": [ + "ir", + "ii", + "p", + "q", + "delta", + "omega", + "speed", + "Eqp", + "psidp", + "psiqpp", + "psidpp", + "vd", + "vq", + "te", + "id", + "iq", + ], + "GenClassical": ["ir", "ii", "p", "q", "delta", "omega"], + "Tgov1": [], + "Ieeet1": ["efd", "ksat"], + "SexsPti": ["efd"], + "Ieeest": ["vss"], + "BusFault": ["state", "ir", "ii"], +} + + +# --------------------------------------------------------------------------- +# UQ workflow helpers +# --------------------------------------------------------------------------- + + +def generate_samples(param_specs, N, seed=None, method="lhs"): + """ + Sample parameters according to per-spec distributions. + + Parameters + ---------- + param_specs : list of dict, each with keys: + "class" - device class string (e.g. "Genrou") + "id" - device id string (e.g. "genrou_2_1") + "param" - parameter name (e.g. "H") + "dist" - one of "uniform" (default) or "normal" + + For dist="uniform": + Option A — percent: "nominal" + "pct" + lo = nominal * (1 - pct), hi = nominal * (1 + pct) + Option B — fixed: "lo" + "hi" + + For dist="normal": + "mean" and "std" + + N : int, number of samples + seed : int or None, for reproducibility + method : "lhs" — Latin Hypercube (better space-filling, recommended) + "random" — independent random draws per parameter + + Returns + ------- + pd.DataFrame shape (N, len(param_specs)) + columns named "{class}_{id}_{param}", e.g. "Genrou_genrou_2_1_H" + + Notes + ----- + LHS works by drawing N uniform samples stratified across [0,1] per + dimension, then mapping each through the inverse CDF (ppf) of the + requested distribution. This gives correct marginals for both uniform + and normal distributions while maximising coverage of the input space. + """ + cols = [f"{s['class']}_{s['id']}_{s['param']}" for s in param_specs] + d = len(param_specs) + + if method == "lhs": + sampler = qmc.LatinHypercube(d=d, seed=seed) + unit = sampler.random(N) # N × d uniform in (0, 1) + else: + rng = np.random.default_rng(seed) + unit = rng.uniform(0, 1, size=(N, d)) + + data = {} + for j, spec in enumerate(param_specs): + dist = spec.get("dist", "uniform") + u = unit[:, j] + + if dist == "uniform": + if "lo" in spec and "hi" in spec: + lo, hi = spec["lo"], spec["hi"] + else: + lo = spec["nominal"] * (1 - spec["pct"]) + hi = spec["nominal"] * (1 + spec["pct"]) + # uniform ppf: lo + u*(hi-lo) + data[cols[j]] = lo + u * (hi - lo) + + elif dist == "normal": + mean, std = spec["mean"], spec["std"] + # normal ppf via scipy + data[cols[j]] = sp_norm.ppf(u, loc=mean, scale=std) + + else: + raise ValueError(f"Unknown dist '{dist}' in param_spec for {cols[j]}") + + return pd.DataFrame(data) + + +def make_run_dir( + base_case_dir, + run_root, + i, + sample_row, + param_specs, + monitors_by_class, + case_fn=None, + solver_fn=None, + solver_overrides=None, +): + """ + Create run_root/run_{i:03d}/, copy case+solver JSON, patch params and monitors. + + Parameters + ---------- + base_case_dir : str, directory containing the base case files + run_root : str, parent directory for all runs + i : int, run index + sample_row : pd.Series, one row from generate_samples() DataFrame + param_specs : list of dict (same structure as generate_samples) + monitors_by_class : dict mapping device/bus class (lowercase ok) -> list of var names + e.g. {"bus": ["Vm","Va"], "infinite_bus": ["Vm","Va"], + "genrou": ["delta","omega"]} + Classes not listed get their mon arrays removed. + case_fn : str or None; if None, auto-detected (single .case.json in dir) + solver_fn : str or None; if None, auto-detected (single .solver.json in dir) + solver_overrides : dict or None; key-value pairs to patch into the solver JSON. + Supports top-level keys (e.g. "tmax") and "events" as a list. + For "events", each entry must have "index" (0-based position in + the events array) plus whichever fields to overwrite, e.g.: + {"tmax": 20.0, + "events": [ + {"index": 0, "time": 2.0}, + {"index": 1, "time": 2.1}, + ]} + + Returns + ------- + str : path to the new run directory + """ + if case_fn is None: + candidates = [f for f in os.listdir(base_case_dir) if f.endswith(".case.json")] + if len(candidates) == 0: + # fall back: any .json that is not a .solver.json (e.g. hawaii.json) + candidates = [ + f + for f in os.listdir(base_case_dir) + if f.endswith(".json") and not f.endswith(".solver.json") + ] + if len(candidates) != 1: + raise ValueError( + f"Expected 1 case .json in {base_case_dir}, found: {candidates}" + ) + case_fn = candidates[0] + if solver_fn is None: + candidates = [ + f for f in os.listdir(base_case_dir) if f.endswith(".solver.json") + ] + if len(candidates) != 1: + raise ValueError( + f"Expected 1 .solver.json in {base_case_dir}, found: {candidates}" + ) + solver_fn = candidates[0] + + run_dir = os.path.join(run_root, f"run_{i:03d}") + os.makedirs(run_dir, exist_ok=True) + + # load solver, strip testing keys; remove output_file so the simulator + # uses its default (always writes mon.csv; output_file only controls a symlink alias) + with open(os.path.join(base_case_dir, solver_fn)) as f: + solver = json.load(f) + solver.pop("reference_file", None) + solver.pop("error_tolerance", None) + solver.pop("output_file", None) + + # apply solver_overrides (tmax, events, dt, ...) + if solver_overrides: + for key, val in solver_overrides.items(): + if key == "events": + # val is a list of {index, ...fields}; patch each event in-place + for patch in val: + idx = patch["index"] + for k, v in patch.items(): + if k != "index": + solver["events"][idx][k] = v + else: + solver[key] = val + + with open(os.path.join(run_dir, solver_fn), "w") as f: + json.dump(solver, f, indent=2) + + # load and deep-copy case + with open(os.path.join(base_case_dir, case_fn)) as f: + case = copy.deepcopy(json.load(f)) + + # --- patch device params --- + spec_lookup = {(s["class"].lower(), s["id"]): s for s in param_specs} + col_lookup = { + (s["class"].lower(), s["id"]): f"{s['class']}_{s['id']}_{s['param']}" + for s in param_specs + } + for dev in case.get("devices", []): + key = (dev.get("class", "").lower(), dev.get("id", "")) + if key in spec_lookup: + spec = spec_lookup[key] + dev.setdefault("params", {})[spec["param"]] = float( + sample_row[col_lookup[key]] + ) + + # --- patch monitors --- + mon_lower = {k.lower(): v for k, v in monitors_by_class.items()} + for bus in case.get("buses", []): + cls = bus.get("class", "").lower() + mon_vars = mon_lower.get(cls, []) + if mon_vars: + bus["mon"] = mon_vars + elif "mon" in bus: + del bus["mon"] + for dev in case.get("devices", []): + cls = dev.get("class", "").lower() + mon_vars = mon_lower.get(cls, []) + if mon_vars: + dev["mon"] = mon_vars + elif "mon" in dev: + del dev["mon"] + + with open(os.path.join(run_dir, case_fn), "w") as f: + json.dump(case, f, indent=2) + + return run_dir + + +def run_sample(run_dir, runner, solver_fn=None, timeout=300): + """ + Run DynamicSimulation in run_dir. + + Parameters + ---------- + run_dir : str + runner : str, full path to DynamicSimulation executable + solver_fn : str or None; if None auto-detected + timeout : int, seconds + + Returns + ------- + subprocess.CompletedProcess + """ + if solver_fn is None: + candidates = [f for f in os.listdir(run_dir) if f.endswith(".solver.json")] + if len(candidates) != 1: + raise ValueError( + f"Expected 1 .solver.json in {run_dir}, found: {candidates}" + ) + solver_fn = candidates[0] + return subprocess.run( + [runner, solver_fn], + cwd=run_dir, + capture_output=True, + text=True, + timeout=timeout, + ) + + +def collect_and_save(run_root, samples_df, out_path, mon_fn="mon.csv", mode="stacked"): + """ + Read all run_i/mon.csv and save as Parquet in one of two formats. + + Parameters + ---------- + run_root : str, parent dir containing run_000/, run_001/, ... + samples_df : pd.DataFrame from generate_samples(); index = run index 0..N-1 + out_path : str + mode="stacked" → path to a single output .parquet file + mode="per_run" → path to an output *directory*; writes run_NNN.parquet per run + mon_fn : str, monitor output filename (default "mon.csv") + mode : "stacked" (default) or "per_run" + "stacked" — single flat file with run_id + time + signals + param cols + "per_run" — one file per run, wide format (rows=timesteps, cols=signals only); + param values are NOT duplicated — they live in samples.csv + + Returns + ------- + mode="stacked" → pd.DataFrame : run_id | time | | + mode="per_run" → list of str : paths to written run_NNN.parquet files + """ + import pyarrow as pa + import pyarrow.parquet as pq + + if mode == "stacked": + frames = [] + for i, row in samples_df.iterrows(): + mon_path = os.path.join(run_root, f"run_{i:03d}", mon_fn) + if not os.path.exists(mon_path) or os.path.islink(mon_path): + print(f" WARNING: missing {mon_path}, skipping run {i}") + continue + df = pd.read_csv(mon_path) + df = df.rename(columns={df.columns[0]: "time"}) + df.insert(0, "run_id", i) + for col, val in row.items(): + df[col] = val + frames.append(df) + + if not frames: + raise RuntimeError("No mon.csv files found — did the runs complete?") + + combined = pd.concat(frames, ignore_index=True) + os.makedirs(os.path.dirname(os.path.abspath(out_path)), exist_ok=True) + pq.write_table(pa.Table.from_pandas(combined, preserve_index=False), out_path) + print(f"Saved {len(combined)} rows ({len(frames)} runs) -> {out_path}") + return combined + + elif mode == "per_run": + os.makedirs(out_path, exist_ok=True) + written = [] + for i, _row in samples_df.iterrows(): + mon_path = os.path.join(run_root, f"run_{i:03d}", mon_fn) + if not os.path.exists(mon_path) or os.path.islink(mon_path): + print(f" WARNING: missing {mon_path}, skipping run {i}") + continue + df = pd.read_csv(mon_path) + df = df.rename(columns={df.columns[0]: "time"}) + out_file = os.path.join(out_path, f"run_{i:03d}.parquet") + pq.write_table(pa.Table.from_pandas(df, preserve_index=False), out_file) + written.append(out_file) + + if not written: + raise RuntimeError("No mon.csv files found — did the runs complete?") + + print(f"Saved {len(written)} run files -> {out_path}/run_NNN.parquet") + return written + + else: + raise ValueError(f"Unknown mode '{mode}': must be 'stacked' or 'per_run'") + + +# --------------------------------------------------------------------------- +# Dispatch editing helpers (aleatoric UQ prep) +# --------------------------------------------------------------------------- + + +def read_genrou_dispatch(case_json_path): + """ + Read all Genrou default dispatch values from a case JSON. + + Parameters + ---------- + case_json_path : str + Path to case JSON file (e.g. hawaii.json). + + Returns + ------- + pd.DataFrame with columns: + gen_id, bus, p0, q0 + """ + with open(case_json_path, "r") as f: + case = json.load(f) + + rows = [] + for dev in case.get("devices", []): + if dev.get("class") != "Genrou": + continue + params = dev.get("params", {}) + ports = dev.get("ports", {}) + rows.append( + { + "gen_id": dev.get("id"), + "bus": ports.get("bus"), + "p0": params.get("p0", np.nan), + "q0": params.get("q0", np.nan), + } + ) + + df = pd.DataFrame(rows) + if not df.empty: + df = df.sort_values(["bus", "gen_id"]).reset_index(drop=True) + return df + + +def plot_genrou_dispatch(dispatch_df, top_n_labels=None, case_name=None): + """ + Plot Genrou dispatch defaults (p0 and q0). + + Parameters + ---------- + dispatch_df : pd.DataFrame + Output of read_genrou_dispatch(). + top_n_labels : int or None + Deprecated compatibility argument. Labels are now shown by default. + case_name : str or None + Optional case label (e.g. "hawaii", "illinois") shown in plot titles. + """ + import plotly.express as px + + if dispatch_df.empty: + print("No Genrou devices found.") + return + + case_prefix = f"[{case_name}] " if case_name else "" + + fig1 = px.bar( + dispatch_df, + x="gen_id", + y="p0", + title=f"{case_prefix}Genrou default dispatch p0 ({len(dispatch_df)} generators)", + hover_data=["bus", "q0"], + ) + _ = fig1.update_layout(xaxis_title="Genrou id", yaxis_title="p0") + _ = fig1.update_xaxes(type="category", tickangle=70) + _ = fig1.show() + + fig2 = px.bar( + dispatch_df, + x="gen_id", + y="q0", + title=f"{case_prefix}Genrou default dispatch q0 ({len(dispatch_df)} generators)", + hover_data=["bus", "p0"], + ) + _ = fig2.update_layout(xaxis_title="Genrou id", yaxis_title="q0") + _ = fig2.update_xaxes(type="category", tickangle=70) + _ = fig2.show() + + fig3 = px.scatter( + dispatch_df, + x="p0", + y="q0", + text="gen_id", + title=f"{case_prefix}Genrou dispatch operating points (p0, q0)", + hover_data=["bus"], + ) + _ = fig3.update_traces(textposition="top center") + _ = fig3.show() + + +def patch_genrou_dispatch(case_json_path, updates, output_case_path=None): + """ + Patch selected Genrou ids with new p0/q0 and write updated case JSON. + + Parameters + ---------- + case_json_path : str + Input case file path. + updates : dict or list[dict] + Either mapping by id: + {"2_1": {"p0": 0.03, "q0": 0.01}, ...} + or list form: + [{"id": "2_1", "p0": 0.03, "q0": 0.01}, ...] + output_case_path : str or None + Output path. If None, overwrite case_json_path. + + Returns + ------- + pd.DataFrame + before/after table with columns: + gen_id, p0_before, p0_after, q0_before, q0_after + """ + with open(case_json_path, "r") as f: + case = copy.deepcopy(json.load(f)) + + if isinstance(updates, dict): + normalized = [] + for gid, vals in updates.items(): + row = {"id": gid} + row.update(vals or {}) + normalized.append(row) + else: + normalized = list(updates) + + seen = set() + duplicates = set() + for row in normalized: + gid = row.get("id") + if gid in seen: + duplicates.add(gid) + seen.add(gid) + if duplicates: + raise ValueError(f"Duplicate ids in updates: {sorted(duplicates)}") + + by_id = { + d.get("id"): d for d in case.get("devices", []) if d.get("class") == "Genrou" + } + + change_rows = [] + missing_ids = [] + for row in normalized: + gid = row.get("id") + if gid not in by_id: + missing_ids.append(gid) + continue + + dev = by_id[gid] + dev.setdefault("params", {}) + + p0_new = row.get("p0", dev["params"].get("p0")) + q0_new = row.get("q0", dev["params"].get("q0")) + + if p0_new is not None and not isinstance(p0_new, Number): + raise TypeError(f"p0 for {gid} must be numeric, got {type(p0_new)}") + if q0_new is not None and not isinstance(q0_new, Number): + raise TypeError(f"q0 for {gid} must be numeric, got {type(q0_new)}") + + old_p0 = dev["params"].get("p0", np.nan) + old_q0 = dev["params"].get("q0", np.nan) + + if "p0" in row: + dev["params"]["p0"] = float(p0_new) + if "q0" in row: + dev["params"]["q0"] = float(q0_new) + + change_rows.append( + { + "gen_id": gid, + "p0_before": old_p0, + "p0_after": dev["params"].get("p0", np.nan), + "q0_before": old_q0, + "q0_after": dev["params"].get("q0", np.nan), + } + ) + + out_path = output_case_path or case_json_path + with open(out_path, "w") as f: + json.dump(case, f, indent=2) + + changes_df = pd.DataFrame(change_rows) + if not changes_df.empty: + changes_df = changes_df.sort_values("gen_id").reset_index(drop=True) + + print(f"Patched {len(change_rows)} Genrou devices in {out_path}") + if missing_ids: + print(f"WARNING: ids not found: {sorted(missing_ids)}") + + return changes_df + + +def attach_json_ids(case_data, json_path) -> None: + """Attach GridKit JSON identifiers to case_data dataframes in-place. + + Adds columns: + - bus_df: ``json_bus_number`` (int, same as BUS_I — confirms correspondence) + - gen_df: ``json_gen_id`` (str like "2_1"; None for offline/absent gens) + - branch_df: ``json_branch_id`` (str like "BR_1_2_1") + + Parameters + ---------- + case_data : MatpowerCaseData + Object whose ``.bus``, ``.gen``, ``.branch`` DataFrames will be mutated. + json_path : str or Path + Path to a GridKit JSON case file (hawaii.json style). + """ + with open(json_path) as fh: + case = json.load(fh) + + # --- bus: json number == BUS_I (trivially identical, added for hover clarity) --- + case_data.bus["json_bus_number"] = case_data.bus["BUS_I"].astype(int) + + # --- generators: derive json_gen_id from rank within GEN_BUS group --- + # Build set of ids that actually exist in the JSON (for validation / offline marking) + json_gen_ids = {d["id"] for d in case["devices"] if d.get("class") == "Genrou"} + + gen_df = case_data.gen + has_status = "GEN_STATUS" in gen_df.columns + + json_gen_id_col = [] + # rank counter per bus + rank: dict[int, int] = {} + for row in gen_df.itertuples(): + bus = int(row.GEN_BUS) + rank[bus] = rank.get(bus, 0) + 1 + candidate = f"{bus}_{rank[bus]}" + # Mark as None if offline in .m (GEN_STATUS=0) or absent from JSON + if has_status and int(getattr(row, "GEN_STATUS", 1)) == 0: + json_gen_id_col.append(None) + elif candidate not in json_gen_ids: + # absent from JSON even though status=1 — mark None + json_gen_id_col.append(None) + else: + json_gen_id_col.append(candidate) + + case_data.gen["json_gen_id"] = json_gen_id_col + + # --- gen: also store the JSON bus "number" for each gen (== GEN_BUS == BUS_I) --- + # Only set when the gen is present in JSON; None for offline/absent gens. + case_data.gen["json_gen_bus_number"] = [ + int(row.GEN_BUS) if gid is not None else None + for row, gid in zip(gen_df.itertuples(), json_gen_id_col) + ] + + # --- branches: derive json_branch_id from (F_BUS, T_BUS) rank --- + parallel_count: dict[tuple, int] = {} + json_branch_id_col = [] + for row in case_data.branch.itertuples(): + key = (int(row.F_BUS), int(row.T_BUS)) + parallel_count[key] = parallel_count.get(key, 0) + 1 + json_branch_id_col.append(f"BR_{key[0]}_{key[1]}_{parallel_count[key]}") + + case_data.branch["json_branch_id"] = json_branch_id_col + + +def get_case_path_for_editing(test_run_dir, example_dir, case_name): + """ + Resolve case path for edits. + + Prefer the copied case file in test_run_dir. If it does not exist, + fall back to the source/example directory. + """ + run_case = os.path.join(test_run_dir, case_name) + if os.path.exists(run_case): + return run_case + return os.path.join(example_dir, case_name) diff --git a/uq-usecase/py-utils/m_viz_utils.py b/uq-usecase/py-utils/m_viz_utils.py new file mode 100644 index 000000000..ad569770f --- /dev/null +++ b/uq-usecase/py-utils/m_viz_utils.py @@ -0,0 +1,1220 @@ +"""Utilities for loading MATPOWER .m cases in m_viz.ipynb.""" + +from __future__ import annotations + +from dataclasses import dataclass +from io import StringIO +from pathlib import Path +import re +import shlex +from typing import Any + +import numpy as np +import pandas as pd +import plotly.graph_objects as go + +from matpowercaseframes import CaseFrames + + +@dataclass +class MatpowerCaseData: + """Normalized MATPOWER case tables and lightweight metadata.""" + + case_name: str + source_path: str + bus: pd.DataFrame + gen: pd.DataFrame + branch: pd.DataFrame + genfuel: pd.Index | None = None + gencost: pd.DataFrame | None = None + bus_name: pd.Index | None = None # mpc.bus_name cell array, one entry per bus row + base_mva: float | None = None + version: str | None = None + + +@dataclass +class GeoPrepResult: + """Result bundle for geo enrichment and split-point preprocessing.""" + + case_data: MatpowerCaseData + bus_geo_df: pd.DataFrame + split_applied: bool + n_unique_bus_locations_before: int + n_unique_bus_locations_after: int + + +def read_matpower_case(m_file_path: str | Path) -> MatpowerCaseData: + """Load a MATPOWER .m file and return normalized case tables. + + Uses the standard ``matpowercaseframes`` package for the numeric matrices + (bus, gen, branch, gencost) and a built-in curly-brace parser for the + string cell-arrays that the standard package does not handle + (``mpc.genfuel``, ``mpc.gentype``, ``mpc.bus_name``). + """ + + source_path = Path(m_file_path).expanduser().resolve() + if not source_path.exists(): + raise FileNotFoundError(f"MATPOWER case file not found: {source_path}") + + case = CaseFrames(str(source_path)) + + bus_df = normalize_bus_df(case.bus) + gen_df = normalize_gen_df(case.gen) + branch_df = normalize_branch_df(case.branch) + + # Parse curly-brace cell arrays directly from the .m file text. + m_text = source_path.read_text(encoding="utf-8", errors="replace") + genfuel = _parse_m_cell_array(m_text, "genfuel") + bus_name_list = _parse_m_cell_array(m_text, "bus_name") + + genfuel_index = pd.Index(genfuel, name="genfuel") if genfuel else None + bus_name_index = pd.Index(bus_name_list, name="bus_name") if bus_name_list else None + + # Attach bus_name as a column on bus_df so name-based lookups work. + if bus_name_index is not None and len(bus_name_index) == len(bus_df): + bus_df = bus_df.copy() + bus_df["bus_name"] = bus_name_index.values + + return MatpowerCaseData( + case_name=getattr(case, "name", source_path.stem), + source_path=str(source_path), + bus=bus_df, + gen=gen_df, + branch=branch_df, + genfuel=genfuel_index, + gencost=_copy_df_attr(case, "gencost"), + bus_name=bus_name_index, + base_mva=_coerce_optional_float(getattr(case, "baseMVA", None)), + version=_coerce_optional_str(getattr(case, "version", None)), + ) + + +def normalize_bus_df(bus_df: pd.DataFrame) -> pd.DataFrame: + """Return bus table indexed by BUS_I with stable numeric typing.""" + + df = bus_df.copy() + if "BUS_I" not in df.columns: + df = df.reset_index() + if "BUS_I" not in df.columns: + raise ValueError("MATPOWER bus table is missing required column 'BUS_I'") + + df["BUS_I"] = _to_int_series(df["BUS_I"], "BUS_I") + df = df.set_index("BUS_I", drop=False).sort_index() + return df + + +def normalize_gen_df(gen_df: pd.DataFrame) -> pd.DataFrame: + """Return generator table with numeric GEN_BUS and a stable row index.""" + + df = gen_df.copy().reset_index(drop=False) + if "GEN_BUS" not in df.columns: + raise ValueError( + "MATPOWER generator table is missing required column 'GEN_BUS'" + ) + + df["GEN_BUS"] = _to_int_series(df["GEN_BUS"], "GEN_BUS") + df = df.rename(columns={df.columns[0]: "gen_row"}) + # Name the GEN_STATUS column (col 8 in MATPOWER gen, index 7 after gen_row prepend) + if "GEN_STATUS" not in df.columns and len(df.columns) > 8: + df = df.rename(columns={df.columns[8]: "GEN_STATUS"}) + df = df.set_index("gen_row", drop=False) + return df + + +def normalize_branch_df(branch_df: pd.DataFrame) -> pd.DataFrame: + """Return branch table with numeric endpoint columns and a stable row index.""" + + df = branch_df.copy().reset_index(drop=False) + required_cols = ["F_BUS", "T_BUS"] + missing = [col for col in required_cols if col not in df.columns] + if missing: + raise ValueError( + f"MATPOWER branch table is missing required columns: {', '.join(missing)}" + ) + + df["F_BUS"] = _to_int_series(df["F_BUS"], "F_BUS") + df["T_BUS"] = _to_int_series(df["T_BUS"], "T_BUS") + df = df.rename(columns={df.columns[0]: "branch_row"}) + df = df.set_index("branch_row", drop=False) + return df + + +def summarize_case(case_data: MatpowerCaseData) -> pd.Series: + """Return a concise summary for quick notebook display.""" + + return pd.Series( + { + "case_name": case_data.case_name, + "source_path": case_data.source_path, + "base_mva": case_data.base_mva, + "version": case_data.version, + "n_buses": len(case_data.bus), + "n_generators": len(case_data.gen), + "n_branches": len(case_data.branch), + "n_unique_gen_buses": int(case_data.gen["GEN_BUS"].nunique()), + "n_unique_branch_from_buses": int(case_data.branch["F_BUS"].nunique()), + "n_unique_branch_to_buses": int(case_data.branch["T_BUS"].nunique()), + } + ) + + +def validate_case_tables(case_data: MatpowerCaseData) -> dict[str, Any]: + """Validate key MATPOWER table relationships for notebook sanity checks.""" + + bus_ids = set(case_data.bus["BUS_I"].tolist()) + gen_bus_ids = set(case_data.gen["GEN_BUS"].tolist()) + from_bus_ids = set(case_data.branch["F_BUS"].tolist()) + to_bus_ids = set(case_data.branch["T_BUS"].tolist()) + + missing_gen_buses = sorted(gen_bus_ids - bus_ids) + missing_from_buses = sorted(from_bus_ids - bus_ids) + missing_to_buses = sorted(to_bus_ids - bus_ids) + + required_checks = { + "bus_has_BUS_I": "BUS_I" in case_data.bus.columns, + "gen_has_GEN_BUS": "GEN_BUS" in case_data.gen.columns, + "branch_has_F_BUS": "F_BUS" in case_data.branch.columns, + "branch_has_T_BUS": "T_BUS" in case_data.branch.columns, + } + + return { + "required_columns": required_checks, + "all_required_columns_present": all(required_checks.values()), + "missing_gen_bus_ids": missing_gen_buses, + "missing_branch_from_bus_ids": missing_from_buses, + "missing_branch_to_bus_ids": missing_to_buses, + "gen_buses_are_valid": len(missing_gen_buses) == 0, + "branch_from_buses_are_valid": len(missing_from_buses) == 0, + "branch_to_buses_are_valid": len(missing_to_buses) == 0, + } + + +def validation_report_df(case_data: MatpowerCaseData) -> pd.DataFrame: + """Render validation output as a two-column DataFrame.""" + + report = validate_case_tables(case_data) + rows = [] + for key, value in report.items(): + rows.append({"check": key, "value": value}) + return pd.DataFrame(rows) + + +def create_bus_geo_df_from_gic(file_path: str | Path) -> pd.DataFrame: + """Parse TAMU .gic and return bus-indexed geographic table with lat/lon.""" + + gic_path = Path(file_path).expanduser().resolve() + if not gic_path.exists(): + raise FileNotFoundError(f"GIC file not found: {gic_path}") + + with open(gic_path, "r") as f: + lines = f.readlines() + + data_lines = [] + for line in lines[1:]: + if line.strip().startswith("0"): + break + data_lines.append(line) + + name_str = "\n".join([line.split("'")[1] for line in data_lines]) + data_str = "".join( + [line.split("'")[0] + line.split("'")[-1] for line in data_lines] + ) + + substation_df = pd.DataFrame( + {"index": range(1, len(data_lines) + 1), "name": name_str.split("\n")} + ) + substation_df.set_index("index", inplace=True) + + df_data = pd.read_csv( + StringIO(data_str), + sep=r"\s+", + header=None, + names=["unknown", "lat", "lon", "value"], + ) + substation_df = pd.concat([substation_df, df_data], axis=1) + + start_idx = None + for i, line in enumerate(lines): + if "Begin Bus Substation Data" in line: + start_idx = i + 1 + break + + if start_idx is None: + raise ValueError( + "Could not find 'Begin Bus Substation Data' section in GIC file" + ) + + data_lines = [] + for line in lines[start_idx:]: + if line.strip().startswith("0"): + break + data_lines.append(line) + + bus_substation_map = {} + for line in data_lines: + parts = line.split() + if len(parts) >= 2: + bus_id = int(parts[0]) + substation_idx = int(parts[1]) + bus_substation_map[bus_id] = substation_idx + + bus_geo_df = pd.DataFrame( + substation_df.loc[[v for v in bus_substation_map.values()], :] + ) + bus_geo_df["sub_id"] = bus_geo_df.index + bus_geo_df.index = list(bus_substation_map.keys()) + bus_geo_df.index.name = "BUS_I" + return bus_geo_df + + +def create_bus_geo_df_from_aux(file_path: str | Path) -> pd.DataFrame: + """Parse PowerWorld AUX and return bus-indexed geographic table with lat/lon. + + Expected mapping path: + - `Substation` table provides substation `Number`, `Latitude`, `Longitude` + - `Bus` table provides bus `Number` and `SubNumber` + """ + + aux_path = Path(file_path).expanduser().resolve() + if not aux_path.exists(): + raise FileNotFoundError(f"AUX file not found: {aux_path}") + + lines = aux_path.read_text().splitlines() + substation_df = _parse_aux_table(lines, "Substation") + bus_df = _parse_aux_table(lines, "Bus") + + required_sub_cols = {"Number", "Latitude", "Longitude"} + required_bus_cols = {"Number", "SubNumber"} + missing_sub = sorted(required_sub_cols - set(substation_df.columns)) + missing_bus = sorted(required_bus_cols - set(bus_df.columns)) + if missing_sub: + raise ValueError( + f"AUX Substation table is missing required columns: {', '.join(missing_sub)}" + ) + if missing_bus: + raise ValueError( + f"AUX Bus table is missing required columns: {', '.join(missing_bus)}" + ) + + sub = substation_df[["Number", "Name", "Latitude", "Longitude"]].copy() + sub["Number"] = pd.to_numeric(sub["Number"], errors="coerce") + sub["Latitude"] = pd.to_numeric(sub["Latitude"], errors="coerce") + sub["Longitude"] = pd.to_numeric(sub["Longitude"], errors="coerce") + sub = sub.dropna(subset=["Number", "Latitude", "Longitude"]) + sub["sub_id"] = sub["Number"].astype(int) + + bus = bus_df[["Number", "SubNumber"]].copy() + bus["Number"] = pd.to_numeric(bus["Number"], errors="coerce") + bus["SubNumber"] = pd.to_numeric(bus["SubNumber"], errors="coerce") + bus = bus.dropna(subset=["Number", "SubNumber"]) + bus["BUS_I"] = bus["Number"].astype(int) + bus["sub_id"] = bus["SubNumber"].astype(int) + + merged = bus.merge( + sub[["sub_id", "Name", "Latitude", "Longitude"]], + on="sub_id", + how="left", + ) + if merged[["Latitude", "Longitude"]].isna().any().any(): + missing_bus_ids = ( + merged[merged[["Latitude", "Longitude"]].isna().any(axis=1)]["BUS_I"] + .astype(int) + .tolist() + ) + raise ValueError( + "AUX bus->substation merge has missing coordinates for BUS_I " + f"(first 20): {missing_bus_ids[:20]}" + ) + + bus_geo_df = merged.rename( + columns={"Name": "name", "Latitude": "lat", "Longitude": "lon"} + )[["BUS_I", "name", "sub_id", "lat", "lon"]] + bus_geo_df = ( + bus_geo_df.drop_duplicates(subset=["BUS_I"]).set_index("BUS_I").sort_index() + ) + return bus_geo_df + + +def create_bus_geo_df(file_path: str | Path) -> pd.DataFrame: + """Parse supported geo source (.gic or .aux) into bus-indexed lat/lon.""" + + geo_path = Path(file_path).expanduser().resolve() + suffix = geo_path.suffix.lower() + + if suffix == ".gic": + return create_bus_geo_df_from_gic(geo_path) + if suffix == ".aux": + return create_bus_geo_df_from_aux(geo_path) + + # Fallback for uncommon extensions: try GIC first, then AUX parser. + try: + return create_bus_geo_df_from_gic(geo_path) + except Exception: + return create_bus_geo_df_from_aux(geo_path) + + +def split_points(bus_df: pd.DataFrame, radius: float = 0.006) -> pd.DataFrame: + """Return bus_df with colocated lat/lon points spread on small circles.""" + + out_df = bus_df.copy() + lat_name, lon_name = _infer_lat_lon_cols(out_df) + + unique_locs = ( + out_df[[lat_name, lon_name]] + .groupby([lat_name, lon_name]) + .size() + .reset_index(name="n") + ) + + for _, loc_row in unique_locs.iterrows(): + if int(loc_row["n"]) <= 1: + continue + + lat = loc_row[lat_name] + lon = loc_row[lon_name] + colocated_idx = out_df[ + (out_df[lat_name] == lat) & (out_df[lon_name] == lon) + ].index + n_buses = len(colocated_idx) + theta = 2 * np.pi / n_buses + + for k, bus_idx in enumerate(colocated_idx): + x_new = lon + radius * np.cos(k * theta) + y_new = lat + radius * np.sin(k * theta) + out_df.loc[bus_idx, lat_name] = y_new + out_df.loc[bus_idx, lon_name] = x_new + + return out_df + + +def attach_geo_to_case( + case_data: MatpowerCaseData, + geo_file_path: str | Path, + split_colocated: bool = True, + split_radius: float = 0.006, + gen_fanout: bool = False, + gen_fanout_radius: float = 0.003, + gen_fanout_singles: bool = False, +) -> GeoPrepResult: + """Attach geo lat/lon to bus/gen/branch tables, optionally applying split-points.""" + + bus_geo_df = create_bus_geo_df(geo_file_path) + + bus_df = case_data.bus.copy() + if "BUS_I" not in bus_df.columns: + raise ValueError("case_data.bus must include BUS_I") + + bus_df["lat"] = bus_df["BUS_I"].map(bus_geo_df["lat"]) + bus_df["lon"] = bus_df["BUS_I"].map(bus_geo_df["lon"]) + + if bus_df[["lat", "lon"]].isna().any().any(): + missing_ids = bus_df[bus_df[["lat", "lon"]].isna().any(axis=1)][ + "BUS_I" + ].tolist() + raise ValueError( + f"Missing GIC coordinates for BUS_I ids (first 20): {missing_ids[:20]}" + ) + + n_before = int(bus_df[["lat", "lon"]].drop_duplicates().shape[0]) + if split_colocated: + bus_df = split_points(bus_df, radius=split_radius) + n_after = int(bus_df[["lat", "lon"]].drop_duplicates().shape[0]) + + gen_df = case_data.gen.copy() + bus_indexed = bus_df.set_index("BUS_I") + gen_df["bus_lat"] = gen_df["GEN_BUS"].map(bus_indexed["lat"]) + gen_df["bus_lon"] = gen_df["GEN_BUS"].map(bus_indexed["lon"]) + gen_df["lat"] = gen_df["bus_lat"] + gen_df["lon"] = gen_df["bus_lon"] + if gen_fanout: + gen_df = split_generator_points( + gen_df, radius=gen_fanout_radius, fan_singles=gen_fanout_singles + ) + + branch_df = case_data.branch.copy() + bus_lat = bus_df.set_index("BUS_I")["lat"] + bus_lon = bus_df.set_index("BUS_I")["lon"] + branch_df["from_lat"] = branch_df["F_BUS"].map(bus_lat) + branch_df["from_lon"] = branch_df["F_BUS"].map(bus_lon) + branch_df["to_lat"] = branch_df["T_BUS"].map(bus_lat) + branch_df["to_lon"] = branch_df["T_BUS"].map(bus_lon) + + if all(col in branch_df.columns for col in ["PF", "QF", "PT", "QT", "RATE_A"]): + branch_df["loading_pct"] = compute_branch_loading(branch_df) + + geo_case = MatpowerCaseData( + case_name=case_data.case_name, + source_path=case_data.source_path, + bus=bus_df, + gen=gen_df, + branch=branch_df, + genfuel=case_data.genfuel.copy() if case_data.genfuel is not None else None, + gencost=case_data.gencost.copy() if case_data.gencost is not None else None, + base_mva=case_data.base_mva, + version=case_data.version, + ) + + return GeoPrepResult( + case_data=geo_case, + bus_geo_df=bus_geo_df, + split_applied=split_colocated, + n_unique_bus_locations_before=n_before, + n_unique_bus_locations_after=n_after, + ) + + +def compute_branch_loading(branch_df: pd.DataFrame) -> pd.Series: + """Compute branch loading (%) from PF/QF/PT/QT and RATE_A columns.""" + + sf = np.sqrt(branch_df["PF"] ** 2 + branch_df["QF"] ** 2) + st = np.sqrt(branch_df["PT"] ** 2 + branch_df["QT"] ** 2) + smax = np.maximum(sf, st) + loading_pct = 100 * smax / branch_df["RATE_A"] + loading_pct = loading_pct.where(branch_df["RATE_A"] != 0, np.nan) + return loading_pct + + +def lookup_fault_bus( + bus_df: pd.DataFrame, + fault_bus: int | str, +) -> pd.DataFrame: + """Return the row(s) of *bus_df* matching *fault_bus*. + + - ``int``: exact match on ``BUS_I``. + - ``str``: exact match on ``bus_name``, then case-insensitive substring fallback. + + Returns an empty DataFrame when no match is found. + """ + if isinstance(fault_bus, int): + return bus_df[bus_df["BUS_I"] == fault_bus] + + if "bus_name" not in bus_df.columns: + print( + "WARNING: bus_name column not available — cannot look up fault bus by name." + ) + return pd.DataFrame() + + row = bus_df[bus_df["bus_name"] == fault_bus] + if row.empty: + row = bus_df[bus_df["bus_name"].str.contains(fault_bus, case=False, na=False)] + return row + + +def plot_grid( + case_data: MatpowerCaseData, + map_style: str = "carto-positron", + zoom: float | None = None, + branch_bins: int = 10, + n_hover_points: int = 10, + show_loading: bool = True, + show_gen_connectors: bool = False, + marker_scale: float = 1.0, +) -> go.Figure: + """Geo plot: branches, buses (PD), generators (fuel-colored when available). + + Parameters + ---------- + show_loading: + True — color branches by loading_pct bins (viridis). Requires ``loading_pct`` + on case_data.branch (set by attach_geo_to_case when PF/QF/PT/QT/RATE_A + are present). Falls back silently to gray when the column is missing. + False — uniform gray branches regardless of whether loading_pct exists. + """ + + required_bus_cols = ["lat", "lon", "PD", "BUS_I"] + missing_bus_cols = [c for c in required_bus_cols if c not in case_data.bus.columns] + if missing_bus_cols: + raise ValueError(f"Missing required bus columns: {missing_bus_cols}") + + required_branch_cols = [ + "from_lat", + "from_lon", + "to_lat", + "to_lon", + "F_BUS", + "T_BUS", + ] + missing_branch_cols = [ + c for c in required_branch_cols if c not in case_data.branch.columns + ] + if missing_branch_cols: + raise ValueError(f"Missing required branch columns: {missing_branch_cols}") + + required_gen_cols = ["lat", "lon", "GEN_BUS"] + missing_gen_cols = [c for c in required_gen_cols if c not in case_data.gen.columns] + if missing_gen_cols: + raise ValueError(f"Missing required gen columns: {missing_gen_cols}") + + fig = go.Figure() + + has_loading = show_loading and "loading_pct" in case_data.branch.columns + + if has_loading: + loading_pct = case_data.branch["loading_pct"] + colorscale = _viridis_colorscale() + + for i in range(branch_bins): + lo = i * 100.0 / branch_bins + hi = (i + 1) * 100.0 / branch_bins + idx = loading_pct[(loading_pct >= lo) & (loading_pct < hi)].index + if len(idx) == 0: + continue + + bin_df = case_data.branch.loc[idx] + lat_vals = _branch_lat_array(bin_df, n_hover_points) + lon_vals = _branch_lon_array(bin_df, n_hover_points) + hover_vals = _branch_hover_array(bin_df, n_hover_points) + frac = (i + 0.5) / branch_bins + color = _sample_viridis(colorscale, frac) + + fig.add_trace( + go.Scattermap( + mode="lines", + lat=lat_vals, + lon=lon_vals, + hovertext=hover_vals, + hoverinfo="text", + line=dict(color=color, width=2.5), + showlegend=False, + ) + ) + + # branches with nan loading (RATE_A == 0 or missing) + nan_idx = loading_pct[loading_pct.isna()].index + if len(nan_idx) > 0: + nan_df = case_data.branch.loc[nan_idx] + lat_vals = _branch_lat_array(nan_df, n_hover_points) + lon_vals = _branch_lon_array(nan_df, n_hover_points) + hover_vals = _branch_hover_array(nan_df, n_hover_points) + fig.add_trace( + go.Scattermap( + mode="lines", + lat=lat_vals, + lon=lon_vals, + hovertext=hover_vals, + hoverinfo="text", + line=dict(color="rgba(120, 120, 120, 0.4)", width=1.5), + name="branches (no rating)", + showlegend=True, + ) + ) + + # dummy colorbar trace for loading scale + tickvals = [i / branch_bins for i in range(branch_bins + 1)] + ticktext = [f"{int(v * 100)}%" for v in tickvals] + dcolorsc = _discrete_colorscale( + tickvals, + [ + _sample_viridis(colorscale, (i + 0.5) / branch_bins) + for i in range(branch_bins) + ], + ) + fig.add_trace( + go.Scattermap( + lat=[None], + lon=[None], + mode="markers", + marker=dict( + colorscale=dcolorsc, + showscale=True, + cmin=0, + cmax=1, + colorbar=dict( + title=dict(text="Line
loading"), + thickness=15, + tickvals=tickvals, + ticktext=ticktext, + outlinewidth=0, + x=1.17, + ), + ), + hoverinfo="none", + showlegend=False, + ) + ) + else: + # fallback: no loading data — uniform gray branches + lat_vals = _branch_lat_array(case_data.branch, n_hover_points) + lon_vals = _branch_lon_array(case_data.branch, n_hover_points) + hover_vals = _branch_hover_array(case_data.branch, n_hover_points) + fig.add_trace( + go.Scattermap( + mode="lines", + lat=lat_vals, + lon=lon_vals, + hovertext=hover_vals, + hoverinfo="text", + line=dict(color="rgba(80, 80, 80, 0.45)", width=1.5), + name="branches", + showlegend=True, + ) + ) + + # buses + bus_sizes = ( + (np.sqrt(case_data.bus["PD"].clip(lower=0.01).values) * 2.2 + 3) * marker_scale + ).tolist() + has_bus_name = "bus_name" in case_data.bus.columns + has_json_bus = "json_bus_number" in case_data.bus.columns + bus_hover = [ + ( + f"BUS_I: {row.BUS_I}
" + + (f"{row.bus_name}
" if has_bus_name else "") + + f"PD: {row.PD:.2f} MW
QD: {row.QD:.2f} MVAr" + + (f"
json: number={row.json_bus_number}" if has_json_bus else "") + ) + for row in case_data.bus.itertuples() + ] + fig.add_trace( + go.Scattermap( + name="buses (PD)", + mode="markers", + lat=case_data.bus["lat"].tolist(), + lon=case_data.bus["lon"].tolist(), + hovertext=bus_hover, + hoverinfo="text", + marker=dict( + color=case_data.bus["PD"].values, + colorscale="Blues", + colorbar=dict( + title=dict(text="Load
(PD in MW)"), + x=1.02, + thickness=15, + outlinewidth=0, + ), + size=bus_sizes, + opacity=0.85, + ), + showlegend=True, + ) + ) + + # optional generator fan-out connectors + if show_gen_connectors: + gen_conn_lat, gen_conn_lon, gen_conn_hover = _build_gen_connector_line_arrays( + case_data.gen + ) + if gen_conn_lat: + fig.add_trace( + go.Scattermap( + mode="lines", + lat=gen_conn_lat, + lon=gen_conn_lon, + hoverinfo="none", + line=dict(color="rgba(60, 60, 60, 0.35)", width=1.0), + name="fake gen connections", + showlegend=True, + ) + ) + + # generators + palette = [ + "#1f77b4", + "#d62728", + "#2ca02c", + "#ff7f0e", + "#9467bd", + "#8c564b", + "#e377c2", + "#7f7f7f", + "#bcbd22", + "#17becf", + ] + fuel_color_map = { + "coal": "#222222", + "nuclear": "#a50026", + "hydro": "#4393c3", + "wind": "#4dac26", + "solar": "#f4a582", + "ng": "#fd8d3c", + "gas": "#fd8d3c", + "oil": "#8c6d31", + "geothermal": "#762a83", + "biomass": "#74c476", + "other": "#969696", + "sync_cond": "#d4b9da", + } + if case_data.genfuel is not None and len(case_data.genfuel) == len(case_data.gen): + gen_plot_df = case_data.gen.copy() + gen_plot_df["genfuel"] = case_data.genfuel.values + palette_fuels = [ + f + for f in sorted(gen_plot_df["genfuel"].unique()) + if f not in fuel_color_map + ] + for i, fuel in enumerate(sorted(gen_plot_df["genfuel"].unique())): + gdf = gen_plot_df[gen_plot_df["genfuel"] == fuel] + color = ( + fuel_color_map[fuel] + if fuel in fuel_color_map + else palette[palette_fuels.index(fuel) % len(palette)] + ) + gen_sizes = ( + (np.sqrt(gdf["PG"].clip(lower=0.01).values) * 1.8 + 4) * marker_scale + if "PG" in gdf.columns + else np.full(len(gdf), 6.0 * marker_scale) + ) + _has_json_gen = "json_gen_id" in gdf.columns + _has_json_gen_bus = "json_gen_bus_number" in gdf.columns + gen_hover = [ + f"GEN_BUS: {row.GEN_BUS}
" + f"gen_row: {getattr(row, 'gen_row', 'n/a')}
" + f"PG: {row.PG:.2f} MW
QG: {row.QG:.2f} MVAr
Fuel: {fuel}" + + ( + ( + f"
json: number={getattr(row, 'json_gen_bus_number', row.GEN_BUS)}, device_id={row.json_gen_id}" + if row.json_gen_id is not None + else "
\u26a0 offline \u2014 not in JSON (GEN_STATUS=0)" + ) + if _has_json_gen + else "" + ) + for row in gdf.itertuples() + ] + fig.add_trace( + go.Scattermap( + name=f"gen: {fuel}", + mode="markers", + lat=gdf["lat"].tolist(), + lon=gdf["lon"].tolist(), + hovertext=gen_hover, + hoverinfo="text", + marker=dict(size=gen_sizes.tolist(), color=color, opacity=0.9), + showlegend=True, + ) + ) + else: + fig.add_trace( + go.Scattermap( + name="generators", + mode="markers", + lat=case_data.gen["lat"].tolist(), + lon=case_data.gen["lon"].tolist(), + hovertext=[ + ( + ( + f"GEN_BUS: {row.GEN_BUS}
gen_row: {row.gen_row}" + if "gen_row" in case_data.gen.columns + else f"GEN_BUS: {row.GEN_BUS}" + ) + + ( + ( + f"
json: number={getattr(row, 'json_gen_bus_number', row.GEN_BUS)}, device_id={row.json_gen_id}" + if row.json_gen_id is not None + else "
\u26a0 offline \u2014 not in JSON (GEN_STATUS=0)" + ) + if "json_gen_id" in case_data.gen.columns + else "" + ) + ) + for row in case_data.gen.itertuples() + ], + hoverinfo="text", + marker=dict(size=6 * marker_scale, color="#d62728", opacity=0.9), + showlegend=True, + ) + ) + + center_lon = float(case_data.bus["lon"].mean()) + center_lat = float(case_data.bus["lat"].mean()) + if zoom is None: + zoom = 7 if len(case_data.bus) <= 400 else 5 + + fig.update_layout( + title=f"{case_data.case_name}: grid map" + + (" (loading)" if has_loading else ""), + map={ + "style": map_style, + "center": {"lon": center_lon, "lat": center_lat}, + "zoom": zoom, + }, + legend={"orientation": "h", "y": -0.15, "x": 0, "itemsizing": "constant"}, + margin=dict(r=120), + height=650, + width=780, + ) + return fig + + +def plot_first_geo( + case_data: MatpowerCaseData, + map_style: str = "carto-positron", + zoom: float | None = None, + show_gen_connectors: bool = False, +) -> go.Figure: + """Alias for plot_grid(show_loading=False). Kept for backward compatibility.""" + + return plot_grid( + case_data, + map_style=map_style, + zoom=zoom, + show_loading=False, + show_gen_connectors=show_gen_connectors, + ) + + +def _infer_lat_lon_cols(df: pd.DataFrame) -> tuple[str, str]: + if "lat" in df.columns and "lon" in df.columns: + return "lat", "lon" + if "latitude" in df.columns and "longitude" in df.columns: + return "latitude", "longitude" + if "Latitude" in df.columns and "Longitude" in df.columns: + return "Latitude", "Longitude" + raise ValueError("Could not find latitude/longitude columns in dataframe") + + +def _parse_aux_table(lines: list[str], table_name: str) -> pd.DataFrame: + """Parse a PowerWorld AUX table into a DataFrame. + + Supports both canonical table headers (`Substation (...)`) and a lightweight + `DATA (Substation)` marker style when present. + """ + + table_start_idx = _find_aux_table_header(lines, table_name) + header_line = lines[table_start_idx] + + header_buffer = [header_line] + idx = table_start_idx + 1 + while idx < len(lines) and ")" not in " ".join(header_buffer): + header_buffer.append(lines[idx]) + idx += 1 + + header_text = " ".join(header_buffer) + col_match = re.search(r"\((.*)\)", header_text) + if col_match is None: + raise ValueError(f"Could not parse column header for AUX table '{table_name}'") + + columns = [c.strip() for c in col_match.group(1).split(",") if c.strip()] + if not columns: + raise ValueError(f"AUX table '{table_name}' has no parsed columns") + + while idx < len(lines) and "{" not in lines[idx]: + idx += 1 + if idx >= len(lines): + raise ValueError(f"AUX table '{table_name}' is missing opening '{{'") + idx += 1 + + rows: list[list[str]] = [] + while idx < len(lines): + line = lines[idx].strip() + if not line or line.startswith("//"): + idx += 1 + continue + if line.startswith("}"): + break + + tokens = shlex.split(line) + if len(tokens) < len(columns): + tokens.extend([""] * (len(columns) - len(tokens))) + elif len(tokens) > len(columns): + tokens = tokens[: len(columns)] + rows.append(tokens) + idx += 1 + + if not rows: + raise ValueError(f"AUX table '{table_name}' has no data rows") + + return pd.DataFrame(rows, columns=columns) + + +def _find_aux_table_header(lines: list[str], table_name: str) -> int: + canonical_regex = re.compile(rf"^\s*{re.escape(table_name)}\s*\(") + data_regex = re.compile( + rf"^\s*DATA\s*\(\s*{re.escape(table_name)}\s*\)", re.IGNORECASE + ) + + for idx, line in enumerate(lines): + if canonical_regex.search(line) or data_regex.search(line): + return idx + raise ValueError( + f"Could not find AUX table '{table_name}' (expected '{table_name} (...)' or 'DATA ({table_name})')" + ) + + +def split_generator_points( + gen_df: pd.DataFrame, + radius: float = 0.003, + fan_singles: bool = False, +) -> pd.DataFrame: + """Spread generators on the same bus around the bus center for visibility.""" + + out_df = gen_df.copy() + required_cols = ["GEN_BUS", "lat", "lon"] + missing_cols = [c for c in required_cols if c not in out_df.columns] + if missing_cols: + raise ValueError( + f"Generator fan-out requires columns: {', '.join(required_cols)}" + ) + + for _, grp in out_df.groupby("GEN_BUS", sort=False): + row_idx = grp.index.tolist() + n_gen = len(row_idx) + + center_lat = float(out_df.loc[row_idx[0], "lat"]) + center_lon = float(out_df.loc[row_idx[0], "lon"]) + # single generator: optionally offset east so bus marker underneath is visible + if n_gen == 1: + if fan_singles: + out_df.loc[row_idx[0], "lon"] = center_lon + radius + continue + + theta = 2 * np.pi / n_gen + for k, idx in enumerate(row_idx): + out_df.loc[idx, "lon"] = center_lon + radius * np.cos(k * theta) + out_df.loc[idx, "lat"] = center_lat + radius * np.sin(k * theta) + + return out_df + + +def _viridis_colorscale() -> list[tuple[float, str]]: + """Return viridis colorscale as [(frac, 'rgb(...)'), ...] list.""" + + import plotly.colors as pc + + colors = pc.sequential.Viridis # list of hex/rgb strings + n = len(colors) + return [(i / (n - 1), c) for i, c in enumerate(colors)] + + +def _sample_viridis(colorscale: list[list], frac: float) -> str: + """Sample a color from a plotly colorscale at position frac in [0, 1].""" + + if frac <= 0: + return colorscale[0][1] + if frac >= 1: + return colorscale[-1][1] + for i in range(len(colorscale) - 1): + lo_frac, lo_col = colorscale[i] + hi_frac, hi_col = colorscale[i + 1] + if lo_frac <= frac <= hi_frac: + # just snap to nearest stop for simplicity + t = (frac - lo_frac) / (hi_frac - lo_frac) + return lo_col if t < 0.5 else hi_col + return colorscale[-1][1] + + +def _discrete_colorscale(bvals: list[float], colors: list[str]) -> list[list]: + """Build a plotly discrete colorscale from bin boundaries and per-bin colors.""" + + n = len(colors) + nvals = [(v - bvals[0]) / (bvals[-1] - bvals[0]) for v in bvals] + dcolorscale: list[list] = [] + for k in range(n): + dcolorscale.extend([[nvals[k], colors[k]], [nvals[k + 1], colors[k]]]) + return dcolorscale + + +def _branch_lat_array(branch_df: pd.DataFrame, n_points: int) -> list: + """Interpolate n_points between from_lat and to_lat for each branch, NaN-separated. + + Parallel branches (same F_BUS/T_BUS) are deduplicated — only the first occurrence + is rendered, matching the grouping in _branch_hover_array. + """ + deduped = branch_df.drop_duplicates(subset=["F_BUS", "T_BUS"], keep="first") + arr = np.pad( + np.linspace( + deduped["from_lat"].values, deduped["to_lat"].values, n_points, axis=1 + ), + [(0, 0), (0, 1)], + constant_values=np.nan, + ) + return arr.reshape(-1).tolist() + + +def _branch_lon_array(branch_df: pd.DataFrame, n_points: int) -> list: + """Interpolate n_points between from_lon and to_lon for each branch, NaN-separated. + + Parallel branches (same F_BUS/T_BUS) are deduplicated — only the first occurrence + is rendered, matching the grouping in _branch_hover_array. + """ + deduped = branch_df.drop_duplicates(subset=["F_BUS", "T_BUS"], keep="first") + arr = np.pad( + np.linspace( + deduped["from_lon"].values, deduped["to_lon"].values, n_points, axis=1 + ), + [(0, 0), (0, 1)], + constant_values=np.nan, + ) + return arr.reshape(-1).tolist() + + +def _branch_hover_detail( + row: object, has_loading: bool, has_branch_row: bool, has_json_branch: bool = False +) -> str: + """Build the per-branch detail line (no F_BUS/T_BUS header — used inside groups).""" + parts = [] + if has_branch_row: + parts.append(f"branch_row: {row.branch_row}") + if has_json_branch: + parts.append(f"json: device_id={row.json_branch_id}") + if has_loading and not pd.isna(row.loading_pct): + parts.append(f"Loading: {row.loading_pct:.1f}%") + pf = getattr(row, "PF", None) + if pf is not None: + parts.append(f"PF: {pf:.1f} MW") + rate_a = getattr(row, "RATE_A", None) + if rate_a is not None: + parts.append(f"RATE_A: {rate_a:.1f} MVA") + return "
".join(parts) if parts else "" + + +def _branch_hover_array(branch_df: pd.DataFrame, n_points: int) -> list: + """Build hover text repeated n_points times per branch, with empty string separator. + + Parallel branches (same F_BUS/T_BUS) are grouped under one header line that + states the count, followed by per-branch detail lines separated by ----. + """ + has_loading = "loading_pct" in branch_df.columns + has_branch_row = "branch_row" in branch_df.columns + has_json_branch = "json_branch_id" in branch_df.columns + + groups: dict[tuple, list[str]] = {} + order: list[tuple] = [] + for row in branch_df.itertuples(): + key = (int(row.F_BUS), int(row.T_BUS)) + detail = _branch_hover_detail(row, has_loading, has_branch_row, has_json_branch) + if key not in groups: + groups[key] = [detail] + order.append(key) + else: + groups[key].append(detail) + + result = [] + for key in order: + f_bus, t_bus = key + details = groups[key] + n = len(details) + count = f" ({n} parallel)" if n > 1 else "" + header = f"F_BUS: {f_bus} \u2192 T_BUS: {t_bus}{count}" + body = "
---
".join(d for d in details if d) + hover = f"{header}
{body}" if body else header + result.extend([hover] * n_points + [""]) + return result + + +def _build_branch_line_arrays( + branch_df: pd.DataFrame, +) -> tuple[list[float | None], list[float | None], list[str]]: + lat_vals: list[float | None] = [] + lon_vals: list[float | None] = [] + hover_vals: list[str] = [] + has_loading = "loading_pct" in branch_df.columns + has_branch_row = "branch_row" in branch_df.columns + has_json_branch = "json_branch_id" in branch_df.columns + + for row in branch_df.itertuples(): + lat_vals.extend([row.from_lat, row.to_lat, None]) + lon_vals.extend([row.from_lon, row.to_lon, None]) + + branch_id = f"
branch_row: {row.branch_row}" if has_branch_row else "" + json_bid = ( + f"
json: device_id={row.json_branch_id}" if has_json_branch else "" + ) + if has_loading and not pd.isna(row.loading_pct): + rate_a = getattr(row, "RATE_A", None) + pf = getattr(row, "PF", None) + rate_str = f"
RATE_A: {rate_a:.1f} MVA" if rate_a is not None else "" + pf_str = f"
PF: {pf:.1f} MW" if pf is not None else "" + htxt = ( + f"F_BUS: {row.F_BUS} \u2192 T_BUS: {row.T_BUS}" + f"{branch_id}{json_bid}
" + f"Loading: {row.loading_pct:.1f}%" + f"{pf_str}{rate_str}" + ) + else: + htxt = f"F_BUS: {row.F_BUS} \u2192 T_BUS: {row.T_BUS}{branch_id}{json_bid}" + hover_vals.extend([htxt, htxt, ""]) + + return lat_vals, lon_vals, hover_vals + + +def _build_gen_connector_line_arrays( + gen_df: pd.DataFrame, +) -> tuple[list[float | None], list[float | None], list[str]]: + required_cols = ["bus_lat", "bus_lon", "lat", "lon", "GEN_BUS"] + if any(c not in gen_df.columns for c in required_cols): + return [], [], [] + + lat_vals: list[float | None] = [] + lon_vals: list[float | None] = [] + hover_vals: list[str] = [] + + for row in gen_df.itertuples(): + if ( + pd.isna(row.bus_lat) + or pd.isna(row.bus_lon) + or pd.isna(row.lat) + or pd.isna(row.lon) + ): + continue + + lat_vals.extend([row.bus_lat, row.lat, None]) + lon_vals.extend([row.bus_lon, row.lon, None]) + + gen_row = getattr(row, "gen_row", "n/a") + htxt = f"GEN_BUS: {row.GEN_BUS}
gen_row: {gen_row}" + hover_vals.extend([htxt, htxt, ""]) + + return lat_vals, lon_vals, hover_vals + + +def _parse_m_cell_array(m_text: str, attr_name: str) -> list[str]: + """Extract a MATPOWER curly-brace string cell array from .m file text. + + Handles patterns like:: + + mpc.genfuel = { + 'coal'; + 'wind'; + }; + + Returns a list of strings (quotes and semicolons stripped), or [] if not found. + """ + + pattern = re.compile( + r"mpc\." + re.escape(attr_name) + r"\s*=\s*\{([^}]*)\}\s*;", + re.DOTALL, + ) + m = pattern.search(m_text) + if m is None: + return [] + + entries = [] + for line in m.group(1).splitlines(): + line = line.strip().rstrip(";").strip().strip("'\"") + if line and not line.startswith("%"): + entries.append(line) + return entries + + +def _copy_df_attr(case: CaseFrames, attr_name: str) -> pd.DataFrame | None: + value = getattr(case, attr_name, None) + if value is None: + return None + if not isinstance(value, pd.DataFrame): + return None + return value.copy() + + +def _to_int_series(series: pd.Series, column_name: str) -> pd.Series: + try: + return pd.to_numeric(series, errors="raise").astype(int) + except Exception as exc: + raise ValueError( + f"Column '{column_name}' could not be converted to int" + ) from exc + + +def _coerce_optional_float(value: Any) -> float | None: + if value is None: + return None + try: + return float(value) + except (TypeError, ValueError): + return None + + +def _coerce_optional_str(value: Any) -> str | None: + if value is None: + return None + return str(value)