HydroWaterFactorModel model tutorial

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_energy")

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

reservoir = only(get_components(HydroReservoir, sys))

HydroReservoir: HydroEnergyReservoir__reservoir:
   name: HydroEnergyReservoir__reservoir
   available: true
   storage_level_limits: (min = 0.0, max = 5000.0)
   initial_level: 0.5
   spillage_limits: nothing
   inflow: 4.0
   outflow: 0.0
   level_targets: 1.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: true

Set the storage level limits using set_storage_level_limits!. Here the limits are set between $4000 and 6000 m^3$ .

set_storage_level_limits!(reservoir, (min = 4000, max = 6000))

(min = 4000, max = 6000)

Set the lower bound of reservoir volume using set_level_targets!. Here the level target is set at $0.9 \cdot 6000 = 5400 m^3$.

set_level_targets!(reservoir, 0.9)

0.9

Set the headtovolume factor to 1.0 using set_head_to_volume_factor!, since is an energy model and no conversion is needed.

set_head_to_volume_factor!(reservoir, LinearCurve(1.0))

InfrastructureSystems.LinearCurve (a type of InfrastructureSystems.InputOutputCurve) where function is: f(x) = 1.0 x + 0.0

Decision Model

Setting up the formulations based on PowerSimulations.jl:

template = ProblemTemplate(CopperPlatePowerModel)
set_device_model!(template, ThermalStandard, ThermalBasicDispatch)
set_device_model!(template, PowerLoad, StaticPowerLoad)

but, now we also include the HydroReservoir and HydroTurbine using HydroWaterFactorModel:

set_device_model!(template, HydroReservoir, HydroWaterFactorModel)
set_device_model!(template, HydroTurbine, HydroWaterFactorModel)

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:

results = 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:

power =
    read_variable(
        results,
        "ActivePowerVariable__HydroTurbine";
        table_format = TableFormat.WIDE,
    )

or the volume level of the reservoir

volume =
    read_variable(
        results,
        "HydroReservoirVolumeVariable__HydroReservoir";
        table_format = TableFormat.WIDE,
    )

Note that the final reservoir level is between the set level target and the maximum storage level.