Commit 4e161ca3 authored by Turner, Sean's avatar Turner, Sean
Browse files

Add time-varying constraints and max operating elevation

parent f6e3d564

.claude/settings.local.json

0 → 100644
+82 −0
Original line number Diff line number Diff line 
{

  "permissions": {

    "allow": [

      "Bash(dir:*)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" log --oneline -20)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" config --list)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" status)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" diff)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" log --oneline -5)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" add CLAUDE.md)",

      "Bash(git -C \"C:\\\\Users\\\\snk\\\\powersheds\" commit -m \"$\\(cat <<''EOF''\nadd CLAUDE.md for Claude Code guidance\n\nProvides build commands, architecture overview, and key concepts\nfor future Claude Code instances working in this repository.\n\nCo-Authored-By: Claude Opus 4.5 <noreply@anthropic.com>\nEOF\n\\)\")",

      "Bash(.venv/Scripts/activate)",

      "Bash(maturin develop:*)",

      "Bash(.venv/Scripts/python:*)",

      "Bash(.venv/Scripts/pip install:*)",

      "Bash(.venv/Scripts/maturin develop:*)",

      "Bash(uv run maturin:*)",

      "Bash(uv run python:*)",

      "Bash(uv run pytest:*)",

      "Bash(uv add:*)",

      "Bash(uv sync)",

      "Bash(uv pip install:*)",

      "Bash(git stash:*)",

      "Bash(git checkout:*)",

      "Bash(git add:*)",

      "Bash(git commit:*)",

      "Bash(cargo clean:*)",

      "Bash(ls:*)",

      "Bash(cargo build:*)",

      "Bash(git -C C:/Users/snk/powersheds status)",

      "Bash(git -C C:/Users/snk/powersheds log --oneline -5)",

      "Bash(git -C C:/Users/snk/powersheds diff src/helpers.rs)",

      "Bash(git -C C:/Users/snk/powersheds diff src/lib.rs)",

      "Bash(git -C C:/Users/snk/powersheds add src/helpers.rs src/lib.rs tests/test_hpf.py)",

      "Bash(git -C C:/Users/snk/powersheds commit -m \"$\\(cat <<''EOF''\nHandle infeasible power targets with symmetric clamping\n\nWhen target power falls outside the feasible range at a given head,\nthe HPF interpolation now clamps to the nearest feasible boundary\n\\(min or max\\) rather than returning NaN. This maintains roundtrip\nconsistency between P→Q and Q→P functions.\n\nChanges:\n- src/helpers.rs: Implement symmetric clamping for infeasible power,\n  use linear scan for P index lookup \\(binary search unreliable with\n  approx_eq\\), remove duplicate actual_power calculation in lib.rs\n- tests/test_hpf.py: Update tolerances for flat regions \\(~2 MW corner,\n  ~1.5 MW edge\\), filter clamping cases from roundtrip tests, increase\n  performance thresholds for linear scan algorithm\n\nFixes NaN values in target_release for Cumberland example simulation.\n\nCo-Authored-By: Claude Opus 4.5 <noreply@anthropic.com>\nEOF\n\\)\")",

      "Bash(git -C C:/Users/snk/powersheds push)",

      "Bash(python:*)",

      "Bash(python3:*)",

      "Bash(find:*)",

      "Bash(.venvScriptsactivate.ps1)",

      "Bash(.\\\\.venv\\\\Scripts\\\\python.exe:*)",

      "Bash(.venv/Scripts/python.exe:*)",

      "Bash(git rm:*)",

      "Bash(git push:*)",

      "Bash(.venv\\\\Scripts\\\\python.exe:*)",

      "Bash(C:Userssnkpowersheds.venvScriptspython.exe examples/plot_network.py examples/Cumberland/cascade_config.yaml)",

      "Bash(cmd /c \"C:\\\\Users\\\\snk\\\\powersheds\\\\.venv\\\\Scripts\\\\activate.bat && maturin develop\")",

      "Bash(cmd /c \"C:\\\\Users\\\\snk\\\\powersheds\\\\.venv\\\\Scripts\\\\activate.bat && maturin develop 2>&1\")",

      "Bash(C:Userssnkpowersheds.venvScriptsmaturin.exe develop)",

      "Bash(powershell -Command:*)",

      "WebSearch",

      "Bash(..venvScriptsactivate.ps1)",

      "Bash(source:*)",

      "Bash(jupyter nbconvert:*)",

      "Bash(head:*)",

      "Bash(. .venv/Scripts/activate)",

      "Bash(done)",

      "Bash(grep:*)",

      "Bash(\"../../../.venv/Scripts/python.exe\" -c:*)",

      "Bash(.venv/Scripts/maturin.exe:*)",

      "Bash(tasklist:*)",

      "Bash(taskkill:*)",

      "Bash(cmd //c \"taskkill /F /PID 48996\")",

      "Bash(cmd //c:*)",

      "Bash(wc:*)",

      "Bash(for f in /c/Users/snk/powersheds/examples/private/Kennebec/storage_HWelevation_tables/*.csv)",

      "Bash(do)",

      "Bash(basename=\"$f##*/\")",

      "Bash(basename=\"$basename%.csv\")",

      "Bash(if [ ! -f \"/c/Users/snk/powersheds/examples/private/Kennebec/outflow_TWelevation_tables/$basename.csv\" ])",

      "Bash(then)",

      "Bash(echo:*)",

      "Bash(fi:*)",

      "Bash(md5sum:*)",

      "Bash(xargs:*)",

      "Bash(\"C:/Users/snk/powersheds/.venv/Scripts/python\" -m jupyter nbconvert --to notebook --execute experiments.ipynb --output executed_experiments.ipynb --ExecutePreprocessor.timeout=600)",

      "Bash(\"C:/Users/snk/powersheds/.venv/Scripts/python\" \"C:/Users/snk/powersheds/examples/private/Kennebec/_test_exp.py\")",

      "Bash(\"C:/Users/snk/powersheds/.venv/Scripts/python\" -c \"import powersheds; print\\(powersheds.__file__\\)\")",

      "Bash(Rscript:*)"

    ]

  }

}

