HydroPumpTurbine with Energy Model

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

Note

HydroPowerSimulations.jl is an extension library of PowerSimulations.jl for modeling hydro units. Users are encouraged to review the tutorial in PowerSimulations.jl on Running a Single-Step Problem before this tutorial.

Load packages

using PowerSystems
using PowerSimulations
using HydroPowerSimulations
using PowerSystemCaseBuilder
using HiGHS ## solver

Data

Note

PowerSystemCaseBuilder.jl is a helper library that makes it easier to reproduce examples in the documentation and tutorials. Normally you would pass your local files to create the system data instead of calling the function PowerSystemCaseBuilder.build_system.

sys = build_system(PSITestSystems, "c_sys5_hydro_pump_energy")
System
Property Value
Name
Description
System Units Base SYSTEM_BASE
Base Power 100.0
Base Frequency 60.0
Num Components 31
Static Components
Type Count
ACBus 5
Arc 6
HydroPumpTurbine 1
HydroReservoir 2
Line 6
PowerLoad 3
RenewableDispatch 3
ThermalStandard 5
Forecast Summary
owner_type owner_category name time_series_type initial_timestamp resolution count horizon interval window_count
String String String String String Dates.CompoundPeriod Int64 Dates.CompoundPeriod Dates.CompoundPeriod Int64
HydroPumpTurbine Component capacity Deterministic 2024-01-01T00:00:00 1 hour 1 1 day 1 day 2
HydroPumpTurbine Component max_active_power Deterministic 2024-01-01T00:00:00 1 hour 1 1 day 1 day 2
PowerLoad Component max_active_power Deterministic 2024-01-01T00:00:00 1 hour 3 1 day 1 day 2
RenewableDispatch Component max_active_power Deterministic 2024-01-01T00:00:00 1 hour 3 1 day 1 day 2

With a single PowerSystems.HydroPumpTurbine connected to two PowerSystems.HydroReservoir (head and tail reservoirs of the turbine):

hy = only(get_components(HydroPumpTurbine, sys))

reservoir_head = get_component(HydroReservoir, sys, "Bat_head_reservoir")
reservoir_tail = get_component(HydroReservoir, sys, "Bat_tail_reservoir")
HydroReservoir: Bat_tail_reservoir:
   name: Bat_tail_reservoir
   available: true
   storage_level_limits: (min = 0.0, max = 0.0)
   initial_level: 0.0
   spillage_limits: nothing
   inflow: 0.0
   outflow: 0.0
   level_targets: nothing
   intake_elevation: 0.0
   head_to_volume_factor: InfrastructureSystems.LinearCurve(0.0, 0.0)
   upstream_turbines: 1-element Vector{PowerSystems.HydroUnit}
   downstream_turbines: 0-element Vector{PowerSystems.HydroUnit}
   upstream_reservoirs: 1-element Vector{PowerSystems.Device}
   operation_cost: 
   level_data_type: PowerSystems.ReservoirDataTypeModule.ReservoirDataType.USABLE_VOLUME = 1
   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: false

Note that the reservoirs has a level_data_type of ENERGY, that implies its storage level limits data are in MWh. That means that the available capacity of the head reservoir is between 5.0 and 400 MWh, while the tail reservoir is set to zero, implying an infinite tail reservoir.

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 HydroPumpEnergyDispatch:

pump_model = DeviceModel(
    HydroPumpTurbine,
    HydroPumpEnergyDispatch;
    attributes = Dict{String, Any}(
        "reservation" => true,
        "energy_target" => false,
    ),
    time_series_names = Dict(),
)
set_device_model!(template, pump_model)

The HydroPumpEnergyDispatch(@ref) is a closed model for turbine and must be connected to independent reservoirs (not connected with other HydroTurbine). For that reason it is not needed to include a model of HydroEnergyModelReservoir. When the attribute reservation is set-up to true it does not allow to simultaneously use the pump and turbine, forcing one of those variables to zero. The energy_target attributes allow to include a final target for the head reservoir based on its level_targets field.

In addition, the time_series_names is set-up to an empty dictionary. By default, the HydroPumpEnergyDispatch(@ref) model allows to include limits on the capacity and max_active_power at each time step if the user need it by properly setting up those time series (similar to a HydroDispatchRunOfRiver model)

time_series_names = Dict(
    ActivePowerTimeSeriesParameter => "max_active_power",
    EnergyCapacityTimeSeriesParameter => "capacity",
)

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

model = DecisionModel(template, sys; optimizer = HiGHS.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

PowerSimulations Problem Expressions Results
ActivePowerBalance__System
ProductionCostExpression__ThermalStandard
PTDFBranchFlow__Line
ActivePowerBalance__ACBus
PowerSimulations Problem Parameters Results
ActivePowerTimeSeriesParameter__PowerLoad
PowerSimulations Problem Variables Results
ActivePowerVariable__ThermalStandard
ActivePowerPumpVariable__HydroPumpTurbine
ReservationVariable__HydroPumpTurbine
ActivePowerVariable__HydroPumpTurbine

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

var = read_variable(
    res,
    "ActivePowerVariable__HydroPumpTurbine";
    table_format = TableFormat.WIDE,
)
14 rows omitted
DateTime Bat_pump
Dates.DateTime Float64?
2024-01-01T00:00:00 200.0
2024-01-01T01:00:00 200.0
2024-01-01T02:00:00 200.0
2024-01-01T03:00:00 200.0
2024-01-01T04:00:00 200.0
2024-01-01T05:00:00 200.0
2024-01-01T06:00:00 200.0
2024-01-01T07:00:00 200.0
2024-01-01T08:00:00 200.0
2024-01-01T09:00:00 200.0

its pump usage

var = read_variable(
    res,
    "ActivePowerPumpVariable__HydroPumpTurbine";
    table_format = TableFormat.WIDE,
)
14 rows omitted
DateTime Bat_pump
Dates.DateTime Float64?
2024-01-01T00:00:00 0.0
2024-01-01T01:00:00 0.0
2024-01-01T02:00:00 0.0
2024-01-01T03:00:00 0.0
2024-01-01T04:00:00 0.0
2024-01-01T05:00:00 0.0
2024-01-01T06:00:00 0.0
2024-01-01T07:00:00 0.0
2024-01-01T08:00:00 0.0
2024-01-01T09:00:00 0.0

and reserve commitment (enabled with reservation = true):

var =
    read_variable(
        res,
        "ReservationVariable__HydroPumpTurbine";
        table_format = TableFormat.WIDE,
    )
14 rows omitted
DateTime Bat_pump
Dates.DateTime Float64?
2024-01-01T00:00:00 1.0
2024-01-01T01:00:00 1.0
2024-01-01T02:00:00 1.0
2024-01-01T03:00:00 1.0
2024-01-01T04:00:00 1.0
2024-01-01T05:00:00 1.0
2024-01-01T06:00:00 1.0
2024-01-01T07:00:00 1.0
2024-01-01T08:00:00 1.0
2024-01-01T09:00:00 1.0

HydroPumpEnergyDispatch does not track stored energy in reservoirs. If you need an MWh trajectory as an optimization state, add a HydroEnergyModelReservoir formulation for the connected reservoirs.