XCSP3 competition

cpmpy/expressions/core.py

@@ -829,7 +829,7 @@ def _wsum_should(arg):
      all substractions are transformed into less readable wsums)
     """
     return isinstance(arg, Operator) and \
-           (arg.name == 'wsum' or \
+           (arg.name == 'sum' or arg.name == 'wsum' or \
             (arg.name == 'mul' and len(arg.args) == 2 and \
              any(is_num(a) for a in arg.args)
             ) )

cpmpy/expressions/globalconstraints.py

@@ -194,7 +194,18 @@ def __init__(self, *args):
     def decompose(self):
         """Returns the decomposition
         """
-        return [var1 != var2 for var1, var2 in all_pairs(self.args)], []
+        if False:
+            return [var1 != var2 for var1, var2 in all_pairs(self.args)], []
+        else:
+            # DO NOT COMMIT IN MAINLINE yet
+            # switch to ILP friendly decomposition
+            cons = []
+            lbs, ubs = get_bounds(self.args)
+            M = [cp.MapDomain(var) for var in self.args if isinstance(var, Expression)]
+            for val in range(min(lbs), max(ubs)+1):
+                # each value can be taken at most once (not necessarily exactly once)
+                cons.append(cp.sum(x == val for x in self.args) <= 1)
+            return cons, M
 
     def value(self):
         return len(set(argvals(self.args))) == len(self.args)
@@ -210,8 +221,20 @@ def __init__(self, arr, n):
         super().__init__("alldifferent_except_n", [flatarr, n])
 
     def decompose(self):
-        # equivalent to (var1 == n) | (var2 == n) | (var1 != var2)
-        return [(var1 == var2).implies(cp.any(var1 == a for a in self.args[1])) for var1, var2 in all_pairs(self.args[0])], []
+        if False:
+            # equivalent to (var1 == n) | (var2 == n) | (var1 != var2)
+            return [(var1 == var2).implies(cp.any(var1 == a for a in self.args[1])) for var1, var2 in all_pairs(self.args[0])], []
+        else:
+            # DO NOT COMMIT IN MAINLINE yet
+            # switch to ILP friendly decomposition
+            cons = []
+            arr, n = self.args
+            lbs, ubs = get_bounds(arr)
+            for val in range(min(lbs), max(ubs)+1):
+                if val != n:
+                    # each value can be taken at most once (not necessarily exactly once)
+                    cons.append(cp.sum(x == val for x in arr) <= 1)
+            return cons, [MapDomain(x) for x in arr]
 
     def value(self):
         vals = [argval(a) for a in self.args[0] if argval(a) not in argvals(self.args[1])]
@@ -246,8 +269,19 @@ def __init__(self, *args):
     def decompose(self):
         """Returns the decomposition
         """
-        # arg0 == arg1, arg1 == arg2, arg2 == arg3... no need to post n^2 equalities
-        return [var1 == var2 for var1, var2 in zip(self.args[:-1], self.args[1:])], []
+        # Not sure we need the Boolean version here...
+        if False:
+            # arg0 == arg1, arg1 == arg2, arg2 == arg3... no need to post n^2 equalities
+            return [var1 == var2 for var1, var2 in zip(self.args[:-1], self.args[1:])], []
+        else:
+            # DO NOT COMMIT IN MAINLINE yet
+            # switch to ILP friendly decomposition
+            cons = []
+            lbs, ubs = get_bounds(self.args)
+            for val in range(min(lbs), max(ubs)+1):
+                # each value can be taken at most once (not necessarily exactly once)
+                cons.extend((x1 == val) == (x2 == val) for x1,x2 in all_pairs(self.args))
+            return cons, [MapDomain(x) for x in self.args]
 
     def value(self):
         return len(set(argvals(self.args))) == 1
@@ -265,8 +299,21 @@ def __init__(self, arr, n):
         super().__init__("allequal_except_n", [flatarr, n])
 
     def decompose(self):
-        return [(cp.any(var1 == a for a in self.args[1]) | (var1 == var2) | cp.any(var2 == a for a in self.args[1]))
-                for var1, var2 in all_pairs(self.args[0])], []
+        # Not sure we need the Boolean version here...
+        if False:
+            return [(cp.any(var1 == a for a in self.args[1]) | (var1 == var2) | cp.any(var2 == a for a in self.args[1]))
+                    for var1, var2 in all_pairs(self.args[0])], []
+        else:
+            # DO NOT COMMIT IN MAINLINE yet
+            # switch to ILP friendly decomposition
+            cons = []
+            arr, n = self.args
+            lbs, ubs = get_bounds(arr)
+            for val in range(min(lbs), max(ubs)+1):
+                if val != n:
+                    # each value can be taken at most once (not necessarily exactly once)
+                    cons.append((x1 == val) == (x2 == val) for x1,x2 in all_pairs(arr))
+            return cons, [MapDomain(x) for x in arr]
 
     def value(self):
         vals = [argval(a) for a in self.args[0] if argval(a) not in argvals(self.args[1])]
@@ -303,6 +350,7 @@ def decompose(self):
                 MiniZinc has slightly different one:
                 https://github.com/MiniZinc/libminizinc/blob/master/share/minizinc/std/fzn_circuit.mzn
         """
