odow / SDDP.jl

A JuMP extension for Stochastic Dual Dynamic Programming
https://sddp.dev
Other
305 stars 61 forks source link

Clarify documentation of non monotonic bound #746

Closed dannyopts closed 4 months ago

dannyopts commented 4 months ago

Hi,

I have noticed that the bound in the iteration log is sometimes not monotonic increasing, although the docs state "bound: this is a lower bound (upper if maximizing) for the value of the optimal policy. This should be monotonically improving (increasing if minimizing, decreasing if maximizing)." https://sddp.dev/stable/tutorial/first_steps/. I have attached an iteration log at the bottom of this issue which demonstrates this.

The policy graph I am using is cyclical, and the bound is mostly monotonicly increasing, but not always.

Does the bound not being monotonic suggest numerical errors, or that I have some how violated one of the assumptions of the model?

I have generated the lp files for the first two iterations and attached below and it seems that the cut made in the first iteration is subsequently removed. I dont understand the logic for cut deletion, but it seems strange to me that a cut would be deleted between iteration 1 and 2, particularly when the the new one doesnt dominate the previous one?

example iteration log:

iteration simulation bound time (s) solves pid

     1   1.026923e+03  1.209591e+02  5.936684e+00      1272   1
     2   1.310886e+03  3.369082e+01  2.047813e+01     14407   1
     3   1.060614e+03  1.686153e+02  2.192366e+01     15679   1
     4   6.421468e+02  1.777345e+02  2.616833e+01     16951   1

LP file for single child of root node after iteration 1:


minimize
obj: 1950 wind_size_out + 1650 solar_size_out + 366 battery_size_out + 3740 _electro_size_out + 1 x11
subject to
c1: 11784.978488769622 _electro_size_out + 1 x12 >= 381.1497345414217
c2: 1 x11 - 1 x12 >= 0
Bounds
wind_size_in = 0
solar_size_in = 0
battery_size_in = 0
battery_state_in = 0
_electro_size_in = 0
0 <= wind_size_out <= 8
0 <= solar_size_out <= 8
0 <= battery_size_out <= 4
0 <= battery_state_out <= 4
0 <= _electro_size_out <= 10
x11 >= 0
x12 >= 0
End

LP file for single child of root node after iteration 2:


minimize
obj: 1950 wind_size_out + 1650 solar_size_out + 366 battery_size_out + 3740 _electro_size_out + 1 x11
subject to
c1: 1 x11 - 1 x12 >= 0
c2: 17808.63265656682 wind_size_out + 4179.291345533028 solar_size_out + 4566.697925341333 battery_size_out - 669.6530173846298 _electro_size_out + 1 x12 >= 419.9917767751505

Bounds
wind_size_in = 0
solar_size_in = 0
battery_size_in = 0
battery_state_in = 0
_electro_size_in = 0
0 <= wind_size_out <= 8
0 <= solar_size_out <= 8
0 <= battery_size_out <= 4
0 <= battery_state_out <= 4
0 <= _electro_size_out <= 10
x11 >= 0
x12 >= 0
End

Just incase I have done something stupid I have included the steps I went through to get these files below.

The policy graph is generated as

graph = SDDP.Graph(
     0,
     [i for i in 1:54],
     vcat(
          [(0 => 1, 1.0)],
          [(1 => 2, 1.0)],
          [(2 => i, 1/52) for i in 3:54],
          [(i => i+1, 0.99) for i in 3:53],
          [(54 => 3, 0.99)],
     )
)

construct the model

model = SDDP.PolicyGraph(
    graph,
    sense = :Min,
    lower_bound = 0.0,
    optimizer = HiGHS.Optimizer,
) do sp, node
....

Train for a single iteration

SDDP.train(model; cut_type=SDDP.MULTI_CUT, log_every_iteration = true, iteration_limit=1, sampling_scheme = SDDP.InSampleMonteCarlo(;max_depth = 5, terminate_on_dummy_leaf = false))

Write out subproblem for single child of root node

SDDP.write_subproblem_to_file(model[1], "subproblem_1_1.lp")

Train for another iteration

SDDP.train(model; cut_type=SDDP.MULTI_CUT, log_every_iteration = true, iteration_limit=1, sampling_scheme = SDDP.InSampleMonteCarlo(;max_depth = 5, terminate_on_dummy_leaf = false))

Write out subproblem for single child of root node

SDDP.write_subproblem_to_file(model[1], "subproblem_1_2.lp")

Thank you so much for all your hard work with this awesome library and also for investing so much time into writing some brilliant documentation 😍

dannyopts commented 4 months ago

Having a poke around in https://github.com/odow/SDDP.jl/blob/496c3c30d66c2961485946095f5f0e2015250db2/src/plugins/bellman_functions.jl#L135 but not yet having fully parsed it yet, it would seem that a cut is deleted if it is not the best lower bound for one atleast one of the sampled states?

But if my thinking is correct, the fact it is not providing a lower bound on the sampled points doesnt mean it is not providing a lower bound for other points in the domain which could mean we end up with a previously infeasible point becoming feasible (as we have seen in the case above).

This seems fine and the comment in the documentation should be caveated to include the fact the bound can decrease due to cut deletion?

odow commented 4 months ago

Yes, as you figured out, we implement a form of cut pruning that--in rare cases, particularly with cyclic graphs--can result in the lower bound being non-monotonic.

You can pass SDDP.train(model; cut_deletion_minimum = 1_000) (or some other large number) to disable this. https://github.com/odow/SDDP.jl/blob/496c3c30d66c2961485946095f5f0e2015250db2/src/algorithm.jl#L906

the documentation should be caveated to include the fact the bound can decrease due to cut deletion?

Okay :smile:

odow commented 4 months ago

Thanks!