StochDynamicProgramming: the multistock problem
This tutorial was generated using Literate.jl. Download the source as a .jl file. Download the source as a .ipynb file.
This example comes from StochDynamicProgramming.jl.
using SDDP, HiGHS, Testfunction test_multistock_example() model = SDDP.LinearPolicyGraph(; stages = 5, lower_bound = -5.0, optimizer = HiGHS.Optimizer, ) do subproblem, stage @variable( subproblem, 0 <= stock[i=1:3] <= 1, SDDP.State, initial_value = 0.5 ) @variables(subproblem, begin 0 <= control[i=1:3] <= 0.5 ξ[i=1:3] # Dummy for RHS noise. end) @constraints( subproblem, begin sum(control) - 0.5 * 3 <= 0 [i = 1:3], stock[i].out == stock[i].in + control[i] - ξ[i] end ) Ξ = collect( Base.product((0.0, 0.15, 0.3), (0.0, 0.15, 0.3), (0.0, 0.15, 0.3)), )[:] SDDP.parameterize(subproblem, Ξ) do ω fix.(ξ, ω) return end @stageobjective(subproblem, (sin(3 * stage) - 1) * sum(control)) end SDDP.train( model; iteration_limit = 100, cut_type = SDDP.SINGLE_CUT, log_frequency = 10, ) @test SDDP.calculate_bound(model) ≈ -4.349 atol = 0.01 simulation_results = SDDP.simulate(model, 5000) @test length(simulation_results) == 5000 μ = SDDP.Statistics.mean( sum(data[:stage_objective] for data in simulation) for simulation in simulation_results ) @test μ ≈ -4.349 atol = 0.1 returnendtest_multistock_example()-------------------------------------------------------------------
SDDP.jl (c) Oscar Dowson and contributors, 2017-26
-------------------------------------------------------------------
problem
nodes : 5
state variables : 3
scenarios : 1.43489e+07
existing cuts : false
options
solver : serial mode
risk measure : SDDP.Expectation()
sampling scheme : SDDP.InSampleMonteCarlo
subproblem structure
VariableRef : [13, 13]
AffExpr in MOI.EqualTo{Float64} : [3, 3]
AffExpr in MOI.LessThan{Float64} : [1, 1]
VariableRef in MOI.GreaterThan{Float64} : [7, 7]
VariableRef in MOI.LessThan{Float64} : [6, 7]
numerical stability report
matrix range [1e+00, 1e+00]
objective range [3e-01, 2e+00]
bounds range [5e-01, 5e+00]
rhs range [2e+00, 2e+00]
-------------------------------------------------------------------
iteration simulation bound time (s) solves pid
-------------------------------------------------------------------
10 -4.724951e+00 -4.443273e+00 2.210419e-01 1400 1
20 -3.477232e+00 -4.396582e+00 3.476410e-01 2800 1
30 -4.071432e+00 -4.381910e+00 4.813440e-01 4200 1
40 -4.425331e+00 -4.370204e+00 6.228900e-01 5600 1
50 -4.238128e+00 -4.364287e+00 7.666900e-01 7000 1
60 -3.991595e+00 -4.358356e+00 9.158289e-01 8400 1
70 -4.464815e+00 -4.354899e+00 1.064723e+00 9800 1
80 -3.912294e+00 -4.352960e+00 1.224767e+00 11200 1
90 -2.784732e+00 -4.350166e+00 1.385526e+00 12600 1
100 -2.789585e+00 -4.349355e+00 1.547386e+00 14000 1
-------------------------------------------------------------------
status : iteration_limit
total time (s) : 1.547386e+00
total solves : 14000
best bound : -4.349355e+00
numeric issues : 0
-------------------------------------------------------------------