Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

fix benders when the master has integral variables #866

Merged
merged 1 commit into from
May 4, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion src/Algorithm/basic/solveipform.jl
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ _optimizer_params(::Formulation, algo::SolveIpForm, ::UserOptimizer) = algo.user
_optimizer_params(form::Formulation, algo::SolveIpForm, ::CustomOptimizer) = getinner(getoptimizer(form, algo.optimizer_id))
_optimizer_params(::Formulation, ::SolveIpForm, ::NoOptimizer) = nothing

function run!(algo::SolveIpForm, env::Env, form::Formulation, input::OptimizationState)
function run!(algo::SolveIpForm, env::Env, form::Formulation, input::OptimizationState, optimizer_id = 1)
opt = getoptimizer(form, algo.optimizer_id)
params = _optimizer_params(form, algo, opt)
if params !== nothing
Expand Down
26 changes: 22 additions & 4 deletions src/Algorithm/benders/default.jl
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,7 @@ Benders.get_reform(ctx::BendersContext) = ctx.reform
Benders.get_master(ctx::BendersContext) = getmaster(ctx.reform)
Benders.get_benders_subprobs(ctx::BendersContext) = get_benders_sep_sps(ctx.reform)
struct BendersMasterResult{F,S}
ip_solver::Bool
result::OptimizationState{F,S}
end

Expand All @@ -45,9 +46,15 @@ function Benders.is_unbounded(master_res::BendersMasterResult)
return status == ClB.UNBOUNDED
end

Benders.get_primal_sol(master_res::BendersMasterResult) = get_best_lp_primal_sol(master_res.result)
function Benders.get_primal_sol(master_res::BendersMasterResult)
if master_res.ip_solver
return get_best_ip_primal_sol(master_res.result)
end
return get_best_lp_primal_sol(master_res.result)
end

Benders.get_dual_sol(master_res::BendersMasterResult) = get_best_lp_dual_sol(master_res.result)
Benders.get_obj_val(master_res::BendersMasterResult) = getvalue(get_best_lp_primal_sol(master_res.result))
Benders.get_obj_val(master_res::BendersMasterResult) = getvalue(Benders.get_primal_sol(master_res))


function _reset_second_stage_cost_var_inc_vals(ctx::BendersContext)
Expand All @@ -69,14 +76,25 @@ end
function Benders.optimize_master_problem!(master, ctx::BendersContext, env)
rm_input = OptimizationState(master)
opt_state = run!(ctx.restr_master_solve_alg, env, master, rm_input, ctx.restr_master_optimizer_id)
master_res = BendersMasterResult(opt_state)
ip_solver = typeof(ctx.restr_master_solve_alg) <: SolveIpForm
master_res = BendersMasterResult(ip_solver, opt_state)
_update_second_stage_cost_var_inc_vals(ctx, master_res)
return master_res
end

function Benders.treat_unbounded_master_problem!(master, ctx::BendersContext, env)
mast_result = nothing
certificate = false
ip_solver = typeof(ctx.restr_master_solve_alg) <: SolveIpForm

# In the unbounded case, to get a dual infeasibility certificate, we need to relax the
# integrality and solve the master again. (at least with GLPK & the current implementation of SolveIpForm)
if ip_solver
relax_integrality!(master)
rm_input = OptimizationState(master)
opt_state = run!(SolveLpForm(get_dual_solution = true), env, master, rm_input, ctx.restr_master_optimizer_id)
enforce_integrality!(master)
end

