cut generators can also generate lazy constraints (initially for gurobi)

Former-commit-id: 25e31e4 [formerly 3e81e0f]
Former-commit-id: 54c66ef

Author: h-g-s

examples/queens-lazy.py

@@ -0,0 +1,61 @@
+"""Example of a solver to the n-queens problem:
+   n chess queens should be placed in a n x n
+   chess board so that no queen can attack another,
+   i.e., just one queen per line, column and diagonal.
+"""
+
+from sys import stdout
+from mip.model import Model, xsum
+from mip.constants import MAXIMIZE, BINARY
+from mip.callbacks import CutsGenerator
+
+
+class DiagonalCutGenerator(CutsGenerator):
+
+    def generate_cuts(self, model: Model):
+        def row(vname: str) -> str:
+            return int(vname.split('(')[1].split(',')[0].split(')')[0])
+
+        def col(vname: str) -> str:
+            return int(vname.split('(')[1].split(',')[1].split(')')[0])
+
+        x = {(row(v.name), col(v.name)): v for v in model.vars}
+        for p, k in enumerate(range(2 - n, n - 2 + 1)):
+            cut = xsum(x[i, j] for i in range(n) for j in range(n)
+                       if i - j == k) <= 1
+            if cut.violation > 0.001:
+                model.add_cut(cut)
+
+        for p, k in enumerate(range(3, n + n)):
+            cut = xsum(x[i, j] for i in range(n) for j in range(n)
+                       if i + j == k) <= 1
+            if cut.violation > 0.001:
+                model.add_cut(cut)
+
+
+# number of queens
+n = 8
+
+queens = Model('queens', MAXIMIZE)
+
+x = [[queens.add_var('x({},{})'.format(i, j), var_type=BINARY)
+      for j in range(n)] for i in range(n)]
+
+# one per row
+for i in range(n):
+    queens += xsum(x[i][j] for j in range(n)) == 1, 'row({})'.format(i)
+
+# one per column
+for j in range(n):
+    queens += xsum(x[i][j] for i in range(n)) == 1, 'col({})'.format(j)
+
+
+queens.cuts_generator = DiagonalCutGenerator()
+queens.cuts_generator.lazy_constraints = True
+queens.optimize()
+
+stdout.write('\n')
+for i, v in enumerate(queens.vars):
+    stdout.write('O ' if v.x >= 0.99 else '. ')
+    if i % n == n-1:
+        stdout.write('\n')

mip/callbacks.py

@@ -15,7 +15,7 @@ class CutsGenerator:
     """abstract class for implementing cut generators"""
 
     def __init__(self):