+        # XXX Right, so for ILP solvers, we would ideally do lazy subcircuit elimination (otherwise MTZ)...
         succ = cpm_array(self.args)
         n = len(succ)
         order = intvar(0,n-1, shape=n)
@@ -386,6 +434,7 @@ def decompose(self):
 
             lb, ub = get_bounds(x)
             if lb >= 0 and ub < len(rev): # safe, index is within bounds
+                # XXX Both rev and x are variables... is there an ILP friendly decomp?
                 constraining.append(rev[x] == i)
             else: # partial! need safening here
                 is_defined, total_expr, toplevel = cp.transformations.safening._safen_range(rev[x], (0, len(rev)-1), 1)
@@ -620,7 +669,10 @@ def decompose(self):
         if expressions:
             return [cp.any(expr == a for a in arr)], defining
         else:
-            return [expr != val for val in range(lb, ub + 1) if val not in arr], defining
+            # XXX do we properly capture that x!=v with b:=x==v should be ~b?
+            # Can't do MapDomain, would need to check if `expr` is a variable...
+            # TODO is this even efficient when there are many gaps? (similar to NotInDomain)
+            return [expr != val for val in range(lb, ub + 1) if val not in arr] + [expr >= min(arr), expr <= max(arr)], defining
 
 
     def value(self):
@@ -629,6 +681,62 @@ def value(self):
     def __repr__(self):
         return "{} in {}".format(self.args[0], self.args[1])
 
+class MapDomain(GlobalConstraint):
+    """
+        Maps an integer decision variable to an array of Boolean variables,
+        one for each value in the domain, e.g. |ub+1 - lb| Boolean variables.
+
+        No Boolean variables are returned or accessible;
+        after declaring this constraint, just
+        use 'ivar == v' and it will be replaced by the correct Boolean ?in flatten()?.
+
+        Note: the `decompose()` takes an optional `csemap` and `is_supported` argument!
+
+        TODO: what with solvers (e.g. SMT) that do not need flatten?
+        (they also don't need us to do CSE, so might be fine?)
+
+        TODO: multiple (e.g. by global constraints) MapDomain's can be declared
+        for the same integer variable...
+    """
+    def __init__(self, ivar):
+        super().__init__("mapdomain", [ivar])
+
+    def decompose(self, is_supported=False, csemap=None):
+        """
+            is_supported: if True, only the CSE map is filled
+                (yes, ortools declares a global for this, called... MapDomain )
+            csemap: if given, will populate the csemap with the Boolean variables
+        """
+        ivar = self.args[0]
+        lb, ub = get_bounds(ivar)
+
+        bvs = cp.boolvar(shape=(ub+1-lb,), name=f"B#{ivar.name}")
+        all_in_csemap = True
+        if csemap is not None:
+            for i,v in enumerate(range(lb, ub+1)):
+                expr = (ivar == v)
+                if expr in csemap:
+                    bvs[i] = csemap[expr]  # overwrite with pre-created one
+                else:
+                    all_in_csemap = False
+                    csemap[expr] = bvs[i]
+
+        if is_supported:
+            # TRICKY HACK to use 2nd argument as all_in_csemap...
+            return [], all_in_csemap
+
+        # ILP friendly decomposition
+        # TODO: if the 'ivar' is eliminated from the model, no need for 2nd constraint...
+        # TRICKY HACK to use 2nd argument as all_in_csemap...
+        return [cp.sum(bvs) == 1,
+                cp.sum(bvs[i]*v for i,v in enumerate(range(lb, ub+1))) == ivar], all_in_csemap
+
+    def value(self):
+        # not much to say...
+        return True
+
+    def __repr__(self):
+        return f"MapDomain({self.args[0]})"
 
 class Xor(GlobalConstraint):
     """
@@ -648,6 +756,7 @@ def __init__(self, arg_list):
 
     def decompose(self):
         # there are multiple decompositions possible, Recursively using sum allows it to be efficient for all solvers.
+        # and ILP friendly...
         decomp = [sum(self.args[:2]) == 1]
         if len(self.args) > 2:
             decomp = Xor([decomp,self.args[2:]]).decompose()[0]

cpmpy/expressions/globalfunctions.py

@@ -326,6 +326,16 @@ def decompose_comparison(self, cmp_op, cmp_rhs):
         arr, val = self.args
         return [eval_comparison(cmp_op, Operator('sum',[ai==val for ai in arr]), cmp_rhs)], []
 
+    def decompose_numerical(self):
+        """
+            Return a numerical expression to replace the array loopup with in the expression tree
+            XXX DOES NOT work automatically, needs hacking into decompose_global.py
+        """
+        arr, val = self.args
+        expr = cp.sum(ai == val for ai in arr)
+        return expr, []# [cp.expressions.globalconstraints.MapDomain(ai) for ai in arr] <- somehow breaks the constraint in the setting "count(..., ...) not in [...]""
+
+
     def value(self):
         arr, val = self.args
         val = argval(val)

cpmpy/solvers/choco.py

@@ -346,7 +346,7 @@ def objective(self, expr, minimize):
 
         # make objective function non-nested
         obj_var = intvar(*get_bounds(expr))
-        self += obj_var == expr
+        self.add(obj_var == expr, internal=True)
 
         self.obj = obj_var
         self.minimize_obj = minimize  # Choco has as default to maximize
@@ -406,7 +406,7 @@ def transform(self, cpm_expr):
 
         return cpm_cons
 
-    def add(self, cpm_expr):
+    def add(self, cpm_expr, internal:bool=False):
         """
             Eagerly add a constraint to the underlying solver.
 
