HydroTurbine + Reservoir for WaterModel

To follow along, you can download this tutorial as a Julia script (.jl) or Jupyter notebook (.ipynb).

Load packages

using PowerSystems
using PowerSimulations
using HydroPowerSimulations
using PowerSystemCaseBuilder
using Ipopt ## solver

Data

sys = build_system(PSITestSystems, "c_sys5_hy_turbine_head")

With a single PowerSystems.HydroTurbine connected downstream to a PowerSystems.HydroReservoir:

hy = only(get_components(HydroTurbine, sys))

reservoir = only(get_components(HydroReservoir, sys))

HydroReservoir: Water_Reservoir:
   name: Water_Reservoir
   available: true
   storage_level_limits: (min = 463.5, max = 555.5)
   initial_level: 0.9
   spillage_limits: nothing
   inflow: 1.0
   outflow: 0.0
   level_targets: 1.0
   intake_elevation: 463.3
   head_to_volume_factor: InfrastructureSystems.LinearCurve(302376.2, 0.0)
   upstream_turbines: 0-element Vector{PowerSystems.HydroUnit}
   downstream_turbines: 1-element Vector{PowerSystems.HydroUnit}
   upstream_reservoirs: 0-element Vector{PowerSystems.Device}
   operation_cost: 
   level_data_type: PowerSystems.ReservoirDataTypeModule.ReservoirDataType.HEAD = 3
   ext: Dict{String, Any}()
   InfrastructureSystems.SystemUnitsSettings:
      base_value: 100.0
      unit_system: InfrastructureSystems.UnitSystemModule.UnitSystem.SYSTEM_BASE = 0
   has_supplemental_attributes: false
   has_time_series: true

Note that the reservoir has a level_data_type of HEAD, that implies its storage level limits data are in meters (above the sea level) and refer to the hydraulic head levels. That means that its available capacity lies with its hydraulic head being within 463.5 and 555.5 meters, and its intake elevation is at 463.3 meters. In addition note that the elevation of the turbine is on 317.12 meters above the sea level.

Decision Model

Setting up the formulations based on PowerSimulations.jl:

template = ProblemTemplate(PTDFPowerModel)
set_device_model!(template, ThermalStandard, ThermalBasicDispatch)
set_device_model!(template, PowerLoad, StaticPowerLoad)
set_device_model!(template, Line, StaticBranch)

but, now we also include the HydroTurbine using HydroTurbineBilinearDispatch:

set_device_model!(template, HydroTurbine, HydroTurbineBilinearDispatch)

This is a nonlinear model that to compute its output power requires the bilinear term head times water flow. For that purpose the non-convex Ipopt solver will be used to solve this problem.

In addition, we need to use the water model for the HydroReservoir via HydroWaterModelReservoir.

set_device_model!(template, HydroReservoir, HydroWaterModelReservoir)

With the template properly set-up, we construct, build and solve the optimization problem:

model = DecisionModel(template, sys; optimizer = Ipopt.Optimizer)
build!(model; output_dir = mktempdir())
solve!(model)

InfrastructureSystems.Simulation.RunStatusModule.RunStatus.SUCCESSFULLY_FINALIZED = 0

Exploring Results

Results can be explored using:

res = OptimizationProblemResults(model)

Start: 2024-01-01T00:00:00

End: 2024-01-01T23:00:00

Resolution: 60 minutes

Use read_variable to read in the dispatch variable results for the hydro:

var =
    read_variable(res, "ActivePowerVariable__HydroTurbine"; table_format = TableFormat.WIDE)

or the water flowing through the turbine (in m³/s):

var = read_expression(
    res,
    "TotalHydroFlowRateTurbineOutgoing__HydroTurbine";
    table_format = TableFormat.WIDE,
)

and the head level of the reservoir:

hydraulic_head = read_variable(
    res,
    "HydroReservoirHeadVariable__HydroReservoir";
    table_format = TableFormat.WIDE,
)

Note that since the water outflow limit of the turbine is limited on 30 m³/s, the optimal solution decides to flow as much water as possible producing power around 190 MW with that flow and hydraulic head.