# We can derive a cut from the extreme ray
certificates = MathProg.get_dual_infeasibility_certificate(master, getoptimizer(master, ctx.restr_master_optimizer_id))
Expand All @@ -88,7 +106,7 @@ function Benders.treat_unbounded_master_problem!(master, ctx::BendersContext, en
)
set_ip_primal_sol!(opt_state, first(certificates))
set_lp_primal_sol!(opt_state, first(certificates))
mast_result = BendersMasterResult(opt_state)
mast_result = BendersMasterResult(ip_solver, opt_state)
_update_second_stage_cost_var_inc_vals(ctx, mast_result)
certificate = true
else
Expand Down
21 changes: 17 additions & 4 deletions src/Algorithm/benders/printer.jl
Original file line number Diff line number Diff line change
Expand Up @@ -79,7 +79,20 @@ function Benders.optimize_master_problem!(master, ctx::BendersPrinterContext, en
return result
end

Benders.treat_unbounded_master_problem!(master, ctx::BendersPrinterContext, env) = Benders.treat_unbounded_master_problem!(master, ctx.inner, env)
function Benders.treat_unbounded_master_problem!(master, ctx::BendersPrinterContext, env)
result, c = Benders.treat_unbounded_master_problem!(master, ctx.inner, env)
if ctx.debug_print_master || ctx.debug_print_master_primal_solution || ctx.debug_print_master_dual_solution
println(crayon"bold blue", repeat('-', 80), crayon"reset")
println(crayon"bold underline blue", "Treat unbounded master", crayon"reset")
@show master
print(crayon"bold underline blue", "Master primal solution:", crayon"!bold !underline")
@show Benders.get_primal_sol(result)
print(crayon"reset")
println(crayon"bold blue", repeat('-', 80), crayon"reset")
end
return result, c
end

Benders.set_second_stage_var_costs_to_zero!(ctx::BendersPrinterContext) = Benders.set_second_stage_var_costs_to_zero!(ctx.inner)
Benders.reset_second_stage_var_costs!(ctx::BendersPrinterContext) = Benders.reset_second_stage_var_costs!(ctx.inner)
Benders.update_sp_rhs!(ctx::BendersPrinterContext, sp, primal_sol) = Benders.update_sp_rhs!(ctx.inner, sp, primal_sol)
Expand All @@ -90,14 +103,14 @@ function Benders.optimize_separation_problem!(ctx::BendersPrinterContext, sp::Fo
if ctx.debug_print_subproblem || ctx.debug_print_subproblem_primal_solution || ctx.debug_print_subproblem_dual_solution
println(crayon"bold green", repeat('-', 80), crayon"reset")
end
ctx.sp_elapsed_time = @elapsed begin
result = Benders.optimize_separation_problem!(ctx.inner, sp, env)
end
if ctx.debug_print_subproblem
print(crayon"bold underline green", "Separation problem:", crayon"!bold !underline")
@show sp
print(crayon"reset")
end
ctx.sp_elapsed_time = @elapsed begin
result = Benders.optimize_separation_problem!(ctx.inner, sp, env)
end
if ctx.debug_print_subproblem_primal_solution
print(crayon"bold underline green", "Separation problem primal solution:", crayon"!bold !underline")
@show Benders.get_primal_sol(result)
Expand Down
24 changes: 15 additions & 9 deletions test/unit/Benders/benders_default.jl
Original file line number Diff line number Diff line change
Expand Up @@ -206,18 +206,24 @@ function benders_iter_default_A_integer()

alg = Coluna.Algorithm.BendersCutGeneration(
max_nb_iterations = 10,
separation_solve_alg = Coluna.Algorithm.SolveIpForm()
restr_master_solve_alg = Coluna.Algorithm.SolveIpForm()
)
ctx = Coluna.Algorithm.BendersPrinterContext(
reform, alg;
print = true
print = true,
# debug_print_master = true,
# debug_print_master_primal_solution = true,
# debug_print_master_dual_solution = true,
# debug_print_subproblem = true,
# debug_print_subproblem_primal_solution = true,
# debug_print_subproblem_dual_solution = true
)
Coluna.set_optim_start_time!(env)

result = Coluna.Benders.run_benders_loop!(ctx, env)
@test result.mlp ≈ 4.0
end
register!(unit_tests, "benders_default", benders_iter_default_A_integer; x = true)
register!(unit_tests, "benders_default", benders_iter_default_A_integer)


# B with continuous first stage finds optimal solution
Expand Down Expand Up @@ -254,9 +260,9 @@ register!(unit_tests, "benders_default", benders_iter_default_B_continuous)
# B with integer first stage finds optimal solution
# Error occurs during test TODO fix
# expected output:
# mlp = 7.083333333333333
# x1 = 0.0, x2 = 1.0
# y1 = 3.1666666666666665, y2 = 0.3333333333333333
# mlp = 7
# x1 = 0.0, x2 = 2.0
# y1 = 2, y2 = 0
function benders_iter_default_B_integer()
#env, reform = benders_simple_example()
env, reform = benders_form_B()
Expand All @@ -271,7 +277,7 @@ function benders_iter_default_B_integer()

alg = Coluna.Algorithm.BendersCutGeneration(
max_nb_iterations = 10,
separation_solve_alg = Coluna.Algorithm.SolveIpForm()
restr_master_solve_alg = Coluna.Algorithm.SolveIpForm()
)
ctx = Coluna.Algorithm.BendersPrinterContext(
reform, alg;
Expand All @@ -280,6 +286,6 @@ function benders_iter_default_B_integer()
Coluna.set_optim_start_time!(env)

result = Coluna.Benders.run_benders_loop!(ctx, env)
@test result.mlp ≈ 7.083333333333333
@test result.mlp ≈ 7
end
register!(unit_tests, "benders_default", benders_iter_default_B_integer; x = true)
register!(unit_tests, "benders_default", benders_iter_default_B_integer)