@@ -425,7 +425,8 @@ def add(self, cpm_expr):
         :return: self
         """
         # add new user vars to the set
-        get_variables(cpm_expr, collect=self.user_vars)
+        if not internal:
+            get_variables(cpm_expr, collect=self.user_vars)
         # ensure all vars are known to solver
 
         # transform and post the constraints
@@ -596,7 +597,7 @@ def _get_constraint(self, cpm_expr):
                 table = table.astype(float) # nan's require float dtype
                 # Choco requires a wildcard value not present in dom of args,
                 # take value lower than anything else
-                chc_star = min(np.nanmin(table), *get_bounds(array)[0]) -1
+                chc_star = int(min(np.nanmin(table), *get_bounds(array)[0]) -1) # should be an int
                 chc_table = np.nan_to_num(table, nan=chc_star).astype(int).tolist()
                 return self.chc_model.table(self.solver_vars(array), chc_table, universal_value=chc_star, algo="STR2+")
             elif cpm_expr.name == "regular":

cpmpy/solvers/cpo.py

@@ -42,7 +42,6 @@
         CPM_cpo
 """
 
-import shutil
 import time
 import warnings
 
@@ -346,6 +345,8 @@ def objective(self, expr, minimize=True):
 
                 technical side note: any constraints created during conversion of the objective are permanently posted to the solver
         """
+        get_variables(expr, collect=self.user_vars)
+
         dom = self.get_docp().modeler
         if self.has_objective():
             self.cpo_model.remove(self.cpo_model.get_objective_expression())

cpmpy/solvers/exact.py

@@ -216,6 +216,11 @@ def solve(self, time_limit=None, assumptions=None, **kwargs):
         self.cpm_status.runtime = end - start
 
         self.objective_value_ = None
+
+        self._fillVars()
+        if self.has_objective():
+            self.objective_value_ = self.objective_.value()
+
         # translate exit status
         #   see 'toOptimum' documentation:
         #   https://gitlab.com/nonfiction-software/exact/-/blob/main/src/interface/IntProg.cpp#L877
@@ -235,19 +240,13 @@ def solve(self, time_limit=None, assumptions=None, **kwargs):
         elif my_status == "INCONSISTENT": # found inconsistency over assumptions
             self.cpm_status.exitstatus = ExitStatus.UNSATISFIABLE
         elif my_status == "TIMEOUT": # found timeout
-            if self.xct_solver.hasSolution(): # found a (sub-)optimal solution
+            if self.objective_value_ is not None: # found a (sub-)optimal solution
                 self.cpm_status.exitstatus = ExitStatus.FEASIBLE
             else: # no solution found
                 self.cpm_status.exitstatus = ExitStatus.UNKNOWN
         else:
             raise NotImplementedError(my_status)  # a new status type was introduced, please report on github
 
-        self._fillVars()
-        if self.has_objective():
-            if self.objective_is_min_:
-                self.objective_value_ = obj_val
-            else: # maximize, so actually negative value
-                self.objective_value_ = -obj_val
 
         # True/False depending on self.cpm_status
         return self._solve_return(self.cpm_status)
@@ -423,10 +422,10 @@ def objective(self, expr, minimize):
         self.objective_is_min_ = minimize
 
         # make objective function non-nested and with positive BoolVars only
+        get_variables(expr, collect=self.user_vars) # add objvars to vars
         (flat_obj, flat_cons) = flatten_objective(expr)
         flat_obj = only_positive_bv_wsum(flat_obj)  # remove negboolviews
-        self.user_vars.update(get_variables(flat_obj))  # add objvars to vars
-        self += flat_cons  # add potentially created constraints
+        self.add(flat_cons, internal=True) # add potentially created constraints
 
         # make objective function or variable and post
         xct_cfvars,xct_rhs = self._make_numexpr(flat_obj,0)
@@ -517,7 +516,7 @@ def _add_xct_reif_right(self, head, sign, xct_cfvars, xct_rhs):
     def is_multiplication(cpm_expr): # helper function
         return isinstance(cpm_expr, Operator) and cpm_expr.name == 'mul'
 
-    def add(self, cpm_expr_orig):
+    def add(self, cpm_expr_orig, internal:bool=False):
         """
             Eagerly add a constraint to the underlying solver.
 