README.md

+70 −39
Original line number Diff line number Diff line
@@ -81,46 +81,74 @@ Dependencies (pandas, numpy, scipy, matplotlib, PyYAML, pyarrow, etc.) are decla 


```python

import powersheds



# Minimal single-reservoir cascade (toy numbers for illustration)

cascade_data = {

    "reservoirs": {},

    "rivers": {},

    "confluences": {}

}

from dataclasses import dataclass, field



@dataclass

class ReservoirData:

    object_type: str

    simulation_order: int

    downstream_object: str

    capacity: float

    initial_pool_elevation: float

    min_power_pool: list[float]

    tailwater_elevation: float

    max_release: list[float]

    min_release: list[float]

    set_storage: list[float]

    set_elevation: list[float]

    catchment_inflow: list[float]

    target_power: list[float]

    hpf_h: list[float]

    hpf_p: list[float]

    hpf_q: list[float]

    set_tw_outflow: list[float] = field(default_factory=list)

    set_tw_elevation: list[float] = field(default_factory=list)

    max_operating_elevation: list[float] = field(default_factory=list)



@dataclass

class CascadeData:

    reservoirs: dict[str, ReservoirData]

    rivers: dict[str, object]

    confluences: dict[str, object]



# Define a reservoir with 24 hours of inflow and power targets

cascade_data["reservoirs"]["DemoReservoir"] = {

    "object_type": "reservoir",

    "simulation_order": 1,

    "downstream_object": "NA",  # terminal object

    "capacity": 500.0,            # Mm3

    "initial_pool_elevation": 200.0,

    "min_power_pool": 190.0,

    "tailwater_elevation": 150.0,

    "max_release": 2.0,           # Mm3/hr

    "min_release": 0.0,           # Mm3/hr

cascade_data = CascadeData(

    reservoirs={

        "DemoReservoir": ReservoirData(

            object_type="reservoir",

            simulation_order=1,

            downstream_object="NA",  # terminal object

            capacity=500.0,            # Mm3

            initial_pool_elevation=200.0,

            min_power_pool=[190.0] * 24,

            tailwater_elevation=150.0,

            max_release=[2.0] * 24,    # Mm3/hr

            min_release=[0.0] * 24,    # Mm3/hr

            # storage–elevation (SET) curve

    "set_storage":   [0.0, 250.0, 500.0],

    "set_elevation": [180.0, 195.0, 210.0],

            set_storage=[0.0, 250.0, 500.0],

            set_elevation=[180.0, 195.0, 210.0],

            # hourly inputs

    "catchment_inflow": [0.2] * 24,

    "target_power":     [50.0] * 24,  # MW

            catchment_inflow=[0.2] * 24,

            target_power=[50.0] * 24,  # MW

            # head–power–flow (HPF) grid: cartesian product of H and P

            # H (m):    [40, 60]

            # P (MW):   [ 0, 50, 100]

            # Q (m3/s): row-major values for each (H,P) pair

    "hpf_h": [40, 40, 40, 60, 60, 60],

    "hpf_p": [ 0, 50,100,  0, 50,100],

    "hpf_q": [ 0, 60,120,  0, 45, 90],

}

            hpf_h=[40, 40, 40, 60, 60, 60],

            hpf_p=[0, 50, 100, 0, 50, 100],

            hpf_q=[0, 60, 120, 0, 45, 90],

        )

    },

    rivers={},

    confluences={},

)



results = powersheds.simulate_cascade(cascade_data)



res = results["DemoReservoir"]

print("release (Mm3/hr)", res.release[:5])

print("storage (Mm3)",    res.storage[:5])

print("power (MW)",       res.actual_power[:5])

print("release (Mm3/hr)", res["release"][:5])

print("storage (Mm3)",    res["storage"][:5])

print("power (MW)",       res["actual_power"][:5])

```



