I have an MILP problem and use CPLEX (Python interface). I am working on user heuristics for branching in the branch-and-bound procedure. With HSCallback I managed to get the information about the current subproblem (all integer variables are binary), in particular, what the variables that have been previously fixed/branched on are and how they've been fixed.
Now, I would like to do something with the current subproblem. Is there a way to obtain it as a model? Do I have to construct it each time in the callback from the original model by adding the constraints corresponding to the fixed variables? Or is there a more efficient way of doing it?
Below is a working example for a problem encoding Boolean satisfiability of CNF with 5 variables that has exactly one assignment (corresponding to the commented out line, uncommenting it makes the problem infeasible). I disable all possible heauristics and cuts to make sure that branch and bound has to do some work. Here I construct the subproblem by creating a copy of the original model and adding the constraints corresponding to the current branch.
import datetime
import docplex.mp.model as cpx
import cplex.callbacks as CPX_CB
from docplex.mp.relax_linear import LinearRelaxer
class MyHSCallback(CPX_CB.HSCallback):
def call(self, args, *kwargs):
print()
print(datetime.datetime.now(), "MyHSCallback working.............")
print("node ID", self.get_node_ID())
fixed = [(ind, int(self.get_lower_bounds(ind)))
for ind in range(self.var_count) if self.get_lower_bounds(ind) == self.get_upper_bounds(ind)]
print("Fixed variables", ", ".join([str(ind) + ":" + str(value) for ind, value in fixed]))
print("Feasible subproblem", self.solve())
# creating a copy of the subproblem
model_copy = self.original_model.copy()
for ind, value in fixed:
model_copy.add_constraint(model_copy.get_var_by_index(ind) == value)
# doing something with the current subproblem
lp_model = LinearRelaxer().linear_relaxation(mdl=model_copy)
sol = lp_model.solve()
if sol is not None:
print("Relaxation Values", sol.get_values([lp_model.get_var_by_index(ind) for ind in range(self.var_count)]))
model = cpx.Model(name="cnf")
vars = model.binary_var_list(5, name="var")
x = vars[0]
y = vars[1]
z = vars[2]
w = vars[3]
v = vars[4]
model.add_constraint(w + x + y >= 1)
model.add_constraint(v + w - y >= 0)
model.add_constraint(w - x + y + z >= 0)
model.add_constraint(v - x + y - z >= -1)
model.add_constraint(v - w + y + z >= 0)
model.add_constraint(v - w + x - z >= -1)
model.add_constraint(v - w + x - y + z >= -1)
model.add_constraint(- w - x - y >= -2)
model.add_constraint(-v + w - y + z >= -1)
model.add_constraint(-v + w + x - y - z >= -2)
model.add_constraint(-v + w - x - z >= -2)
model.add_constraint(-v - w + x >= -1)
model.add_constraint(-v - w - x + y >= -2)
hs_cb = model.register_callback(MyHSCallback)
hs_cb.var_count = len(vars)
hs_cb.original_model = model
model.parameters.preprocessing.symmetry = 0
model.parameters.preprocessing.presolve = 0
model.parameters.preprocessing.boundstrength = 0
model.parameters.preprocessing.coeffreduce = 0
model.parameters.mip.strategy.search = 1 #traditional
model.parameters.mip.strategy.probe = -1
model.parameters.mip.strategy.heuristicfreq = -1
model.parameters.mip.strategy.fpheur = -1
model.parameters.mip.limits.cutpasses = -1
model.parameters.mip.cuts.bqp = -1
model.parameters.mip.cuts.cliques = -1
model.parameters.mip.cuts.covers = -1
model.parameters.mip.cuts.disjunctive = -1
model.parameters.mip.cuts.flowcovers = -1
model.parameters.mip.cuts.gomory = -1
model.parameters.mip.cuts.gubcovers = -1
model.parameters.mip.cuts.implied = -1
model.parameters.mip.cuts.liftproj = -1
model.parameters.mip.cuts.localimplied = -1
model.parameters.mip.cuts.mcfcut = -1
model.parameters.mip.cuts.mircut = -1
model.parameters.mip.cuts.pathcut = -1
model.parameters.mip.cuts.rlt = -1
model.parameters.mip.cuts.zerohalfcut = -1
sol = model.solve()
print()
if sol is None:
print("Infeasible")
else:
model.print_solution()
values = sol.get_values(vars)
print(values)
```