@@ -537,7 +536,8 @@ def add(self, cpm_expr_orig):
         """
 
         # add new user vars to the set
-        get_variables(cpm_expr_orig, collect=self.user_vars)
+        if not internal:
+            get_variables(cpm_expr_orig, collect=self.user_vars)
 
         # transform and post the constraints
         for cpm_expr in self.transform(cpm_expr_orig):
@@ -552,7 +552,7 @@ def add(self, cpm_expr_orig):
                         assert pkg_resources.require("exact>=2.1.0"), f"Multiplication constraint {cpm_expr} " \
                                                                       f"only supported by Exact version 2.1.0 and above"
                         if is_num(rhs): # make dummy var
-                            rhs = intvar(rhs, rhs)
+                            rhs = cp.intvar(rhs, rhs)
                         xct_rhs = self.solver_var(rhs)
                         assert all(isinstance(v, _IntVarImpl) for v in lhs.args), "constant * var should be " \
                                                                                   "rewritten by linearize"

cpmpy/solvers/gcs.py

@@ -364,11 +364,11 @@ def objective(self, expr, minimize=True):
         """
         # make objective function non-nested
         (flat_obj, flat_cons) = flatten_objective(expr)
-        self += flat_cons # add potentially created constraints
+        self.add(flat_cons, internal=True) # add potentially created constraints
         self.user_vars.update(get_variables(flat_obj)) # add objvars to vars
 
         (obj, obj_cons) = get_or_make_var(flat_obj, csemap=self._csemap)
-        self += obj_cons
+        self.add(obj_cons, internal=True)
 
         self.objective_var = obj
 
@@ -471,7 +471,7 @@ def verify(self, name=None, location=".", time_limit=None, display_output=False,
 
         return self.veripb_return_code
 
-    def add(self, cpm_cons):
+    def add(self, cpm_cons, internal:bool=False):
         """
         Post a (list of) CPMpy constraints(=expressions) to the solver
         Note that we don't store the constraints in a cpm_model,
@@ -482,7 +482,8 @@ def add(self, cpm_cons):
         """
         # add new user vars to the set
                 # add new user vars to the set
-        get_variables(cpm_cons, collect=self.user_vars)
+        if not internal:
+            get_variables(cpm_cons, collect=self.user_vars)
 
         for con in self.transform(cpm_cons):
             cpm_expr = con
@@ -551,10 +552,10 @@ def add(self, cpm_cons):
                             # lt == x < y
                             # gt == x > y
                             lt_bool, gt_bool = boolvar(shape=2)
-                            self += (lhs < rhs) == lt_bool
-                            self += (lhs > rhs) == gt_bool
+                            self.add( (lhs < rhs) == lt_bool, internal=True )
+                            self.add( (lhs > rhs) == gt_bool, internal=True )
                             if fully_reify:
-                                self += (~bool_lhs).implies(lhs == rhs)
+                                self.add( (~bool_lhs).implies(lhs == rhs), internal=True )
                             self.gcs.post_or_reif(self.solver_vars([lt_bool, gt_bool]), reif_var, False)
                         else:
                             raise NotImplementedError("Not currently supported by Glasgow Constraint Solver API '{}' {}".format)
@@ -665,7 +666,7 @@ def add(self, cpm_cons):
             elif isinstance(cpm_expr, GlobalConstraint):
                 # GCS also has SmartTable, Regular Language Membership, Knapsack constraints
                 # which could be added in future. 
-                self += cpm_expr.decompose()  # assumes a decomposition exists...
+                self.add(cpm_expr.decompose(), internal=True)  # assumes a decomposition exists...
             else:
                 # Hopefully we don't end up here.
                 raise NotImplementedError(cpm_expr)