diff --git a/src/Algorithm/colgen/default.jl b/src/Algorithm/colgen/default.jl index 8044b220a..d71246d3b 100644 --- a/src/Algorithm/colgen/default.jl +++ b/src/Algorithm/colgen/default.jl @@ -176,8 +176,8 @@ function ColGen.check_primal_ip_feasibility(master_lp_primal_sol, phase, reform) end # Reduced costs calculation -ColGen.get_orig_costs(ctx::ColGenContext) = ctx.reduced_cost_helper.c -ColGen.get_coef_matrix(ctx::ColGenContext) = ctx.reduced_cost_helper.A +ColGen.get_subprob_var_orig_costs(ctx::ColGenContext) = ctx.reduced_cost_helper.dw_subprob_c +ColGen.get_subprob_var_coef_matrix(ctx::ColGenContext) = ctx.reduced_cost_helper.dw_subprob_A function ColGen.update_sp_vars_red_costs!(ctx::ColGenContext, sp::Formulation{DwSp}, red_costs) for (var_id, _) in getvars(sp) @@ -234,6 +234,7 @@ end struct ColGenPricingResult{F,S} result::OptimizationState{F,S} columns::Vector{GeneratedColumn} + best_red_cost::Float64 end function ColGen.is_infeasible(pricing_res::ColGenPricingResult) @@ -247,9 +248,10 @@ function ColGen.is_unbounded(pricing_res::ColGenPricingResult) end ColGen.get_primal_sols(pricing_res) = pricing_res.columns -ColGen.get_dual_bound(pricing_res) = get_ip_dual_bound(pricing_res.result) +ColGen.get_dual_bound(pricing_res) = pricing_res.best_red_cost -has_improving_red_cost(column::GeneratedColumn) = column.red_cost < 0 +is_improving_red_cost(red_cost) = red_cost < 0 +has_improving_red_cost(column::GeneratedColumn) = is_improving_red_cost(column.red_cost) # In our implementation of `push_in_set!`, we keep only columns that have improving reduced # cost. function ColGen.push_in_set!(pool, column) @@ -270,16 +272,20 @@ function ColGen.optimize_pricing_problem!(ctx::ColGenContext, sp::Formulation{Dw # (C) the contribution of the pure master variables. # Master convexity constraints contribution. + # TODO: talk with fv & Ruslan because this way to take into account convexity constraints has + # drawbacks (numerical stability). lb_dual = master_dual_sol[sp.duty_data.lower_multiplicity_constr_id] ub_dual = master_dual_sol[sp.duty_data.upper_multiplicity_constr_id] # Pure master variables contribution. - # TODO + # TODO (only when stabilization is used otherwise already taken into account by master obj val) generated_columns = GeneratedColumn[] for col in get_ip_primal_sols(opt_state) red_cost = getvalue(col) - lb_dual - ub_dual push!(generated_columns, GeneratedColumn(col, red_cost)) end - return ColGenPricingResult(opt_state, generated_columns) + + best_red_cost = getvalue(get_ip_dual_bound(opt_state)) - lb_dual - ub_dual + return ColGenPricingResult(opt_state, generated_columns, best_red_cost) end \ No newline at end of file diff --git a/src/Algorithm/colgen/utils.jl b/src/Algorithm/colgen/utils.jl index c4d766dce..a2cb233ba 100644 --- a/src/Algorithm/colgen/utils.jl +++ b/src/Algorithm/colgen/utils.jl @@ -67,43 +67,67 @@ function _submatrix( end """ -Extracted information to speed-up calculation of reduced costs of master variables. -We extract from the master only the information we need to compute the reduced cost of DW +Extracted information to speed-up calculation of reduced costs of subproblem representatives +and pure master variables. +We extract from the master the information we need to compute the reduced cost of DW subproblem variables: -- `c` contains the perenial cost of DW subproblem representative variables -- `A` is a submatrix of the master coefficient matrix that involves only DW subproblem +- `dw_subprob_c` contains the perenial cost of DW subproblem representative variables +- `dw_subprob_A` is a submatrix of the master coefficient matrix that involves only DW subproblem representative variables. +We also extract from the master the information we need to compute the reduced cost of pure +master variables: +- `pure_master_c` contains the perenial cost of pure master variables +- `pure_master_A` is a submatrix of the master coefficient matrix that involves only pure master + variables. Calculation is `c - transpose(A) * master_lp_dual_solution`. This information is given to the generic implementation of the column generation algorithm through methods: -- ColGen.get_orig_costs +- ColGen.get_subprob_var_orig_costs - ColGen.get_orig_coefmatrix """ struct ReducedCostsCalculationHelper - c::SparseVector{Float64,VarId} - A::DynamicSparseMatrix{ConstrId,VarId,Float64} + dw_subprob_c::SparseVector{Float64,VarId} + dw_subprob_A::DynamicSparseMatrix{ConstrId,VarId,Float64} + master_c::SparseVector{Float64,VarId} + master_A::DynamicSparseMatrix{ConstrId,VarId,Float64} end -function ReducedCostsCalculationHelper(master) - dwspvar_ids = VarId[] +""" +Function `var_duty_func(var_id)` returns `true` if we want to keep the variable `var_id`; `false` otherwise. +Same for `constr_duty_func(constr_id)`. +""" +function _get_costs_and_coeffs(master, var_duty_func, constr_duty_func) + var_ids = VarId[] peren_costs = Float64[] for var_id in Iterators.keys(getvars(master)) - if iscuractive(master, var_id) && getduty(var_id) <= AbstractMasterRepDwSpVar - push!(dwspvar_ids, var_id) + if iscuractive(master, var_id) && var_duty_func(var_id) + push!(var_ids, var_id) push!(peren_costs, getperencost(master, var_id)) end end - dwsprep_costs = sparsevec(dwspvar_ids, peren_costs, Coluna.MAX_NB_ELEMS) - dwsprep_coefmatrix = _submatrix( + costs = sparsevec(var_ids, peren_costs, Coluna.MAX_NB_ELEMS) + coef_matrix = _submatrix(master, constr_duty_func, var_duty_func) + return costs, coef_matrix +end + +function ReducedCostsCalculationHelper(master) + dw_subprob_c, dw_subprob_A = _get_costs_and_coeffs( master, - constr_id -> !(getduty(constr_id) <= MasterConvexityConstr), - var_id -> getduty(var_id) <= AbstractMasterRepDwSpVar + var_id -> getduty(var_id) <= AbstractMasterRepDwSpVar, + constr_id -> !(getduty(constr_id) <= MasterConvexityConstr) + ) + + master_c, master_A = _get_costs_and_coeffs( + master, + var_id -> getduty(var_id) <= AbstractOriginMasterVar, + constr_id -> !(getduty(constr_id) <= MasterConvexityConstr) ) - return ReducedCostsCalculationHelper(dwsprep_costs, dwsprep_coefmatrix) + + return ReducedCostsCalculationHelper(dw_subprob_c, dw_subprob_A, master_c, master_A) end """ diff --git a/src/ColGen/interface.jl b/src/ColGen/interface.jl index cde2e48b5..8c1857695 100644 --- a/src/ColGen/interface.jl +++ b/src/ColGen/interface.jl @@ -149,13 +149,13 @@ LP solution is integer feasible; `nothing` otherwise. Returns the original cost `c` of subproblems variables. to compute reduced cost `̄c = c - transpose(A) * π`. """ -@mustimplement "ColGenReducedCosts " get_orig_costs(ctx::AbstractColGenContext) +@mustimplement "ColGenReducedCosts " get_subprob_var_orig_costs(ctx::AbstractColGenContext) """ Returns the coefficient matrix `A` of subproblem variables in the master to compute reduced cost `̄c = c - transpose(A) * π`. """ -@mustimplement "ColGenReducedCosts" get_coef_matrix(ctx::AbstractColGenContext) +@mustimplement "ColGenReducedCosts" get_subprob_var_coef_matrix(ctx::AbstractColGenContext) "Updates reduced costs of variables of a given subproblem." @mustimplement "ColGenReducedCosts" update_sp_vars_red_costs!(ctx::AbstractColGenContext, sp, red_costs) @@ -184,6 +184,8 @@ end function compute_dual_bound(ctx, phase, master_lp_obj_val, master_dbs) # TODO pure master variables are missing. + @show master_lp_obj_val + @show master_dbs return master_lp_obj_val + mapreduce(((id, val),) -> val, +, master_dbs) end @@ -242,15 +244,14 @@ function run_colgen_iteration!(context, phase, env) # stabcenter is master_dual_sol # return alpha * stab_center + (1 - alpha) * lp_dual_sol - # With stabilization, you solve several times the suproblem because you can have misprice # loop: # - solve all subproblems # - check if misprice # Compute reduced cost (generic operation) by you must support math operations. - c = get_orig_costs(context) - A = get_coef_matrix(context) + c = get_subprob_var_orig_costs(context) + A = get_subprob_var_coef_matrix(context) red_costs = c - transpose(A) * mast_dual_sol # Updates subproblems reduced costs. diff --git a/test/unit/ColGen/colgen_default.jl b/test/unit/ColGen/colgen_default.jl index b718c52fe..7e96eb0d8 100644 --- a/test/unit/ColGen/colgen_default.jl +++ b/test/unit/ColGen/colgen_default.jl @@ -2,11 +2,12 @@ form1() = """ master min - 3x1 + 2x2 + 5x3 + 4y1 + 3y2 + 5y3 + 3x1 + 2x2 + 5x3 + 4y1 + 3y2 + 5y3 + z s.t. - x1 + x2 + x3 + y1 + y2 + y3 >= 10 - x1 + 2x2 + y1 + 2y2 <= 100 - x1 + 3x3 + y1 + + 3y3 == 100 + x1 + x2 + x3 + y1 + y2 + y3 + 2z >= 10 + x1 + 2x2 + y1 + 2y2 + z <= 100 + x1 + 3x3 + y1 + + 3y3 == 100 + z <= 5 dw_sp min @@ -17,6 +18,9 @@ dw_sp integer representatives x1, x2, x3, y1, y2, y3 + + pure + z bounds x1 >= 0 @@ -25,6 +29,7 @@ dw_sp y1 >= 0 y2 >= 0 y3 >= 0 + z >= 0 """ function get_name_to_varids(form) @@ -46,33 +51,52 @@ end # Simple case with only subproblem representatives variables. function test_reduced_costs_calculation_helper() _, master, _, _, _ = reformfromstring(form1()) + @show master vids = get_name_to_varids(master) cids = get_name_to_constrids(master) helper = ClA.ReducedCostsCalculationHelper(master) - @test helper.c[vids["x1"]] == 3 - @test helper.c[vids["x2"]] == 2 - @test helper.c[vids["x3"]] == 5 - @test helper.c[vids["y1"]] == 4 - @test helper.c[vids["y2"]] == 3 - @test helper.c[vids["y3"]] == 5 - - @test helper.A[cids["c1"], vids["x1"]] == 1 - @test helper.A[cids["c1"], vids["x2"]] == 1 - @test helper.A[cids["c1"], vids["x3"]] == 1 - @test helper.A[cids["c1"], vids["y1"]] == 1 - @test helper.A[cids["c1"], vids["y2"]] == 1 - @test helper.A[cids["c1"], vids["y3"]] == 1 - - @test helper.A[cids["c2"], vids["x1"]] == 1 - @test helper.A[cids["c2"], vids["x2"]] == 2 - @test helper.A[cids["c2"], vids["y1"]] == 1 - @test helper.A[cids["c2"], vids["y2"]] == 2 - - @test helper.A[cids["c3"], vids["x1"]] == 1 - @test helper.A[cids["c3"], vids["x3"]] == 3 - @test helper.A[cids["c3"], vids["y1"]] == 1 - @test helper.A[cids["c3"], vids["y3"]] == 3 + @test helper.dw_subprob_c[vids["x1"]] == 3 + @test helper.dw_subprob_c[vids["x2"]] == 2 + @test helper.dw_subprob_c[vids["x3"]] == 5 + @test helper.dw_subprob_c[vids["y1"]] == 4 + @test helper.dw_subprob_c[vids["y2"]] == 3 + @test helper.dw_subprob_c[vids["y3"]] == 5 + @test helper.dw_subprob_c[vids["z"]] == 0 + + @test helper.dw_subprob_A[cids["c1"], vids["x1"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["x2"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["x3"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["y1"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["y2"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["y3"]] == 1 + @test helper.dw_subprob_A[cids["c1"], vids["z"]] == 0 # z is not in the subproblem. + + @test helper.dw_subprob_A[cids["c2"], vids["x1"]] == 1 + @test helper.dw_subprob_A[cids["c2"], vids["x2"]] == 2 + @test helper.dw_subprob_A[cids["c2"], vids["y1"]] == 1 + @test helper.dw_subprob_A[cids["c2"], vids["y2"]] == 2 + @test helper.dw_subprob_A[cids["c2"], vids["z"]] == 0 # z is not in the subproblem. + + @test helper.dw_subprob_A[cids["c3"], vids["x1"]] == 1 + @test helper.dw_subprob_A[cids["c3"], vids["x3"]] == 3 + @test helper.dw_subprob_A[cids["c3"], vids["y1"]] == 1 + @test helper.dw_subprob_A[cids["c3"], vids["y3"]] == 3 + @test helper.dw_subprob_A[cids["c3"], vids["z"]] == 0 # z is not in the subproblem. + + @test helper.master_c[vids["x1"]] == 0 # x1 is not in the master. + @test helper.master_c[vids["x2"]] == 0 # x2 is not in the master. + @test helper.master_c[vids["x3"]] == 0 # x3 is not in the master. + @test helper.master_c[vids["y1"]] == 0 # y1 is not in the master. + @test helper.master_c[vids["y2"]] == 0 # y2 is not in the master. + @test helper.master_c[vids["y3"]] == 0 # y3 is not in the master. + @test helper.master_c[vids["z"]] == 1 + + @test helper.master_A[cids["c1"], vids["x1"]] == 0 # x1 is not in the master. + @test helper.master_A[cids["c1"], vids["z"]] == 2 + @test helper.master_A[cids["c2"], vids["z"]] == 1 + @test helper.master_A[cids["c3"], vids["z"]] == 0 + @test helper.master_A[cids["c4"], vids["z"]] == 1 end register!(unit_tests, "colgen_default", test_reduced_costs_calculation_helper) @@ -486,7 +510,6 @@ function ColGen.optimize_master_lp_problem!(master, ctx::TestColGenIterationCont end dual_sol = ColGen.get_dual_sol(output) - @show dual_sol for (constr_id, constr) in ClMP.getconstrs(master) name = ClMP.getname(master, constr) if !haskey(ctx.master_lp_dual_sol, name) @@ -501,17 +524,13 @@ end ColGen.is_unbounded(ctx::TestColGenIterationContext) = ColGen.is_unbounded(ctx.context) ColGen.is_infeasible(ctx::TestColGenIterationContext) = ColGen.is_infeasible(ctx.context) ColGen.update_master_constrs_dual_vals!(ctx::TestColGenIterationContext, phase, reform, master_lp_dual_sol) = ColGen.update_master_constrs_dual_vals!(ctx.context, phase, reform, master_lp_dual_sol) -ColGen.get_orig_costs(ctx::TestColGenIterationContext) = ColGen.get_orig_costs(ctx.context) -ColGen.get_coef_matrix(ctx::TestColGenIterationContext) = ColGen.get_coef_matrix(ctx.context) +ColGen.get_subprob_var_orig_costs(ctx::TestColGenIterationContext) = ColGen.get_subprob_var_orig_costs(ctx.context) +ColGen.get_subprob_var_coef_matrix(ctx::TestColGenIterationContext) = ColGen.get_subprob_var_coef_matrix(ctx.context) function ColGen.update_sp_vars_red_costs!(ctx::TestColGenIterationContext, sp::Formulation{DwSp}, red_costs) - for i in 1:5 - println("\e[34m ***************** \e[00m") - end ColGen.update_sp_vars_red_costs!(ctx.context, sp, red_costs) for (_, var) in ClMP.getvars(sp) name = ClMP.getname(sp, var) - println(" ---- name = $(name) ---- expected : $(ctx.pricing_var_reduced_costs[name]) ---- actual : $(ClMP.getcurcost(sp, var)) --- cur_cost = $(ClMP.getcurcost(sp, var)))") @test ctx.pricing_var_reduced_costs[name] ≈ ClMP.getcurcost(sp, var) end return @@ -597,7 +616,7 @@ function test_colgen_iteration_min_gap() @test output.infeasible_subproblem == false @test output.unbounded_subproblem == false end -#register!(unit_tests, "colgen_default", test_colgen_iteration_min_gap) +register!(unit_tests, "colgen_default", test_colgen_iteration_min_gap) function test_colgen_iteration_max_gap() env, master, sps, reform = max_toy_gap() @@ -659,10 +678,8 @@ function test_colgen_iteration_max_gap() @test output.unbounded_master == false @test output.infeasible_subproblem == false @test output.unbounded_subproblem == false - - end -#register!(unit_tests, "colgen_default", test_colgen_iteration_max_gap) +register!(unit_tests, "colgen_default", test_colgen_iteration_max_gap) function test_colgen_iteration_pure_master_vars() env, master, sps, reform = toy_gap_with_penalties() @@ -681,13 +698,13 @@ function test_colgen_iteration_pure_master_vars() "c6" => 15.41666667, "c7" => 19.26666667, # fixed "c8" => -10.86666667, - "c10" => -22.55, + "c10" => -22.55, "c12" => -30.83333334 ) master_obj_val = 52.95 pricing_var_reduced_costs = Dict( - "x_11" => 0.26666666999999933, + "x_11" => - 0.26666666999999933, "x_12" => - 12.13333333, "x_13" => - 7.56666667, "x_14" => 0.0, diff --git a/test/unit/ColGen/colgen_iteration.jl b/test/unit/ColGen/colgen_iteration.jl index d2fca5309..5630cdd38 100644 --- a/test/unit/ColGen/colgen_iteration.jl +++ b/test/unit/ColGen/colgen_iteration.jl @@ -153,8 +153,8 @@ function ColGen.optimize_pricing_problem!(ctx::ColGenIterationTestContext, form, end # Reduced costs -ColGen.get_orig_costs(::ColGenIterationTestContext) = [7, 2, 1, 5, 3, 6, 8, 4, 2, 9] -ColGen.get_coef_matrix(::ColGenIterationTestContext) = [ +ColGen.get_subprob_var_orig_costs(::ColGenIterationTestContext) = [7, 2, 1, 5, 3, 6, 8, 4, 2, 9] +ColGen.get_subprob_var_coef_matrix(::ColGenIterationTestContext) = [ 1 1 1 1 0 0 0 0 0 0; 1 0 0 0 1 1 1 0 0 0; 0 1 0 0 1 0 0 1 1 0;