#### Notes
@@ -133,7 +161,7 @@ print("power (MW)", res.actual_power[:5]) 


## 📚 Data Model



Pass a single `cascade_data` mapping with three sections: `reservoirs`, `rivers`, and `confluences`. Each contains a dictionary of named objects.

Pass a single `cascade_data` object with `reservoirs`, `rivers`, and `confluences` attributes. A lightweight Python `dataclass` works well; each section should be a dictionary of named objects with matching attributes.



Reservoir
@@ -141,9 +169,12 @@ Reservoir 
- `simulation_order`: integer; lower runs earlier

- `downstream_object`: name of downstream object or "NA"

- `capacity`: reservoir storage capacity (Mm3)

- `initial_pool_elevation` and `min_power_pool` (m)

- `initial_pool_elevation` (m)

- `min_power_pool`: list[float] length T (m)

- `tailwater_elevation` (m)

- `max_release` and `min_release` (Mm3/hr)

- `max_release`: list[float] length T (Mm3/hr)

- `min_release`: list[float] length T (Mm3/hr)

- `max_operating_elevation`: optional list[float] length T, or empty

- `set_storage` and `set_elevation`: arrays for the storage–elevation curve

- `catchment_inflow`: list[float] length T

- `target_power`: list[float] length T (MW)
@@ -170,12 +201,12 @@ Confluence 
Python entry point



```python

powersheds.simulate_cascade(cascade_data: dict) -> dict[str, Result]

powersheds.simulate_cascade(cascade_data) -> dict[str, Result]

```



Returned per-object results:



- Reservoir: `storage`, `pool_elevation`, `head`, `release`, `spill`, `catchment_inflow`, `total_inflow`, `target_release`, `target_power`, `actual_power`

- Reservoir: `storage`, `pool_elevation`, `head`, `release`, `spill`, `catchment_inflow`, `total_inflow`, `target_release`, `target_power`, `actual_power`, `tailwater_elevation`, `violation_min_power_pool`, `violation_min_release`, `violation_max_release`, `violation_max_operating_elevation`

- River: `inflow`, `outflow` (with travel-time lag)

- Confluence: `inflow`, `outflow` (no lag)

examples/Cumberland/cumberland_diagnostics.html

0 → 100644
+3888 −0

File added.

Preview size limit exceeded, changes collapsed.

examples/Cumberland/notebook.ipynb

+12 −3
Original line number Diff line number Diff line 
%% Cell type:markdown id: tags:



# Powersheds Cumberland Cascade Diagnostics



This notebook loads the 8-reservoir Cumberland River Basin cascade,

runs the Powersheds simulation, and produces:

1. A **static network diagram** showing the cascade topology

2. An **interactive time series diagnostic** (Plotly) with dropdown reservoir selection



%% Cell type:code id: tags:



``` python

from dataclasses import dataclass

from pathlib import Path



import pandas as pd

import yaml



import powersheds

from powersheds.diagnostics import cascade_diagnostics, network_diagram, save_diagnostics_html

```



%% Cell type:markdown id: tags:



## 1. Load cascade configuration and time series



%% Cell type:code id: tags:



``` python

from dataclasses import dataclass, field



@dataclass

class ReservoirData:

    object_type: str

    capacity: float

    initial_pool_elevation: float

    min_power_pool: float

    min_power_pool: list[float]

    set_storage: list[float]

    set_elevation: list[float]

    hpf_h: list[float]

    hpf_p: list[float]

    hpf_q: list[float]

    tailwater_elevation: float

    max_release: float

    min_release: float

    max_release: list[float]

    min_release: list[float]

    catchment_inflow: list[float]

    target_power: list[float]

    simulation_order: int

    downstream_object: str

    set_tw_outflow: list[float] = field(default_factory=list)

    set_tw_elevation: list[float] = field(default_factory=list)

    max_operating_elevation: list[float] = field(default_factory=list)



@dataclass

class RiverData:

    object_type: str

    simulation_order: int

    downstream_object: str

    lag: int

    legacy_flows: list[float]



@dataclass

class ConfluenceData:

    object_type: str

    simulation_order: int

    downstream_object: str



@dataclass

class CascadeData:

    reservoirs: dict[str, ReservoirData]

    rivers: dict[str, RiverData]

    confluences: dict[str, ConfluenceData]

```



