diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index b3154ac..1924612 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -41,6 +41,7 @@ jobs: HOPE/tests PSCAD/tests PLEXOSDB/tests + PSSE/tests surge/test_integration.py -q -rs powerfactory-wheel-catalogs: diff --git a/PSSE/dynamic_stability.py b/PSSE/dynamic_stability.py new file mode 100644 index 0000000..bdf4a48 --- /dev/null +++ b/PSSE/dynamic_stability.py @@ -0,0 +1,124 @@ +"""Pure-Python transient response metrics for PSS/E dynamic output data.""" + +from __future__ import annotations + +import math +from numbers import Real +from typing import Any, Mapping + + +def _finite_real(value: Any, name: str) -> float: + if isinstance(value, bool) or not isinstance(value, Real): + raise ValueError(f"{name} must contain only real numeric values") + result = float(value) + if not math.isfinite(result): + raise ValueError(f"{name} must contain only finite numeric values") + return result + + +def _nonnegative_finite(value: Any, name: str) -> float: + result = _finite_real(value, name) + if result < 0: + raise ValueError(f"{name} must be non-negative") + return result + + +def analyze_channel_series( + time: Any, + channels: Mapping[Any, Any], + *, + tolerance: float = 0.02, + settling_window: float = 0.5, + channel_tolerances: Mapping[str, float] | None = None, +) -> dict[str, dict[str, Any]]: + """Summarize transient channels without returning their full time series. + + tolerance and channel_tolerances are absolute tolerances in each channel's + native units. A channel is considered settled at the earliest sample after + its last excursion outside the final-value band, provided the remaining + simulation duration is at least settling_window seconds. This is a + descriptive signal metric, not a substitute for a study-specific + stability criterion. + """ + default_tolerance = _nonnegative_finite(tolerance, "tolerance") + required_window = _nonnegative_finite(settling_window, "settling_window") + + try: + point_count = len(time) + except (TypeError, AttributeError) as exc: + raise ValueError("time must be a non-empty sized sequence") from exc + if point_count < 2: + raise ValueError("time must contain at least two samples") + + times = [_finite_real(value, "time") for value in time] + if any(right <= left for left, right in zip(times, times[1:])): + raise ValueError("time values must be strictly increasing") + + if not channels: + raise ValueError("at least one dynamic channel is required") + + overrides = dict(channel_tolerances or {}) + channel_ids = {str(key) for key in channels} + unknown = sorted(set(overrides) - channel_ids) + if unknown: + raise ValueError(f"channel_tolerances contains unknown channel(s): {unknown}") + + metrics: dict[str, dict[str, Any]] = {} + for key, raw_values in channels.items(): + channel_id = str(key) + try: + if len(raw_values) != point_count: + raise ValueError( + f"channel {channel_id!r} has {len(raw_values)} samples; " + f"expected {point_count} to match time" + ) + values = [ + _finite_real(value, f"channel {channel_id!r}") + for value in raw_values + ] + except TypeError as exc: + raise ValueError( + f"channel {channel_id!r} must be a sized numeric sequence" + ) from exc + + channel_tolerance = _nonnegative_finite( + overrides.get(channel_id, default_tolerance), + f"tolerance for channel {channel_id!r}", + ) + min_index = min(range(point_count), key=values.__getitem__) + max_index = max(range(point_count), key=values.__getitem__) + final_value = values[-1] + last_outside = max( + ( + index + for index, value in enumerate(values) + if abs(value - final_value) > channel_tolerance + ), + default=-1, + ) + candidate_index = last_outside + 1 + settled = ( + candidate_index < point_count + and times[-1] - times[candidate_index] >= required_window + ) + settling_time = times[candidate_index] if settled else None + + metrics[channel_id] = { + "sample_count": point_count, + "initial_value": values[0], + "final_value": final_value, + "minimum": values[min_index], + "minimum_index": min_index, + "minimum_time": times[min_index], + "maximum": values[max_index], + "maximum_index": max_index, + "maximum_time": times[max_index], + "peak_to_peak": values[max_index] - values[min_index], + "max_abs_deviation_from_final": max( + abs(value - final_value) for value in values + ), + "tolerance": channel_tolerance, + "settled": settled, + "settling_time": settling_time, + } + return metrics diff --git a/PSSE/psse_mcp.py b/PSSE/psse_mcp.py index 1d71b00..ca350ec 100644 --- a/PSSE/psse_mcp.py +++ b/PSSE/psse_mcp.py @@ -2,11 +2,23 @@ import os import json import io +import math from pathlib import Path from mcp.server.mcpserver import MCPServer as FastMCP from typing import Dict, List, Optional, Any +if __package__: + from .dynamic_stability import analyze_channel_series +else: + try: + # Tests and embedded loaders may execute this file under a synthetic + # module name, where the package-relative import is unavailable. + from PSSE.dynamic_stability import analyze_channel_series + except ImportError: + # Preserve direct execution via `python PSSE/psse_mcp.py`. + from dynamic_stability import analyze_channel_series + _repo_root = str(Path(__file__).resolve().parents[1]) _repo_root_added = _repo_root not in sys.path if _repo_root_added: @@ -1033,5 +1045,123 @@ def read_dynamic_output( } +@mcp.tool() +def analyze_dynamic_stability( + outfile: str, + channels: Optional[List[int]] = None, + outvrsn: int = 0, + tolerance: float = 0.02, + settling_window: float = 0.5, + channel_tolerances: Optional[Dict[str, float]] = None, +) -> Dict[str, Any]: + """Compute compact transient-response metrics from a PSS/E output file. + + Tolerances are absolute and expressed in each channel's native units. + The settling time is the first sample after the last excursion outside + the final-value tolerance band, provided the remaining simulation covers + settling_window seconds. This metric is descriptive, not a universal + pass/fail stability criterion. + """ + try: + outfile = checked_path(outfile, purpose="dynamic output") + except PathNotAllowed as exc: + return {"status": "error", "message": str(exc)} + + if not os.path.isfile(outfile): + return { + "status": "error", + "message": f"Dynamic output file does not exist: {os.path.abspath(outfile)}", + } + + for name, value in (("tolerance", tolerance), ("settling_window", settling_window)): + if ( + isinstance(value, bool) + or not isinstance(value, (int, float)) + or not math.isfinite(value) + or value < 0 + ): + return {"status": "error", "message": f"{name} must be a finite non-negative number"} + + if channels is not None: + if not channels: + return {"status": "error", "message": "channels must contain at least one channel number"} + if any( + not isinstance(channel, int) + or isinstance(channel, bool) + or channel < 1 + for channel in channels + ): + return {"status": "error", "message": "channels must contain positive integer channel numbers"} + if len(set(channels)) != len(channels): + return {"status": "error", "message": "channels must not contain duplicates"} + + if channel_tolerances is not None: + if not isinstance(channel_tolerances, dict): + return {"status": "error", "message": "channel_tolerances must be an object keyed by channel number"} + for channel_id, value in channel_tolerances.items(): + if ( + not isinstance(channel_id, str) + or not channel_id.isdecimal() + or isinstance(value, bool) + or not isinstance(value, (int, float)) + or not math.isfinite(value) + or value < 0 + ): + return { + "status": "error", + "message": "channel_tolerances must map numeric channel IDs to finite non-negative numbers", + } + + try: + title, channel_names, data = _dyntools().CHNF( + outfile, outvrsn=outvrsn + ).get_data() + if channels is None: + selected = [ + key for key in channel_names + if key != "time" and key in data + ] + else: + missing = [channel for channel in channels if channel not in channel_names] + if missing: + return {"status": "error", "message": f"Requested channel(s) not found: {missing}"} + selected = list(channels) + + missing_data = [channel for channel in selected if channel not in data] + if missing_data: + return { + "status": "error", + "message": f"Dynamic output is missing data for channel(s): {missing_data}", + } + + metrics = analyze_channel_series( + data.get("time"), + {channel: data[channel] for channel in selected}, + tolerance=tolerance, + settling_window=settling_window, + channel_tolerances=channel_tolerances, + ) + return { + "status": "success", + "outfile": os.path.abspath(outfile), + "title": title, + "sample_count": len(data["time"]), + "time_start": float(data["time"][0]), + "time_end": float(data["time"][-1]), + "settling_window": float(settling_window), + "default_tolerance": float(tolerance), + "channels": { + str(channel): channel_names[channel] + for channel in selected + }, + "metrics": metrics, + } + except Exception as exc: + return { + "status": "error", + "message": f"Failed to analyze dynamic stability: {exc}", + } + + if __name__ == "__main__": mcp.run(transport="stdio") diff --git a/PSSE/tests/test_dynamic_output.py b/PSSE/tests/test_dynamic_output.py index 48129dd..b33fd35 100644 --- a/PSSE/tests/test_dynamic_output.py +++ b/PSSE/tests/test_dynamic_output.py @@ -5,7 +5,8 @@ from PSSE import psse_mcp -def test_dynamic_output_missing_file(): +def test_dynamic_output_missing_file(monkeypatch): + monkeypatch.setattr(psse_mcp, "checked_path", lambda path, purpose: path) result = psse_mcp.list_dynamic_output_channels( os.path.join(os.path.dirname(__file__), "does_not_exist.out") ) @@ -14,7 +15,8 @@ def test_dynamic_output_missing_file(): assert "does not exist" in result["message"] -def test_dynamic_output_read_missing_file(): +def test_dynamic_output_read_missing_file(monkeypatch): + monkeypatch.setattr(psse_mcp, "checked_path", lambda path, purpose: path) result = psse_mcp.read_dynamic_output( os.path.join(os.path.dirname(__file__), "does_not_exist.out") ) diff --git a/PSSE/tests/test_dynamic_stability.py b/PSSE/tests/test_dynamic_stability.py new file mode 100644 index 0000000..c5f54d3 --- /dev/null +++ b/PSSE/tests/test_dynamic_stability.py @@ -0,0 +1,132 @@ +"""Regression tests for PSS/E transient-response analysis.""" + +import math + +import pytest + +from PSSE.dynamic_stability import analyze_channel_series +from PSSE import psse_mcp + + +def test_reports_extrema_and_settling_time_with_irregular_timestamps(): + metrics = analyze_channel_series( + [0.0, 0.1, 0.2, 0.4, 0.6, 1.0], + {"1": [1.0, 1.4, 0.5, 1.02, 1.01, 1.0]}, + tolerance=0.03, + settling_window=0.5, + )["1"] + + assert metrics["minimum"] == 0.5 + assert metrics["minimum_time"] == 0.2 + assert metrics["maximum"] == 1.4 + assert metrics["maximum_time"] == 0.1 + assert metrics["peak_to_peak"] == pytest.approx(0.9) + assert metrics["max_abs_deviation_from_final"] == pytest.approx(0.5) + assert metrics["settled"] is True + assert metrics["settling_time"] == 0.4 + + +def test_late_excursion_prevents_false_settled_result(): + metrics = analyze_channel_series( + [0, 1, 2, 3, 4], + {1: [1.0, 1.01, 1.0, 1.4, 1.0]}, + tolerance=0.02, + settling_window=1.0, + )["1"] + + assert metrics["settled"] is False + assert metrics["settling_time"] is None + + +def test_per_channel_tolerance_overrides_default(): + metrics = analyze_channel_series( + [0, 1, 2, 3], + {1: [1.0, 1.1, 1.02, 1.0], 2: [60.0, 60.1, 60.02, 60.0]}, + tolerance=0.01, + settling_window=1.0, + channel_tolerances={"1": 0.11, "2": 0.03}, + ) + + assert metrics["1"]["tolerance"] == 0.11 + assert metrics["1"]["settled"] is True + assert metrics["2"]["tolerance"] == 0.03 + assert metrics["2"]["settled"] is True + + +@pytest.mark.parametrize( + ("time", "channels", "message"), + [ + ([0, 1, 1], {1: [1, 2, 3]}, "strictly increasing"), + ([0, math.inf], {1: [1, 2]}, "finite numeric"), + ([0, 1], {1: [1]}, "expected 2 to match time"), + ([0, 1], {1: [1, math.nan]}, "finite numeric"), + ([0], {1: [1]}, "at least two samples"), + ([0, 1], {}, "at least one dynamic channel"), + ], +) +def test_rejects_invalid_dynamic_series(time, channels, message): + with pytest.raises(ValueError, match=message): + analyze_channel_series(time, channels) + + +@pytest.mark.parametrize( + ("kwargs", "message"), + [ + ({"tolerance": -1}, "tolerance must be non-negative"), + ({"settling_window": math.inf}, "settling_window must contain only finite"), + ({"channel_tolerances": {"99": 0.1}}, "unknown channel"), + ({"channel_tolerances": {"1": -0.1}}, "tolerance for channel '1' must be non-negative"), + ], +) +def test_rejects_invalid_analysis_options(kwargs, message): + with pytest.raises(ValueError, match=message): + analyze_channel_series([0, 1], {"1": [1.0, 1.0]}, **kwargs) + + +def test_mcp_tool_returns_compact_metrics_without_raw_series(tmp_path, monkeypatch): + outfile = tmp_path / "case.out" + outfile.touch() + + class FakeCHNF: + def __init__(self, path, outvrsn=0): + self.path = path + self.outvrsn = outvrsn + + def get_data(self): + return ( + "test case", + {"time": "Time(s)", 1: "Bus voltage"}, + {"time": [0.0, 0.5, 1.0, 1.5, 2.0], 1: [1.0, 0.8, 1.0, 1.0, 1.0]}, + ) + + class FakeDyntools: + CHNF = FakeCHNF + + monkeypatch.setattr(psse_mcp, "checked_path", lambda path, purpose: path) + monkeypatch.setattr(psse_mcp, "_dyntools", lambda: FakeDyntools) + + result = psse_mcp.analyze_dynamic_stability( + str(outfile), channels=[1], tolerance=0.01, settling_window=0.5 + ) + + assert result["status"] == "success", result + assert result["sample_count"] == 5 + assert result["channels"] == {"1": "Bus voltage"} + assert result["metrics"]["1"]["minimum"] == 0.8 + assert result["metrics"]["1"]["minimum_time"] == 0.5 + assert result["metrics"]["1"]["settled"] is True + assert "data" not in result + + +def test_mcp_tool_rejects_duplicate_channels_before_engine_start(tmp_path, monkeypatch): + outfile = tmp_path / "case.out" + outfile.touch() + calls = [] + monkeypatch.setattr(psse_mcp, "checked_path", lambda path, purpose: path) + monkeypatch.setattr(psse_mcp, "_dyntools", lambda: calls.append("engine")) + + result = psse_mcp.analyze_dynamic_stability(str(outfile), channels=[1, 1]) + + assert result["status"] == "error" + assert "duplicates" in result["message"] + assert calls == []