HydroPumpTurbine with Energy Model

Note

HydroPowerSimulations.jl is an extension library of PowerSimulations.jl for modeling hydro units. Users are encouraged to review the single-step tutorial in PowerSimulations.jl before this tutorial.

Load packages

julia> using PowerSystems
julia> using PowerSimulations
julia> using HydroPowerSimulations
julia> using PowerSystemCaseBuilder
julia> 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 build_system. For more details visit PowerSystemCaseBuilder README

julia> sys = build_system(PSITestSystems, "c_sys5_hydro_pump_energy")┌ Info: Building new system c_sys5_hydro_pump_energy from raw data
  sys_descriptor.raw_data = "/home/runner/.julia/artifacts/edcb5940e84a802a86ad4f2223214d33121ac044/PowerSystemsTestData-4.0.2/psy_data/data_5bus_pu.jl"
[ Info: Serialized time series data to /home/runner/.julia/packages/PowerSystemCaseBuilder/zW01F/data/serialized_system/4e67b70ea6977dbe21c7731d72cdc1494adf072a7f3f08d921db740cf264ce79/c_sys5_hydro_pump_energy_time_series_storage.h5.
[ Info: Serialized System to /home/runner/.julia/packages/PowerSystemCaseBuilder/zW01F/data/serialized_system/4e67b70ea6977dbe21c7731d72cdc1494adf072a7f3f08d921db740cf264ce79/c_sys5_hydro_pump_energy.json
[ Info: Serialized System metadata to /home/runner/.julia/packages/PowerSystemCaseBuilder/zW01F/data/serialized_system/4e67b70ea6977dbe21c7731d72cdc1494adf072a7f3f08d921db740cf264ce79/c_sys5_hydro_pump_energy_metadata.json
              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  i ⋯
│ String             String          String            String            S ⋯
├───────────────────┼────────────────┼──────────────────┼──────────────────┼────
│ HydroPumpTurbine  │ Component      │ capacity         │ Deterministic    │ 2 ⋯
│ HydroPumpTurbine  │ Component      │ max_active_power │ Deterministic    │ 2 ⋯
│ PowerLoad         │ Component      │ max_active_power │ Deterministic    │ 2 ⋯
│ RenewableDispatch │ Component      │ max_active_power │ Deterministic    │ 2 ⋯
└───────────────────┴────────────────┴──────────────────┴──────────────────┴────
                                                               6 columns omitted

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

julia> hy = only(get_components(HydroTurbine, sys))ERROR: ArgumentError: Collection is empty, must contain exactly 1 element
julia> res_head = get_component(HydroReservoir, sys, "Bat_head_reservoir")HydroReservoir: Bat_head_reservoir:
   name: Bat_head_reservoir
   available: true
   storage_level_limits: (min = 5.0, max = 400.0)
   initial_level: 0.5
   spillage_limits: nothing
   inflow: 1.0
   outflow: 0.0
   level_targets: 0.0
   intake_elevation: 0.0
   head_to_volume_factor: InfrastructureSystems.LinearCurve(0.0, 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.ENERGY = 4
   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
julia> res_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:

julia> template = ProblemTemplate(PTDFPowerModel)                     Network Model
┌────────────────────┬─────────────────────────────────┐
│ Network Model      │ PowerSimulations.PTDFPowerModel │
│ Slacks             │ false                           │
│ PTDF               │ false                           │
│ Duals              │ None                            │
│ HVDC Network Model │ None                            │
└────────────────────┴─────────────────────────────────┘

            Device Models
┌─────────────┬─────────────┬────────┐
│ Device Type  Formulation  Slacks │
└─────────────┴─────────────┴────────┘
julia> set_device_model!(template, ThermalStandard, ThermalBasicDispatch)
julia> set_device_model!(template, PowerLoad, StaticPowerLoad)
julia> set_device_model!(template, Line, StaticBranch)

but, now we also include the HydroTurbine using HydroPumpEnergyDispatch:

julia> pump_model = DeviceModel(
           HydroPumpTurbine,
           HydroPumpEnergyDispatch;
           attributes = Dict{String, Any}(
               "reservation" => true,
               "energy_target" => false,
           ),
           time_series_names = Dict(),
       )
                            Device Model
┌───────────────────────────────┬─────────────────────────┬────────┐
│ Device Type                    Formulation              Slacks │
├───────────────────────────────┼─────────────────────────┼────────┤
│ PowerSystems.HydroPumpTurbine │ HydroPumpEnergyDispatch │ false  │
└───────────────────────────────┴─────────────────────────┴────────┘

       Attributes
┌───────────────┬───────┐
│ Name           Value │
├───────────────┼───────┤
│ reservation   │ true  │
│ energy_target │ false │
└───────────────┴───────┘

No FeedForwards Assigned
julia> 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:

julia> model = DecisionModel(template, sys; optimizer = HiGHS.Optimizer)[ Info: Overriding time_series_cache_size because time series is stored in memory
                     Network Model
┌────────────────────┬─────────────────────────────────┐
│ Network Model      │ PowerSimulations.PTDFPowerModel │
│ Slacks             │ false                           │
│ PTDF               │ false                           │
│ Duals              │ None                            │
│ HVDC Network Model │ None                            │
└────────────────────┴─────────────────────────────────┘

                                 Device Models
┌───────────────────────────────┬───────────────────────────────────────┬───────
│ Device Type                    Formulation                            Slac ⋯
├───────────────────────────────┼───────────────────────────────────────┼───────
│ PowerSystems.ThermalStandard  │ PowerSimulations.ThermalBasicDispatch │ fals ⋯
│ PowerSystems.HydroPumpTurbine │ HydroPumpEnergyDispatch               │ fals ⋯
│ PowerSystems.PowerLoad        │ PowerSimulations.StaticPowerLoad      │ fals ⋯
└───────────────────────────────┴───────────────────────────────────────┴───────
                                                                1 column omitted

                        Branch Models
┌───────────────────┬───────────────────────────────┬────────┐
│ Branch Type        Formulation                    Slacks │
├───────────────────┼───────────────────────────────┼────────┤
│ PowerSystems.Line │ PowerSimulations.StaticBranch │ false  │
└───────────────────┴───────────────────────────────┴────────┘
julia> build!(model; output_dir = mktempdir())InfrastructureSystems.Optimization.ModelBuildStatusModule.ModelBuildStatus.BUILT = 0
julia> solve!(model)InfrastructureSystems.Simulation.RunStatusModule.RunStatus.SUCCESSFULLY_FINALIZED = 0

Exploring Results

Results can be explored using:

julia> res = OptimizationProblemResults(model)Start: 2024-01-01T00:00:00
End: 2024-01-01T23:00:00
Resolution: 60 minutes

PowerSimulations Problem Expressions Results
┌───────────────────────────────────────────┐
│ PTDFBranchFlow__Line                      │
│ ActivePowerBalance__ACBus                 │
│ ProductionCostExpression__ThermalStandard │
│ ActivePowerBalance__System                │
└───────────────────────────────────────────┘

 PowerSimulations Problem Parameters Results
┌───────────────────────────────────────────┐
│ ActivePowerTimeSeriesParameter__PowerLoad │
└───────────────────────────────────────────┘

 PowerSimulations Problem Variables Results
┌───────────────────────────────────────────┐
│ ActivePowerVariable__ThermalStandard      │
│ ActivePowerPumpVariable__HydroPumpTurbine │
│ ActivePowerVariable__HydroPumpTurbine     │
│ ReservationVariable__HydroPumpTurbine     │
└───────────────────────────────────────────┘

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

julia> var = read_variable(
           res,
           "ActivePowerVariable__HydroPumpTurbine";
           table_format = TableFormat.WIDE,
       )24×2 DataFrame
 Row  DateTime             Bat_pump 
      DateTime             Float64? 
─────┼───────────────────────────────
   1 │ 2024-01-01T00:00:00     200.0
   2 │ 2024-01-01T01:00:00     200.0
   3 │ 2024-01-01T02:00:00     200.0
   4 │ 2024-01-01T03:00:00     200.0
   5 │ 2024-01-01T04:00:00     200.0
   6 │ 2024-01-01T05:00:00     200.0
   7 │ 2024-01-01T06:00:00     200.0
   8 │ 2024-01-01T07:00:00     200.0
  ⋮  │          ⋮              ⋮
  18 │ 2024-01-01T17:00:00     200.0
  19 │ 2024-01-01T18:00:00     200.0
  20 │ 2024-01-01T19:00:00     200.0
  21 │ 2024-01-01T20:00:00     200.0
  22 │ 2024-01-01T21:00:00     200.0
  23 │ 2024-01-01T22:00:00     200.0
  24 │ 2024-01-01T23:00:00     200.0
                       9 rows omitted

its pump usage

julia> var = read_variable(
           res,
           "ActivePowerPumpVariable__HydroPumpTurbine";
           table_format = TableFormat.WIDE,
       )24×2 DataFrame
 Row  DateTime             Bat_pump 
      DateTime             Float64? 
─────┼───────────────────────────────
   1 │ 2024-01-01T00:00:00       0.0
   2 │ 2024-01-01T01:00:00       0.0
   3 │ 2024-01-01T02:00:00       0.0
   4 │ 2024-01-01T03:00:00       0.0
   5 │ 2024-01-01T04:00:00       0.0
   6 │ 2024-01-01T05:00:00       0.0
   7 │ 2024-01-01T06:00:00       0.0
   8 │ 2024-01-01T07:00:00       0.0
  ⋮  │          ⋮              ⋮
  18 │ 2024-01-01T17:00:00       0.0
  19 │ 2024-01-01T18:00:00       0.0
  20 │ 2024-01-01T19:00:00       0.0
  21 │ 2024-01-01T20:00:00       0.0
  22 │ 2024-01-01T21:00:00       0.0
  23 │ 2024-01-01T22:00:00       0.0
  24 │ 2024-01-01T23:00:00       0.0
                       9 rows omitted

and stored energy:

julia> var =
           read_variable(res, "EnergyVariable__HydroPumpTurbine"; table_format = TableFormat.WIDE)ERROR: EnergyVariable__HydroPumpTurbine is not stored