%% Cell type:code id: tags:



``` python

CUMBERLAND_DIR = Path(".").resolve()



with open(CUMBERLAND_DIR / "cascade_config.yaml") as f:

    config_dict = yaml.safe_load(f)



reservoir_dict = {}

river_dict = {}

confluence_dict = {}



for name, specs in config_dict.items():

    if specs["object_type"] == "reservoir":

        ts = pd.read_csv(CUMBERLAND_DIR / "time_series" / f"{name}.csv")

        se = pd.read_csv(CUMBERLAND_DIR / "storage_HWelevation_tables" / f"{name}.csv")

        hpf = pd.read_parquet(CUMBERLAND_DIR / "head_power_flow_tables" / name)

        tw_path = CUMBERLAND_DIR / "outflow_TWelevation_tables" / f"{name}.csv"

        if tw_path.exists():

            tw = pd.read_csv(tw_path)

            set_tw_outflow = tw["outflow_cumecs"].tolist()

            set_tw_elevation = tw["elevation_m"].tolist()

        else:

            set_tw_outflow, set_tw_elevation = [], []



        n_hours = len(ts)

        # Convert scalar YAML constraint values to time series

        specs = dict(specs)  # copy to avoid mutating config_dict

        specs["min_power_pool"] = [specs["min_power_pool"]] * n_hours

        specs["max_release"] = [specs["max_release"]] * n_hours

        specs["min_release"] = [specs["min_release"]] * n_hours



        reservoir_dict[name] = ReservoirData(

            **specs,

            catchment_inflow=ts["catchment_inflow"].tolist(),

            target_power=ts["target_power"].tolist(),

            set_storage=se["storage_Mm3"].tolist(),

            set_elevation=se["elevation_m"].tolist(),

            hpf_h=hpf["H_m"].astype(float).tolist(),

            hpf_p=hpf["P_MW"].astype(float).tolist(),

            hpf_q=hpf["Q_cumecs"].astype(float).tolist(),

            set_tw_outflow=set_tw_outflow,

            set_tw_elevation=set_tw_elevation,

        )

    elif specs["object_type"] == "river":

        legacy = pd.read_csv(CUMBERLAND_DIR / "time_series" / f"{name}.csv")

        river_dict[name] = RiverData(

            **specs,

            legacy_flows=legacy["legacy_flow"].tolist(),

        )

    elif specs["object_type"] == "confluence":

        confluence_dict[name] = ConfluenceData(**specs)



cascade_data = CascadeData(

    reservoirs=reservoir_dict,

    rivers=river_dict,

    confluences=confluence_dict,

)



print(f"Loaded {len(reservoir_dict)} reservoirs, {len(river_dict)} rivers, {len(confluence_dict)} confluences")

```



%% Cell type:markdown id: tags:



## 2. Run the simulation



%% Cell type:code id: tags:



``` python

results = powersheds.simulate_cascade(cascade_data)



print(f"Simulation complete: {len(results)} objects")

for name, data in results.items():

    if isinstance(data, dict) and "actual_power" in data:

        n = len(data["actual_power"])

        print(f"  {name}: {n} hours")

```



%% Cell type:markdown id: tags:



## 3. Network diagram



Static visualization of the cascade topology. Trapezoid nodes represent

reservoirs (amber bar = hydropower), circles represent confluences, and

edges show river connections with lag times in hours.



%% Cell type:code id: tags:



``` python

network_diagram(cascade_data, title="Cumberland River Basin Cascade");

```



%% Cell type:markdown id: tags:



## 4. Interactive time series diagnostics



Three linked panels (Flows, Elevation, Power) with shared x-axis zoom/pan.

Use the **dropdown** to switch between reservoirs.



- **Flows panel**: release, target release, inflow, spill, with grey band showing min/max release constraints

- **Elevation panel**: pool elevation, tailwater elevation, and head

- **Power panel**: target vs actual power, with red shading for **infeasible generation**



%% Cell type:code id: tags:



``` python

fig = cascade_diagnostics(

    results,

    cascade_data=cascade_data,

    title="Cumberland River Basin \u2014 Cascade Diagnostics",

    height=900,

    width=1100,

)

fig.show()

```



%% Cell type:markdown id: tags:



## 5. Export to standalone HTML



Save the interactive time series as a self-contained HTML file.



%% Cell type:code id: tags:



``` python

output_path = save_diagnostics_html(

    results,

    output_path="cumberland_diagnostics.html",

    title="Cumberland River Basin \u2014 Cascade Diagnostics",

)

print(f"Saved: {output_path.resolve()}")

```

python/powersheds/diagnostics.py

+132 −27

File changed.

Preview size limit exceeded, changes collapsed.