-        self.lazyConstraints = False
+        self.lazy_constraints = False
 
     def generate_cuts(self, model: Model):
         """Method called by the solver engine to generate cuts

mip/constants.py

@@ -2,7 +2,7 @@
 
 from enum import Enum
 
-VERSION = '1.3.10'
+VERSION = '1.3.11'
 
 # epsilon number (practical zero)
 EPS = 10e-6

mip/gurobi.py

@@ -20,11 +20,11 @@
     if 'GUROBI_HOME' in environ:
         if platform.lower().startswith('win'):
             libfile = glob(os.path.join(os.environ['GUROBI_HOME'],
-                                        'bin\gurobi[0-9][0-9].dll'))
+                                        'bin\\gurobi[0-9][0-9].dll'))
         else:
             libfile = glob(os.path.join(os.environ['GUROBI_HOME'],
                                         'lib/libgurobi[0-9][0-9].*'))
-        
+
         if libfile:
             lib_path = libfile[0]
 
@@ -249,6 +249,8 @@
 GRBgetstrattr = grblib.GRBgetstrattr
 GRBsetstrattr = grblib.GRBsetstrattr
 
+GRB_CB_MIPSOL    = 4
+GRB_CB_MIPNODE   = 5
 
 GRB_CB_PRE_COLDEL = 1000
 GRB_CB_PRE_ROWDEL = 1001
@@ -390,7 +392,6 @@ def add_var(self,
     def add_cut(self, lin_expr: LinExpr):
         # int GRBcbcut(void *cbdata, int cutlen, const int *cutind, const double *cutval, char cutsense, double cutrhs);
         # int GRBcbcut(void *cbdata, int cutlen, const int *cutind, const double *cutval, char cutsense, double cutrhs);
-        
 
         return
 
@@ -498,9 +499,12 @@ def callback(p_model: CData,
                                               obj_bound, obj_best))
                                     log.append((sec, (obj_bound, obj_best)))
 
-            # adding cuts
-            if where == 5:  # MIPNODE == 5
-                if self.model.cuts_generator:
+            # adding cuts or lazy constraints
+            if self.model.cuts_generator:
+                if where == GRB_CB_MIPNODE or \
+                    (where == GRB_CB_MIPSOL and
+                     hasattr(self.model.cuts_generator, 'lazy_constraints') and
+                     self.model.cuts_generator.lazy_constraints):
                     mgc = ModelGurobiCB(p_model, p_cbdata, where)
                     self.model.cuts_generator.generate_cuts(mgc)
 
@@ -538,6 +542,12 @@ def callback(p_model: CData,
            self.model.store_search_progress_log:
             GRBsetcallbackfunc(self._model, callback, ffi.NULL)
 
+        if (self.model.cuts_generator is not None and
+            hasattr(self.model.cuts_generator, 'lazy_constraints') and
+            self.model.cuts_generator.lazy_constraints) or \
+                self.model.lazy_constrs_generator is not None:
+            self.set_int_param("LazyConstraints", 1)
+
         if self.__threads >= 1:
             self.set_int_param("Threads", self.__threads)
 
@@ -1141,15 +1151,16 @@ def __init__(self, model: Model, grb_model: CData = ffi.NULL,
         self._obj_value = INF
         self._best_bound = INF
         self._status = OptimizationStatus.LOADED
+        self._where = where
 
         # pre-allocate temporary space to query names
         self.__name_space = ffi.new("char[{}]".format(MAX_NAME_SIZE))
         # in cut generation
         self.__name_spacec = ffi.new("char[{}]".format(MAX_NAME_SIZE))
 
         self.__relaxed = False
-        if where == 5:
-            gstatus = ffi.new('int *')
+        gstatus = ffi.new('int *')
+        if where == 5:  # GRB_CB_MIPNODE
             res = GRBcbget(cb_data, where, GRB_CB_MIPNODE_STATUS, gstatus)
             if res != 0:
                 raise Exception('Error getting status')
@@ -1168,6 +1179,21 @@ def __init__(self, model: Model, grb_model: CData = ffi.NULL,
                     raise Exception('Error getting fractional solution')
             else:
                 self._cb_sol = ffi.NULL
+        elif where == 4:  # GRB_CB_MIPSOL
+            self._status = OptimizationStatus.FEASIBLE
+            ires = ffi.new('int *')
+            st = GRBgetintattr(grb_model, 'NumVars'.encode('utf-8'), ires)
+            if st != 0:
+                raise Exception('Could not query number of variables in Gurobi \
+                                callback')
+            ncols = ires[0]
+
+            self._cb_sol = \
+                ffi.new('double[{}]'.format(ncols))
+            res = GRBcbget(cb_data, where, GRB_CB_MIPSOL_SOL, self._cb_sol)
+            if res != 0:
+                raise Exception('Error getting integer solution in gurobi \
+                                callback')
         else:
             self._cb_sol = ffi.NULL
 
@@ -1182,7 +1208,10 @@ def add_cut(self, cut: LinExpr):
         sense = cut.sense.encode("utf-8")
         rhs = -cut.const
 
-        GRBcbcut(self._cb_data, numnz, cind, cval, sense, rhs)
+        if self._where == GRB_CB_MIPNODE:
+            GRBcbcut(self._cb_data, numnz, cind, cval, sense, rhs)
+        elif self._where == GRB_CB_MIPSOL:
+            GRBcblazy(self._cb_data, numnz, cind, cval, sense, rhs)
 
     def get_status(self):
         return self._status

mip/model.py

@@ -363,6 +363,25 @@ def sense(self, value):
         """
         self.__sense = value
 
+    @property
+    def violation(self):
+        """Amount that current solution violates this constraint
+
+        If a solution is available, than this property indicates how much
+        the current solution violates this constraint.
+        """
+        lhs = sum(coef*var.x for (var, coef) in self.__expr.items())
+        rhs = -self.const
+        viol = 0.0
+        if self.sense == '=':
+            viol = abs(lhs-rhs)
+        elif self.sense == '<':
+            viol = max(lhs-rhs, 0.0)
+        elif self.sense == '>':
+            viol = max(rhs-lhs, 0.0)
+
+        return viol
+
 
 class ProgressLog:
     """Class to store the improvement of lower