Custom Variables and Cuts
Coluna allows users to attach custom data to variables and constraints. This data is useful to store information about the variables or constraints in a custom format much easier to process than extracted information from the formulation (coefficient matrix, bounds, costs, and right-hand side).
In this example, we will show how to attach custom data to variables and constraints and use them to separate non-robust cuts. We will use the Bin Packing problem as an example.
Let us consider a Bin Packing problem with only 3 items such that any pair of items fits into one bin but the 3 items do not. The objective function is to minimize the number of bins being used. Pricing is done by inspection over the 6 combinations of items (3 pairs and 3 singletons). The master LP solution has 1.5 bins at the root node, each 0.5 corresponding to a bin with one of the possible pairs of items.
In this example, we will show you how to use non-robust cuts to improve the master LP solution at the root node. Obviously, Coluna is able to solve this instance by branching on the number of bins but the limit one on the number of nodes prevents it to be solved without cuts.
We define the dependencies:
using JuMP, BlockDecomposition, Coluna, GLPK;We define the solver.
coluna = JuMP.optimizer_with_attributes(
Coluna.Optimizer,
"default_optimizer" => GLPK.Optimizer,
"params" => Coluna.Params(
solver = Coluna.Algorithm.TreeSearchAlgorithm(
conqueralg = Coluna.Algorithm.ColCutGenConquer(
colgen = Coluna.Algorithm.ColumnGeneration(
pricing_prob_solve_alg = Coluna.Algorithm.SolveIpForm(
optimizer_id = 1
))
),
maxnumnodes = 1 # we only treat the root node.
)
)
);Let's define the model. Let's $B$ the set of bins and $I$ the set of items. We introduce variable $y_b$ that is equal to 1 if a bin $b$ is used and 0 otherwise. We introduce variable $x_{b,i}$ that is equal to 1 if item $i$ is put in a bin $b$ and 0 otherwise.
model = BlockModel(coluna);We must assign three items:
I = [1, 2, 3];And we have three bins:
B = [1, 2, 3];Each bin is defining a subproblem, we declare our axis:
@axis(axis, collect(B));We declare subproblem variables y[b]:
@variable(model, y[b in axis], Bin);And x[b,i]:
@variable(model, x[b in axis, i in I], Bin);Each item must be assigned to one bin:
@constraint(model, sp[i in I], sum(x[b,i] for b in axis) == 1);We minimize the number of bins and we declare the decomposition:
@objective(model, Min, sum(y[b] for b in axis))
@dantzig_wolfe_decomposition(model, dec, axis);Custom data for non-robust cuts
As said previously, at the end of the column generation at the root node, the master LP solution has 1.5 bins. It corresponds to three bins, each of them used 0.5 times containing one pair (1,2), (1, 3), or (2, 3) of items. We are going to introduce the following non-robust cut to make the master LP solution integral:
\[\sum\limits_{s \in S~if~length(s) \geq 2} λ_s \leq 1\]
where :
- $S$ is the set of possible bin assignments generated by the pricing problem.
- $length(s)$ the number of items in bin assignment $s \in S$.
This cut means that we cannot have more than one bin with at least two items.
But the problem is that the cut is expressed over the master column and we don't have access to these variables from the JuMP model. To address this problem, Coluna offers a way to compute the coefficient of a column in a constraint by implementing the following method:
Coluna.MathProg.computecoeff — Functioncomputecoeff(var_custom_data, constr_custom_data) -> Float64Dispatches on the type of custom data attached to the variable and the constraint to compute the coefficient of the variable in the constraint.
We therefore need to attach custom data to the master columns and the non-robust cut to use the method compute_coeff.
For every subproblem solution $s$, we define custom data with the number of items in the bin.
struct MyCustomVarData <: BlockDecomposition.AbstractCustomVarData
nb_items::Int
end
BlockDecomposition.customvars!(model, MyCustomVarData);We define custom data for the cut that will contain the minimum number of items in a bin that can be used. The value will be 2 in this example.
struct MyCustomCutData <: BlockDecomposition.AbstractCustomConstrData
min_items::Int
end
BlockDecomposition.customconstrs!(model, MyCustomCutData);We implement the computecoeff method for the custom data we defined.
function Coluna.MathProg.computecoeff(
var_custom_data::MyCustomVarData, constr_custom_data::MyCustomCutData
)
return (var_custom_data.nb_items >= constr_custom_data.min_items) ? 1.0 : 0.0
endPricing callback
We define the pricing callback that will generate the bin with best-reduced cost. Be careful, when using non-robust cuts, you must take into account the contribution of the non-robust cuts to the reduced cost of your solution.
function my_pricing_callback(cbdata)
# Get the reduced costs of the original variables.
I = [1, 2, 3]
b = BlockDecomposition.callback_spid(cbdata, model)
rc_y = BlockDecomposition.callback_reduced_cost(cbdata, y[b])
rc_x = [BlockDecomposition.callback_reduced_cost(cbdata, x[b, i]) for i in I]
# Get the dual values of the custom cuts (to calculate contributions of
# non-robust cuts to the cost of the solution).
custduals = Tuple{Int, Float64}[]
for (_, constr) in Coluna.MathProg.getconstrs(cbdata.form.parent_formulation)
if typeof(constr.custom_data) == MyCustomCutData
push!(custduals, (
constr.custom_data.min_items,
Coluna.MathProg.getcurincval(cbdata.form.parent_formulation, constr)
))
end
end
# Pricing by inspection.
sols = [[1], [2], [3], [1, 2], [1, 3], [2, 3]]
best_s = Int[]
best_rc = Inf
for s in sols
rc_s = rc_y + sum(rc_x[i] for i in s) # reduced cost of the subproblem variables
if !isempty(custduals)
# contribution of the non-robust cuts
rc_s -= sum((length(s) >= minits) ? dual : 0.0 for (minits, dual) in custduals)
end
if rc_s < best_rc
best_rc = rc_s
best_s = s
end
end
@show best_s
# build the best one and submit
solcost = best_rc
solvars = JuMP.VariableRef[]
solvarvals = Float64[]
for i in best_s
push!(solvars, x[b, i])
push!(solvarvals, 1.0)
end
push!(solvars, y[b])
push!(solvarvals, 1.0)
# submit the solution
MOI.submit(
model, BlockDecomposition.PricingSolution(cbdata),
solcost,
solvars,
solvarvals,
MyCustomVarData(length(best_s)) # attach a custom data to the column
)
MOI.submit(model, BlockDecomposition.PricingDualBound(cbdata), solcost)
return
endmy_pricing_callback (generic function with 1 method)The pricing callback is done, we define it as the solver of our pricing problem.
subproblems = BlockDecomposition.getsubproblems(dec)
BlockDecomposition.specify!.(
subproblems,
solver = my_pricing_callback
);Non-robust cut separation callback.
We now define the cut separation callback for our non-robust cut. This is the same callback as the one used for robust cuts. There is just one slight difference when you submit the non-robust cut. Since cuts are expressed over the master variables and these variables are inaccessible from the JuMP model, you'll submit a constraint with an empty left-hand side and you'll leave Coluna populate the left-hand side with the values returned by Coluna.MathProg.computecoeff.
So let's define the callback. Basically, if the solution uses more than one bin with two items, The cut is added to the model.
function custom_cut_sep(cbdata)
# Compute the constraint violation by iterating over the master solution.
viol = -1.0
for (varid, varval) in cbdata.orig_sol
var = Coluna.MathProg.getvar(cbdata.form, varid)
if !isnothing(var.custom_data)
if var.custom_data.nb_items >= 2
viol += varval
end
end
end
# Add the cut (at most one variable with 2 or more of the 3 items) if violated.
if viol > 0.001
MOI.submit(
model, MOI.UserCut(cbdata),
JuMP.ScalarConstraint(
JuMP.AffExpr(0.0), # We cannot express the left-hand side so we push 0.
MOI.LessThan(1.0)
),
MyCustomCutData(2) # Cut custom data.
)
end
return
end
MOI.set(model, MOI.UserCutCallback(), custom_cut_sep)
JuMP.optimize!(model)Coluna
Version 0.8.1 | https://github.com/atoptima/Coluna.jl
***************************************************************************************
**** B&B tree root node
**** Local DB = -Inf, global bounds: [ -Inf , Inf ], time = 0.00 sec.
***************************************************************************************
best_s = [1, 2]
best_s = [1, 2]
best_s = [1, 2]
<st= 1> <it= 1> <et= 0.79> <mst= 0.01> <sp= 0.00> <cols= 3> <al= 0.00> <DB=-29997.0000> <mlp=60000.0000> <PB=Inf>
best_s = [1, 3]
best_s = [1, 3]
best_s = [1, 3]
<st= 1> <it= 2> <et= 0.81> <mst= 0.02> <sp= 0.00> <cols= 3> <al= 0.00> <DB=-29999.0000> <mlp=30001.0000> <PB=Inf>
best_s = [2, 3]
best_s = [2, 3]
best_s = [2, 3]
<st= 1> <it= 3> <et= 0.81> <mst= 0.00> <sp= 0.00> <cols= 3> <al= 0.00> <DB=-10001.0000> <mlp=20002.0000> <PB=Inf>
best_s = [1]
best_s = [1]
best_s = [1]
<st= 1> <it= 4> <et= 0.81> <mst= 0.00> <sp= 0.00> <cols= 3> <al= 0.00> <DB= 3.0000> <mlp=15001.5000> <PB=Inf>
best_s = [3]
best_s = [3]
best_s = [3]
<st= 1> <it= 5> <et= 0.81> <mst= 0.00> <sp= 0.00> <cols= 3> <al= 0.00> <DB=-19995.0000> <mlp=10002.0000> <PB=Inf>
best_s = [2]
best_s = [2]
best_s = [2]
<st= 1> <it= 6> <et= 0.81> <mst= 0.00> <sp= 0.00> <cols= 3> <al= 0.00> <DB=-19995.0000> <mlp=10002.0000> <PB=Inf>
Cut separation callback adds 0 new essential cuts and 0 new facultative cuts.
[ Info: Improving primal solution with value 3.0 is found during column generation
best_s = [1]
best_s = [1]
best_s = [1]
<st= 1> <it= 7> <et= 0.81> <mst= 0.00> <sp= 0.00> <cols= 0> <al= 0.00> <DB= 3.0000> <mlp= 3.0000> <PB=3.0000>
───────────────────────────────────────────────────────────────────────────────
Time Allocations
─────────────────────── ────────────────────────
Tot / % measured: 3.78s / 21.5% 2.29GiB / 3.4%
Section ncalls time %tot avg alloc %tot avg
───────────────────────────────────────────────────────────────────────────────
Coluna 1 814ms 100.0% 814ms 80.8MiB 100.0% 80.8MiB
Pricing callback 21 275ms 33.8% 13.1ms 27.0MiB 33.4% 1.29MiB
SolveLpForm 7 25.7ms 3.2% 3.67ms 1.68MiB 2.1% 246KiB
───────────────────────────────────────────────────────────────────────────────
[ Info: Terminated
[ Info: Primal bound: 3.0
[ Info: Dual bound: 3.0We see on the output that the algorithm has converged a first time before a cut is added. Coluna then starts a new iteration taking into account the cut. We notice here an improvement of the value of the dual bound: before the cut, we converge towards 1.5. After the cut, we reach 2.0.
This page was generated using Literate.jl.