Infinite horizon hydro-thermal

This tutorial was generated using Literate.jl. Download the source as a .jl file. Download the source as a .ipynb file.

using SDDP, HiGHS, Test, Statisticsfunction infinite_hydro_thermal(; cut_type)    Ω = [        (inflow = 0.0, demand = 7.5),        (inflow = 5.0, demand = 5),        (inflow = 10.0, demand = 2.5),    ]    graph = SDDP.Graph(        :root_node,        [:week],        [(:root_node => :week, 1.0), (:week => :week, 0.9)],    )    model = SDDP.PolicyGraph(        graph;        lower_bound = 0,        optimizer = HiGHS.Optimizer,    ) do subproblem, node        @variable(            subproblem,            5.0 <= reservoir <= 15.0,            SDDP.State,            initial_value = 10.0        )        @variables(subproblem, begin            thermal_generation >= 0            hydro_generation >= 0            spill >= 0            inflow            demand        end)        @constraints(            subproblem,            begin                reservoir.out == reservoir.in - hydro_generation - spill + inflow                hydro_generation + thermal_generation == demand            end        )        @stageobjective(subproblem, 10 * spill + thermal_generation)        SDDP.parameterize(subproblem, Ω) do ω            fix(inflow, ω.inflow)            fix(demand, ω.demand)            return        end    end    SDDP.train(        model;        cut_type = cut_type,        log_frequency = 100,        sampling_scheme = SDDP.InSampleMonteCarlo(; terminate_on_cycle = true),        parallel_scheme = SDDP.Serial(),        cycle_discretization_delta = 0.1,    )    @test SDDP.calculate_bound(model)  119.167 atol = 0.1    results = SDDP.simulate(model, 500)    objectives =        [sum(s[:stage_objective] for s in simulation) for simulation in results]    sample_mean = round(Statistics.mean(objectives); digits = 2)    sample_ci = round(1.96 * Statistics.std(objectives) / sqrt(500); digits = 2)    println("Confidence_interval = $(sample_mean) ± $(sample_ci)")    @test sample_mean - sample_ci <= 119.167 <= sample_mean + sample_ci    returnendinfinite_hydro_thermal(; cut_type = SDDP.SINGLE_CUT)infinite_hydro_thermal(; cut_type = SDDP.MULTI_CUT)
-------------------------------------------------------------------
         SDDP.jl (c) Oscar Dowson and contributors, 2017-26
-------------------------------------------------------------------
problem
  nodes           : 1
  state variables : 1
  scenarios       : Inf
  existing cuts   : false
options
  solver          : serial mode
  risk measure    : SDDP.Expectation()
  sampling scheme : SDDP.InSampleMonteCarlo
subproblem structure
  VariableRef                             : [8, 8]
  AffExpr in MOI.EqualTo{Float64}         : [2, 2]
  VariableRef in MOI.GreaterThan{Float64} : [5, 5]
  VariableRef in MOI.LessThan{Float64}    : [1, 1]
numerical stability report
  matrix range     [1e+00, 1e+00]
  objective range  [1e+00, 1e+01]
  bounds range     [5e+00, 2e+01]
  rhs range        [0e+00, 0e+00]
-------------------------------------------------------------------
 iteration    simulation      bound        time (s)     solves  pid
-------------------------------------------------------------------
       100   2.500000e+01  1.189023e+02  3.425889e-01      1946   1
       200   2.500000e+01  1.191635e+02  4.185739e-01      3920   1
       300   0.000000e+00  1.191666e+02  5.002570e-01      5902   1
       330   2.500000e+01  1.191667e+02  5.167680e-01      6224   1
-------------------------------------------------------------------
status         : simulation_stopping
total time (s) : 5.167680e-01
total solves   : 6224
best bound     :  1.191667e+02
numeric issues : 0
-------------------------------------------------------------------

Confidence_interval = 120.33 ± 12.06
-------------------------------------------------------------------
         SDDP.jl (c) Oscar Dowson and contributors, 2017-26
-------------------------------------------------------------------
problem
  nodes           : 1
  state variables : 1
  scenarios       : Inf
  existing cuts   : false
options
  solver          : serial mode
  risk measure    : SDDP.Expectation()
  sampling scheme : SDDP.InSampleMonteCarlo
subproblem structure
  VariableRef                             : [8, 8]
  AffExpr in MOI.EqualTo{Float64}         : [2, 2]
  VariableRef in MOI.GreaterThan{Float64} : [5, 5]
  VariableRef in MOI.LessThan{Float64}    : [1, 1]
numerical stability report
  matrix range     [1e+00, 1e+00]
  objective range  [1e+00, 1e+01]
  bounds range     [5e+00, 2e+01]
  rhs range        [0e+00, 0e+00]
-------------------------------------------------------------------
 iteration    simulation      bound        time (s)     solves  pid
-------------------------------------------------------------------
       100   0.000000e+00  1.191283e+02  1.173790e-01      2874   1
       200   2.500000e+00  1.191666e+02  2.081449e-01      4855   1
       282   7.500000e+00  1.191667e+02  2.612770e-01      5733   1
-------------------------------------------------------------------
status         : simulation_stopping
total time (s) : 2.612770e-01
total solves   : 5733
best bound     :  1.191667e+02
numeric issues : 0
-------------------------------------------------------------------

Confidence_interval = 116.06 ± 13.65