From b694018d7c5d41950ee1fce124390a93ad82b0ee Mon Sep 17 00:00:00 2001 From: culi122 <40186418+cjuli1@users.noreply.github.com> Date: Tue, 13 Aug 2024 10:20:23 +0200 Subject: [PATCH 1/3] test benders with warm start in sub problems --- src/prob_stoch.jl | 74 +++++++++++++++++++++++++++++++++-------------- 1 file changed, 52 insertions(+), 22 deletions(-) diff --git a/src/prob_stoch.jl b/src/prob_stoch.jl index e995522..564cc11 100644 --- a/src/prob_stoch.jl +++ b/src/prob_stoch.jl @@ -164,6 +164,11 @@ function save_storagevalues(prob, cuts, storagevalues) return end +""" +Two possibilities to warm start benders solve: +- Start with solving master problem if cuts from previous step are available. NB! we reuse cuts although scenarios can have changed +- Start with solving subproblems with start reservoirs from db.startstates if cuts from previous step are NOT available +""" function solve_benders(stepnr, subix) db = get_local_db() settings = get_settings(db) @@ -182,7 +187,6 @@ function solve_benders(stepnr, subix) while !((abs((ub-lb)/ub) < reltol) || abs(ub-lb) < 1) maintiming[4] += @elapsed begin - count == 0 && TuLiPa.setwarmstart!(mp.prob, false) if (count == 1 && cutreuse) TuLiPa.updatecuts!(mp.prob, mp.cuts) elseif count != 0 @@ -191,43 +195,52 @@ function solve_benders(stepnr, subix) end maintiming[2] += @elapsed begin - if cutreuse # try to reuse cuts from last time step, NB! we reuse cuts although scenarios can have changed - try + if getnumcuts(cuts) != 0 + count == 0 && TuLiPa.setwarmstart!(mp.prob, false) + if cutreuse + try + TuLiPa.solve!(mp.prob) + count == 0 && TuLiPa.clearcuts!(mp.cuts) + catch + count == 0 && println("Retrying first iteration without cuts from last time step") + count > 0 && println("Restarting iterations without cuts from last time step") + TuLiPa.clearcuts!(mp.prob, mp.cuts) + cutreuse = false + count = 0 + end + else TuLiPa.solve!(mp.prob) - catch - count == 0 && println("Retrying first iteration without cuts from last time step") - count > 0 && println("Restarting iterations without cuts from last time step") - TuLiPa.clearcuts!(mp.prob, mp.cuts) - TuLiPa.solve!(mp.prob) - cutreuse = false end - else - TuLiPa.solve!(mp.prob) + count == 0 && TuLiPa.setwarmstart!(mp.prob, true) end end maintiming[4] += @elapsed begin - lb = TuLiPa.getvarvalue(mp.prob, TuLiPa.getfuturecostvarid(mp.cuts), 1) - ub = 0.0 - - count == 0 && TuLiPa.setwarmstart!(mp.prob, true) - (count == 0 && cutreuse) && TuLiPa.clearcuts!(mp.cuts) - TuLiPa.getoutgoingstates!(mp.prob, mp.states) - cutix = TuLiPa.getcutix(mp.cuts) + 1 - if cutix > TuLiPa.getmaxcuts(mp.cuts) - cutix = 1 + if getnumcuts(cuts) != 0 + lb = TuLiPa.getvarvalue(mp.prob, TuLiPa.getfuturecostvarid(mp.cuts), 1) + TuLiPa.getoutgoingstates!(mp.prob, mp.states) + count += 1 end end futures = [] @sync for (_scenix, _subix, _core) in db.dist_sp if _subix == subix - f = @spawnat _core solve_sp(_scenix, _subix, mp.states) + if getnumcuts(cuts) != 0 + f = @spawnat _core solve_sp(_scenix, _subix, mp.states) + else + f = @spawnat _core solve_sp_with_startreservoirs(_scenix, _subix) + end push!(futures, f) end end maintiming[4] += @elapsed begin + ub = 0.0 + cutix = TuLiPa.getcutix(mp.cuts) + 1 + if cutix > TuLiPa.getmaxcuts(mp.cuts) + cutix = 1 + end for future in futures scenix, objectivevalue, scenslopes, scenconstant = fetch(future) @@ -237,7 +250,6 @@ function solve_benders(stepnr, subix) end TuLiPa.updatecutparameters!(mp.prob, mp.cuts) - count += 1 end end maintiming[5] = count @@ -262,6 +274,24 @@ function solve_sp(scenix, subix, states) return (scenix, objectivevalue, scenslopes, scenconstant) end +function solve_sp_with_startreservoirs(scenix, subix) + db = get_local_db() + + sp = db.sp[(scenix, subix)] + maintiming = sp.div[MainTiming] + + maintiming[3] += @elapsed TuLiPa.setstartstates!(sp.prob, getstorages(getobjects(sp.prob)), db.startstates) + maintiming[2] += @elapsed TuLiPa.solve!(sp.prob) + maintiming[3] += @elapsed begin + get_scencutparameters!(sp, states) + + objectivevalue = TuLiPa.getobjectivevalue(sp.prob) + scenslopes = sp.scenslopes + scenconstant = sp.scenconstant + end + return (scenix, objectivevalue, scenslopes, scenconstant) +end + # TODO: Simplify TuLiPa version of getscencutparameters? function get_scencutparameters!(sp::ScenarioProblem, states::Dict{TuLiPa.StateVariableInfo, Float64}) sp.scenconstant = TuLiPa.getobjectivevalue(sp.prob) From 8817f6e0cceb46e72eefc80d94c9760e675e699c Mon Sep 17 00:00:00 2001 From: culi122 <40186418+cjuli1@users.noreply.github.com> Date: Tue, 13 Aug 2024 11:41:08 +0200 Subject: [PATCH 2/3] fix getnumcuts --- src/prob_stoch.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/prob_stoch.jl b/src/prob_stoch.jl index 564cc11..1a6e773 100644 --- a/src/prob_stoch.jl +++ b/src/prob_stoch.jl @@ -195,7 +195,7 @@ function solve_benders(stepnr, subix) end maintiming[2] += @elapsed begin - if getnumcuts(cuts) != 0 + if TuLiPa.getnumcuts(cuts) != 0 count == 0 && TuLiPa.setwarmstart!(mp.prob, false) if cutreuse try @@ -216,7 +216,7 @@ function solve_benders(stepnr, subix) end maintiming[4] += @elapsed begin - if getnumcuts(cuts) != 0 + if TuLiPa.getnumcuts(cuts) != 0 lb = TuLiPa.getvarvalue(mp.prob, TuLiPa.getfuturecostvarid(mp.cuts), 1) TuLiPa.getoutgoingstates!(mp.prob, mp.states) count += 1 @@ -226,7 +226,7 @@ function solve_benders(stepnr, subix) futures = [] @sync for (_scenix, _subix, _core) in db.dist_sp if _subix == subix - if getnumcuts(cuts) != 0 + if TuLiPa.getnumcuts(cuts) != 0 f = @spawnat _core solve_sp(_scenix, _subix, mp.states) else f = @spawnat _core solve_sp_with_startreservoirs(_scenix, _subix) From e2dac9849449820735844ae19d8720a0f3f2c394 Mon Sep 17 00:00:00 2001 From: culi122 <40186418+cjuli1@users.noreply.github.com> Date: Tue, 13 Aug 2024 14:51:47 +0200 Subject: [PATCH 3/3] Use TuLiPa. --- src/prob_stoch.jl | 24 +++++++++++++----------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/src/prob_stoch.jl b/src/prob_stoch.jl index 1a6e773..4e25190 100644 --- a/src/prob_stoch.jl +++ b/src/prob_stoch.jl @@ -185,7 +185,7 @@ function solve_benders(stepnr, subix) lb = mp.cuts.lower_bound reltol = settings["problems"]["stochastic"]["reltol"] # relative tolerance - while !((abs((ub-lb)/ub) < reltol) || abs(ub-lb) < 1) + while !((abs((ub-lb)/ub) < reltol) || abs(ub-lb) < 1) && count < 15 maintiming[4] += @elapsed begin if (count == 1 && cutreuse) TuLiPa.updatecuts!(mp.prob, mp.cuts) @@ -195,12 +195,13 @@ function solve_benders(stepnr, subix) end maintiming[2] += @elapsed begin - if TuLiPa.getnumcuts(cuts) != 0 + if TuLiPa.getnumcuts(mp.cuts) != 0 count == 0 && TuLiPa.setwarmstart!(mp.prob, false) if cutreuse try TuLiPa.solve!(mp.prob) count == 0 && TuLiPa.clearcuts!(mp.cuts) + count += 1 catch count == 0 && println("Retrying first iteration without cuts from last time step") count > 0 && println("Restarting iterations without cuts from last time step") @@ -210,26 +211,27 @@ function solve_benders(stepnr, subix) end else TuLiPa.solve!(mp.prob) + count += 1 end count == 0 && TuLiPa.setwarmstart!(mp.prob, true) + lb = TuLiPa.getvarvalue(mp.prob, TuLiPa.getfuturecostvarid(mp.cuts), 1) + TuLiPa.getoutgoingstates!(mp.prob, mp.states) end end maintiming[4] += @elapsed begin - if TuLiPa.getnumcuts(cuts) != 0 - lb = TuLiPa.getvarvalue(mp.prob, TuLiPa.getfuturecostvarid(mp.cuts), 1) - TuLiPa.getoutgoingstates!(mp.prob, mp.states) - count += 1 + if TuLiPa.getnumcuts(mp.cuts) != 0 + end end futures = [] @sync for (_scenix, _subix, _core) in db.dist_sp if _subix == subix - if TuLiPa.getnumcuts(cuts) != 0 + if TuLiPa.getnumcuts(mp.cuts) != 0 f = @spawnat _core solve_sp(_scenix, _subix, mp.states) else - f = @spawnat _core solve_sp_with_startreservoirs(_scenix, _subix) + f = @spawnat _core solve_sp_with_startreservoirs(_scenix, _subix, mp.states) end push!(futures, f) end @@ -274,16 +276,16 @@ function solve_sp(scenix, subix, states) return (scenix, objectivevalue, scenslopes, scenconstant) end -function solve_sp_with_startreservoirs(scenix, subix) +function solve_sp_with_startreservoirs(scenix, subix, states) db = get_local_db() sp = db.sp[(scenix, subix)] maintiming = sp.div[MainTiming] - maintiming[3] += @elapsed TuLiPa.setstartstates!(sp.prob, getstorages(getobjects(sp.prob)), db.startstates) + maintiming[3] += @elapsed set_startstates!(sp.prob, TuLiPa.getstorages(TuLiPa.getobjects(sp.prob)), db.startstates) maintiming[2] += @elapsed TuLiPa.solve!(sp.prob) maintiming[3] += @elapsed begin - get_scencutparameters!(sp, states) + get_scencutparameters!(sp, states) # TODO: Even better replacing states with db.startstates? objectivevalue = TuLiPa.getobjectivevalue(sp.prob) scenslopes = sp.scenslopes