Starting default implementation of Colgen API

src/Algorithm/colgen.jl

@@ -217,6 +217,7 @@ function clear_before_colgen_iteration!(spinfo::SubprobInfo)
     return
 end
 
+# used by stabiization only
 function compute_subgradient_contribution(
     algo::ColumnGeneration, master::Formulation, puremastervars::Vector{Pair{VarId,Float64}},
     spinfos::Dict{FormId,SubprobInfo}
@@ -241,6 +242,7 @@ function compute_subgradient_contribution(
     return sparsevec(var_ids, var_vals)
 end
 
+# used in reduced_costs_of_solutions
 function compute_reduced_cost(
     stab_is_used, masterform::Formulation,
     spsol::PrimalSolution, lp_dual_sol::DualSolution

src/Algorithm/colgen/default.jl

@@ -1,12 +1,40 @@
 struct ColGenContext <: ColGen.AbstractColGenContext
-    # Information to solve the master
-    master_solve_alg
-    master_optimizer_id
+    reform::Reformulation
+    optim_sense
+    current_ip_primal_bound
 
-    # Memoization to compute reduced costs (this is a precompute)
-    redcost_mem
+    restr_master_solve_alg
+    restr_master_optimizer_id::Int
+
+    pricing_solve_alg
+
+    reduced_cost_helper::ReducedCostsCalculationHelper
+
+    # # Information to solve the master
+    # master_solve_alg
+    # master_optimizer_id
+
+    # # Memoization to compute reduced costs (this is a precompute)
+    # redcost_mem
+    function ColGenContext(reform, alg)
+        rch = ReducedCostsCalculationHelper(getmaster(reform))
+        return new(
+            reform, 
+            getobjsense(reform), 
+            0.0, 
+            alg.restr_master_solve_alg, 
+            alg.restr_master_optimizer_id,
+            alg.pricing_prob_solve_alg,
+            rch
+        )
+    end
 end
 
+ColGen.get_reform(ctx::ColGenContext) = ctx.reform
+ColGen.get_master(ctx::ColGenContext) = getmaster(ctx.reform)
+ColGen.get_pricing_subprobs(ctx::ColGenContext) = get_dw_pricing_sps(ctx.reform)
+
+
 ###############################################################################
 # Sequence of phases
 ###############################################################################
@@ -38,6 +66,8 @@ struct ColGenPhase3 <: ColGen.AbstractColGenPhase end
 ## Implementation of `initial_phase`.
 ColGen.initial_phase(::ColunaColGenPhaseIterator) = ColGenPhase3()
 
+colgen_mast_lp_sol_has_art_vars(ctx) = false
+
 ## Implementation of `next_phase`.
 function ColGen.next_phase(::ColunaColGenPhaseIterator, ::ColGenPhase1, ctx)
     # If master LP solution has no artificial vars, it means that the phase 1 has succeeded.
@@ -63,12 +93,7 @@ function ColGen.next_phase(::ColunaColGenPhaseIterator, ::ColGenPhase3, ctx)
     return nothing
 end
 
-## Methods used in the implementation and that we should mock in tests.
-function colgen_mast_lp_sol_has_art_vars(ctx::ColGenContext)
-
-end
-
-## Implementatation of `setup_reformulation!`
+# Implementatation of `setup_reformulation!`
 ## Phase 1 => non-artifical variables have cost equal to 0
 function ColGen.setup_reformulation!(reform, ::ColGenPhase1)
     master = getmaster(reform)
@@ -106,59 +131,151 @@ function ColGen.setup_reformulation!(reform, ::ColGenPhase3)
     return
 end
 
-######### Pricing strategy
-struct ClassicPricingStrategy <: ColGen.AbstractPricingStrategy 
-    subprobs::Dict{FormId, Formulation{DwSp}}
-end
 
-# function ColGen.get_pricing_strategy(ctx::ColGen.AbstractColGenContext, _)
-#     ClassicPricingStrategy(Dict(i => sp for (i, sp) in ColGen.get_pricing_subprobs(ctx)))
-# end
 
-# ColGen.pricing_strategy_iterate(ps::ClassicPricingStrategy) = iterate(ps.subprobs)
+# Master resolution
+
+"""
+    ColGenMasterResult{F,S}
 
+Contains the solution to the master LP.
+- `F` is the formulation type
+- `S` is the objective sense Type
+"""
+struct ColGenMasterResult{F,S}
+    result::OptimizationState{F,S}
+end
 
-######### Column generation
+# TODO: not type stable !!
+function ColGen.optimize_master_lp_problem!(master, ctx::ColGenContext, env)
+    rm_input = OptimizationState(master, ip_primal_bound=ctx.current_ip_primal_bound)
+    opt_state = run!(ctx.restr_master_solve_alg, env, master, rm_input, ctx.restr_master_optimizer_id)
+    return ColGenMasterResult(opt_state)
+end
 
-# Placeholder methods:  
-ColGen.before_colgen_iteration(::ColGenContext, _, _) = nothing
-ColGen.after_colgen_iteration(::ColGenContext, _, _, _) = nothing
+function ColGen.is_infeasible(master_res::ColGenMasterResult)
+    status = getterminationstatus(master_res.result)
+    return status == ClB.INFEASIBLE || status == ClB.INFEASIBLE_OR_UNBOUNDED
+end
 
-######### Column generation iteration
-function ColGen.optimize_master_lp_problem!(master, context, env)
-    println("\e[31m optimize master lp problem \e[00m")
-    input = OptimizationState(master, ip_primal_bound=0.0) # TODO : ip_primal_bound=get_ip_primal_bound(cg_optstate)
-    return run!(context.master_solve_alg, env, master, input, context.master_optimizer_id)
+function ColGen.is_unbounded(master_res::ColGenMasterResult)
+    status = getterminationstatus(master_res.result)
+    return status == ClB.DUAL_INFEASIBLE || status == ClB.INFEASIBLE_OR_UNBOUNDED
 end
 
-#get_primal_sol(mast_result)
+ColGen.get_primal_sol(master_res::ColGenMasterResult) = get_best_lp_primal_sol(master_res.result)
+ColGen.get_dual_sol(master_res::ColGenMasterResult) = get_best_lp_dual_sol(master_res.result)
+ColGen.get_obj_val(master_res::ColGenMasterResult) = get_lp_primal_bound(master_res.result)
+
+function ColGen.update_master_constrs_dual_vals!(ctx::ColGenContext, phase, reform, master_lp_dual_sol)
 
-function ColGen.check_primal_ip_feasibility(ctx, mast_lp_primal_sol)
-    println("\e[31m check primal ip feasibility \e[00m")
-    return !contains(mast_lp_primal_sol, varid -> isanArtificialDuty(getduty(varid))) &&
-        isinteger(proj_cols_on_rep(mast_lp_primal_sol, getmodel(mast_lp_primal_sol)))
 end
 
-#update_inc_primal_sol!
+function ColGen.check_primal_ip_feasibility(master_lp_primal_sol, phase, reform)
+
+end
 
-#get_dual_sol(mast_result)
+# 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
 
-function ColGen.update_master_constrs_dual_vals!(ctx, master, smooth_dual_sol)
-    println("\e[32m update_master_constrs_dual_vals \e[00m")
-    # Set all dual value of all constraints to 0.
-    for constr in Iterators.values(getconstrs(master))
-        setcurincval!(master, constr, 0.0)
+function ColGen.update_sp_vars_red_costs!(ctx::ColGenContext, sp::Formulation{DwSp}, red_costs)
+    for (var_id, _) in getvars(sp)
+        setcurcost!(sp, var_id, red_costs[var_id])
     end
-    # Update constraints that have non-zero dual values.
-    for (constr_id, val) in smooth_dual_sol
-        setcurincval!(master, constr_id, val)
+    return
+end
+
+# Columns insertion
+function ColGen.insert_columns!(reform, ctx::ColGenContext, phase, columns)
+    for column in columns
+        @show column
     end
+    return 1
 end
 
-function ColGen.update_sp_vars_red_costs!(ctx, sp, red_costs)
-    println("\e[34m update_sp_vars_red_costs \e[00m")
-    for (var_id, _) in getvars(sp)
-        setcurcost!(sp, var_id, red_costs[var_id])
+#############################################################################
+# Pricing strategy
+#############################################################################
+struct ClassicColGenPricingStrategy <: ColGen.AbstractPricingStrategy
+    subprobs::Dict{FormId, Formulation{DwSp}}
+end
+
+ColGen.get_pricing_strategy(ctx::ColGen.AbstractColGenContext, _) = ClassicColGenPricingStrategy(ColGen.get_pricing_subprobs(ctx))
+ColGen.pricing_strategy_iterate(ps::ClassicColGenPricingStrategy) = iterate(ps.subprobs)
+ColGen.pricing_strategy_iterate(ps::ClassicColGenPricingStrategy, state) = iterate(ps.subprobs, state)
+
+#############################################################################
+# Column generation
+#############################################################################
+function ColGen.compute_sp_init_db(ctx::ColGenContext, sp::Formulation{DwSp})
+    return ctx.optim_sense == MinSense ? -Inf : Inf
+end
+
+struct GeneratedColumn
+    column::PrimalSolution{Formulation{DwSp}}
+    red_cost::Float64
+end
+
+"""
+    A structure to store a collection of columns
+"""
+struct ColumnsSet
+    columns::Vector{GeneratedColumn}
+    ColumnsSet() = new(GeneratedColumn[])
+end
+Base.iterate(set::ColumnsSet) = iterate(set.columns)
+Base.iterate(set::ColumnsSet, state) = iterate(set.columns, state)
+
+function ColGen.set_of_columns(ctx::ColGenContext)
+    return ColumnsSet()
+end
+
+struct ColGenPricingResult{F,S}
+    result::OptimizationState{F,S}
+    columns::Vector{GeneratedColumn}
+end
+
+function ColGen.is_infeasible(pricing_res::ColGenPricingResult)
+    status = getterminationstatus(pricing_res.result)
+    return status == ClB.INFEASIBLE || status == ClB.INFEASIBLE_OR_UNBOUNDED
+end
+
+function ColGen.is_unbounded(pricing_res::ColGenPricingResult)
+    status = getterminationstatus(pricing_res.result)
+    return status == ClB.DUAL_INFEASIBLE || status == ClB.INFEASIBLE_OR_UNBOUNDED
+end
+
+ColGen.get_primal_sols(pricing_res) = pricing_res.columns
+ColGen.get_dual_bound(pricing_res) = get_ip_dual_bound(pricing_res.result)
+
+has_improving_red_cost(column::GeneratedColumn) = column.red_cost < 0
+# In our implementation of `push_in_set!`, we keep only columns that have improving reduced 
+# cost.
+function ColGen.push_in_set!(pool, column)
+    println("push_in_set!")
+    # We keep only columns that improve reduced cost
+    if has_improving_red_cost(column)
+        push!(pool.columns, column)
     end
     return
+end
+
+function ColGen.optimize_pricing_problem!(ctx::ColGenContext, sp::Formulation{DwSp}, env, master_dual_sol)
+    input = OptimizationState(sp)
+    opt_state = run!(ctx.pricing_solve_alg, env, sp, input) # master & master dual sol for non robust cuts
+
+    # Reduced cost of a column is composed of the cost of the subproblem variables
+    # and the contribution of the master convex combination constraints.
+    # TODO: find better name
+    # TODO; make sure duty_data is well defined.
+    lb_dual = 0.0 # master_dual_sol[sp.duty_data.lower_multiplicity_constr_id]
+    ub_dual = 0.0 # master_dual_sol[sp.duty_data.upper_multiplicity_constr_id]
+    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)
 end

src/Algorithm/colgen/pricing.jl

@@ -0,0 +1 @@
+

src/Algorithm/colgen/utils.jl

@@ -75,6 +75,11 @@ subproblem variables:
   representative 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_orig_coefmatrix
 """
 struct ReducedCostsCalculationHelper
     c::SparseVector{Float64,VarId}

src/Algorithm/colgenstabilization.jl

@@ -45,7 +45,6 @@ function ClB.restore_from_record!(
     return
 end
 
-
 # function ColGenStabilizationUnit(master::Formulation)
 #     return ColGenStabilizationUnit(
 #         0.5, 0.0, 0, DualBound(master), DualBound(master), nothing, nothing, nothing
@@ -101,7 +100,6 @@ function update_alpha_automatically!(
     unit::ColGenStabilizationUnit, nb_new_col::Int64, lp_dual_sol::DualSolution{M},  
     smooth_dual_sol::DualSolution{M}, h, subgradient_contribution
 ) where {M}    
-
     master = getmodel(lp_dual_sol)
 
     # first we calculate the in-sep direction

src/ColGen/interface.jl

@@ -140,7 +140,7 @@ information at two different places.
 Returns a primal solution expressed in the original problem variables if the current master
 LP solution is integer feasible; `nothing` otherwise.
 """
-@mustimplement "ColGenMaster" check_primal_ip_feasibility(phase, mast_lp_primal_sol, reform)
+@mustimplement "ColGenMaster" check_primal_ip_feasibility(mast_lp_primal_sol, phase, reform)
 
 ############################################################################################
 # Reduced costs calculation.
@@ -173,14 +173,17 @@ if something unexpected happens.
 
 
 function check_master_termination_status(mast_result)
-    # TODO
+    if !is_infeasible(mast_result) && !is_unbounded(mast_result)
+        @assert !isnothing(get_dual_sol(mast_result))
+    end
 end
 
 function check_pricing_termination_status(pricing_result)
     # TODO
 end
 
 function compute_dual_bound(ctx, phase, master_lp_obj_val, master_dbs)
+    # TODO pure master variables are missing.
     return master_lp_obj_val + mapreduce(((id, val),) -> val, +, master_dbs)
 end
 
@@ -224,6 +227,7 @@ function run_colgen_iteration!(context, phase, env)
 
     mast_dual_sol = get_dual_sol(mast_result)
     if isnothing(mast_dual_sol)
+        error("Cannot continue")
         # error or stop? (depends on the context)
     end
 
@@ -233,6 +237,16 @@ function run_colgen_iteration!(context, phase, env)
     update_master_constrs_dual_vals!(context, phase, get_reform(context), mast_dual_sol)
 
     # Stabilization
+    # initialize stabilisation for the iteration
+    # update_stab_after_rm_solve! 
+    # 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)
@@ -262,7 +276,7 @@ function run_colgen_iteration!(context, phase, env)
 
     while !isnothing(sp_to_solve_it)
         (sp_id, sp_to_solve), state = sp_to_solve_it
-        pricing_result = optimize_pricing_problem!(context, sp_to_solve)
+        pricing_result = optimize_pricing_problem!(context, sp_to_solve, env, mast_dual_sol)
 
         # Iteration continues only if the pricing solution is not infeasible nor unbounded.
         if is_infeasible(pricing_result)
@@ -294,9 +308,17 @@ function run_colgen_iteration!(context, phase, env)
     nb_cols_inserted = insert_columns!(get_reform(context), context, phase, generated_columns)
 
     master_lp_obj_val = get_obj_val(mast_result)
-    db = compute_dual_bound(context, phase, master_lp_obj_val, sps_db)
+
+    # compute valid dual bound using the dual bounds returned by the user (cf pricing result).
+    valid_db = compute_dual_bound(context, phase, master_lp_obj_val, sps_db)
+
+    pseudo_db = 0 # same but using primal bound of the pricing result.
+    # pseudo_db used only in the stabilization (update_stability_center!)
+
+    # update_stab_after_gencols!
+
     # check gap
 
-    return ColGenIterationOutput(master_lp_obj_val, db, nb_cols_inserted, false, false, false, false)
+    return ColGenIterationOutput(master_lp_obj_val, valid_db, nb_cols_inserted, false, false, false, false)
 end